Q1 (Quantitative Comparisons)

While the results are suggestive that it does, they are somewhat limited and mostly limited to qualitative comparisons with little quantitative characterization of the methods efficacy.

In the band structure interpolation case, no quantitative errors/discrepancies are provided for the structure. For well separated (e.g., insulating) bands this should be easy to do. Is the method accurate on a fine k-mesh? Even for "entangled" bands this could be done sufficiently far from the largest eigenvalues considered. The plots (especially in Fig. 5) show clear differences, are these meaningful? could the accuracy of the method be improved? what is considered the "ground truth?"

Thank you for these insightful comments. We would like to address this concern by two results:

• Band Gap Comparisons

Our original goal for providing band structures was to demonstrate the capability of our program, however, the reviewer is correct that we have not included quantitative data and this has been an oversight. In the table below we compile the band gaps of diamond and silicon based on the valence band maximum and the conduction band minimum. It can clearly be seen that the band gaps of diamond are reasonably consistent, whereas the band gaps of silicon are not, which quantifies some of the error/discrepancies in bands observed by the reviewer.

The reviewer also identified that an analysis of accuracy with respect to the k-mesh and cut-off energy is absent. We note here that the k-mesh and maximum cut-off energy accessible to our method are currently modest due a need for code optimization. Consequently, calculations were conducted with a k-mesh of 1x1x1 and cut-off energy of 200/400/800 Ha as recorded in the captions of the band structures. These settings were reproduced in the other programs.

The ground truth would be the computed band gap at an infinitely high k-mesh and cut-off energy which should be consistent between all programs at this limit. This is usually approximated by increasing both settings until the band gap or desired observable is converged. Consequently, the primary means to improve the accuracy of our method is to increase the k-mesh and cut-off energy, which we expect to be able to do as we optimize our method.

• Energy differences of different crystal structures of carbon:

Additionally, we provide a table showcasing the energy differences between points A and B from Figure 3. This measure further assesses the accuracy of our computations compared to established methods.

These results indicate that our approach achieves comparable accuracy to conventional SCF methods, providing insight into the reliability of our computations.

Q2 (Comparison with Wannier Interpolation)

On a related note, a common technique to go from the k-mesh to band structure on some path in k-space is Wannier interpolation. Any comparison with this is missing from the present manuscript.

We acknowledge the reviewer's suggestion regarding comparing our approach with Wannier interpolation for deriving band structures along specific k-space paths. However, our manuscript's primary focus is to show the full differentiability of DFT for periodic systems, and comparing it with Wannier interpolation falls outside the intended scope of this work.

We appreciate the importance of such comparisons and will duly consider them in future research endeavors that align more closely with these specific comparative analyses.

Since I am not sure a notification is generated when the review is updated (I did not get one), this comment is to note that I have updated my review after reading the author response.

Q1 (Potential Applications)

I think the work would be more interesting and relevant if it had more discussion of how the work might synergize with other work, in particular materials generation and property prediction. Explicit explanation and examples of impact would be helfpul.

Thank you for the constructive comments. We agree that a more explicit discussion on the potential synergies of our work with other areas in materials science and computational chemistry would greatly enhance its relevance and impact. Here, we outline two examples of how our work could facilitate advanced analyses and the development of novel chemical models:

• Vibrational Analysis of Systems: Our method significantly simplifies the implementation of advanced vibrational analysis techniques. Typically, vibrational analyses use the harmonic approximation, which is limited to quadratic potentials and double derivatives. However, for a more comprehensive understanding, especially in anharmonic systems, higher-order derivatives are required. The inherent capability of our approach to handle such complex derivative calculations more efficiently makes it ideally suited for extended vibrational analyses, including anharmonic vibrations. By enabling easier implementation of methods requiring triple derivatives and beyond, our work could pave the way for more accurate and detailed vibrational studies of complex materials. Other examples of properties include polarisability, dipole, quadrupole and octupole moments, Raman spectroscopy methods, and so on.

• Development of Advanced Chemical Models: Many chemical models, particularly in density functional theory (DFT), rely heavily on derivatives of properties of the chemical system. Standard density functionals typically require the density, its gradient, and Laplacian, among other derivative quantities. Our work provides a robust framework that can facilitate the testing and development of new density functionals that require novel derivatives. The ease of implementing and testing these new functionals within our framework could significantly accelerate the advancement of DFT methods and lead to more accurate models for predicting material properties.

We will expand our discussion in the manuscript to explicitly outline these potential applications.

Q2 (Plane-wave Functions)

This paper may have just replaced these basis functions with the "plane wave" basis functions used in this paper. A switch which hardly seems like a big advance, though nevertheless a useful one.

We appreciate the comments and agree that differentiable DFT methods have been previously conceived in a molecular basis. However, we would like to emphasize that the implementation of these methods in a plane-wave basis, as presented in our work, is not a minor improvement but rather a significant advancement. The differences between molecular and plane-wave bases are substantial and necessitate a fundamentally different approach to the problem. Below, we detail some of these key differences:

1. Assumption of System Boundaries:

• Molecular Bases: These are typically used for finite chemical systems and utilize atom-centered basis functions, such as Gaussian type orbitals (GTOs). Implementing DFT in this context requires specifically written integral packages such as libint, tailored to these finite systems.

• Plane-Wave Bases: In contrast, plane-wave bases assume periodic boundary conditions, allowing us to access and model technically infinite crystal systems. This requires a significantly different mathematical framing, including the definition of energy space (k-space), integration techniques, and band structure analyses.

2. Technical Implementation and Mathematical Framing: The use of a plane-wave basis necessitates a distinct approach in defining energy space, integrating over this space, analyzing band structures, and considering electron occupation across energy bands. This complexity is distinct from the challenges faced in a molecular basis.

3. Efficiency and Accuracy Considerations: While molecular calculations can approximate crystal systems (and vice versa), such approaches are generally less desirable due to efficiency and accuracy concerns. Our method, by focusing on the plane-wave basis, addresses these challenges more directly and effectively for crystal systems.

4. Specialization in Computational Chemistry Modeling: The differences between these bases are so pronounced that most computational chemistry modeling tools specialize in either a molecular basis (e.g., Q-Chem, Gaussian, ORCA) or a plane-wave basis (e.g., VASP, Quantum Espresso). This specialization underscores the distinct challenges and approaches required for each basis type.

We will ensure that our manuscript more clearly articulates these points to emphasize the significant advancements our work represents in the context of differentiable DFT methods using a plane-wave basis.

Q3 (The Venue)

As such, I wonder if it would be more appropriate in a computational physics venue.

We value Reviewer 2's perspective and acknowledgment that our work aligns well with a computational physics forum. While many intended studies are indeed planned in that direction, leveraging this work as a foundation, we present this study as a novel tool for materials discovery rooted in machine learning infrastructure. Hence, we believe that ICLR serves as the appropriate venue for the debut of this work.

Furthermore, the significance of materials discovery within the machine learning community is continually expanding, as evidenced by notable contributions such as works [1-3]. We believe that members of the machine learning community are well-positioned to appreciate and make effective use of the additional functionalities embedded in the method we present. Although we do plan future endeavors targeting the computational physics community, these efforts have yet to materialize.

[1] Schütt, Kristof, et al. "Schnet: A continuous-filter convolutional neural network for modeling quantum interactions." Advances in neural information processing systems 30 (NeurIPS 2017).

[2] Cho, Youngwoo, et al. "Deep-DFT: A Physics-ML Hybrid Approach to Predict Molecular Energy using Transformer." Advances in Neural Information Processing Systems (NeurIPS) Workshop. 2021.

[3] Pope, Phillip, and David Jacobs. "Towards Combinatorial Generalization for Catalysts: A Kohn-Sham Charge-Density Approach." Advances in on Neural Information Processing Systems (NeurIPS 2023).

We would like to thank the reviewer for these insightful comments, which helped us improve the quality of this paper.

Q1 (Potential Applications)

If I am not mistaken, there is nothing much stopping this approach from being applied to other DFT settings, such as small molecules, using some form of Gaussian-harmonic basis set.

We appreciate the reviewer's comments into the potential application of our approach to other DFT settings, such as small molecules employing Gaussian-harmonic basis sets.

We would like to acknowledge that the prior work (Li et al., 2023), which is one of the most important foundation of this work, has indeed explored the extension of our methodology to diverse DFT settings, including applications involving small molecules utilizing Gaussian-harmonic basis sets. The results from this investigation supported the adaptability and effectiveness of our approach across a broader spectrum of electronic structure calculations.

For a comparison of plane-wave basis and other orthogonal basis, we advise the reviewer refer to Reviewer 2 Question 2.

Q2 (More Comparisons)

A lot of related work is mentioned, some of which close enough to allow for direct numerical comparisons. These would have improved the evaluation of the method.

Thank you for the comments. We would like to provide some numerical comparison results against existing methods.

1. Total energy differences.

We present the following tables of the difference in energy between points A and B from Figure 3. The result can measure the accuracy of our computation.

This results shows that our approach can achieve similar results as conventional SCF methods.

2. Band gap comparison.

In the table below we compile the band gaps of diamond and silicon based on the valence band maximum and the conduction band minimum. It can clearly be seen that the band gaps of diamond are reasonably consistent, whereas the band gaps of silicon are not, which quantifies some of the error/discrepancies in bands observed by the reviewer.

We note here that the k-mesh and maximum cut-off energy accessible to our method are currently modest due a need for code optimization. Consequently, calculations were conducted with a k-mesh of 1x1x1 and cut-off energy of 200/400/800 Ha as recorded in the captions of the band structures. These settings were reproduced in the other programs.

We would like to thank the reviewer for these insightful comments, which helped us improve the quality of this paper.

Q1 (The Advantages of Direct Optimization over SCF)

The claim that the iterative nature of SCF obstructs the application of deep learning is not very suitable.

We appreciate the reviewer's perspective regarding the iterative nature of gradient-based optimization and its potential compatibility with deep learning methods.

We want to emphasize that SCF methods can fail to find solutions in some cases (for example near-degeneracies or highly correlated electrons [1]). Direct optimization in our method can avoid these failures.

We also wish to emphasize that the differentiability of our method provides more mechanisms to interface with deep learning architectures than existing SCF methods. For example, the differentiable KS-DFT method proposed in [2] involves a scf loop in the computational graph, where the gradient backpropagation goes through all the scf iterations and which requires a fully converged solution to the scf loop and risks the failure mentioned above. In contrast, the differentiable optimization enables the development of end-to-end trainable models without going back and forth between SCF in another package, and neural networks in deep learning framework.

We apologize for the poor phrasing, and have modified our manuscript to make these points more clear.

[1] Schlegel, H. Bernhard, and J. J. W. McDouall. "Do you have SCF stability and convergence problems?." Computational advances in organic chemistry: Molecular structure and reactivity (1991): 167-185. [2] Kasim, Muhammad F., and Sam M. Vinko. "Learning the exchange-correlation functional from nature with fully differentiable density functional theory." Physical Review Letters 127.12 (2021): 126403.

Q2 (Band Structure)

For the band structure experiment, although several popular packages are compared, the advantages of the proposed method are not clear.

Thank you for this comment. We would like to emphasize that the band structure experiments are a standard application for solid-state DFT computation. The purpose of this experiment is to show the efficacy of our direct optimization approach, which can produce comparable results to SCF methods. Here we are not claiming that we outperformed these popular packages. However, there are some potential advantages that our method might possess due to its full differentiability nature. For instance, we should be able to design pseudo potentials or exchange-correlation functionals with guidance of the band gap, which is more easy to access from experiments and no existing methods have done so, to the best of our knowledge. This is an intriguing application that we will explore further.

For more quantitative comparison of band structure, please refer to Reviewer 4 Question 2.

Q3 (Geometry Optimization)

For the geometry optimization experiment, no baseline methods are compared, while there are mature automation tools such as geomeTRIC and pyberny. …traditionally there can be analytical derivatives available. So I wonder if there is any advantage of using the proposed method?

Thank you for your insightful comments and we acknowledge the importance and widespread use of analytical derivatives in traditional SCF methods for geometry optimization.

The key advantage of our approach is its inherent adaptability and flexibility. By integrating the computation of derivatives directly into the optimization process, our method simplifies the workflow and reduces the complexity of implementing geometry optimization, especially in systems where deriving analytical expressions is challenging or impractical.

Furthermore, our method's capacity to handle derivative calculations internally lends itself to a broader range of applications. This includes complex systems where traditional methods might struggle due to the intricate nature of deriving and implementing analytical derivatives. As detailed in our response to Reviewer 2 Question 1, we foresee numerous potential applications in future work, such as handling complex molecular systems, materials with unconventional properties, and scenarios where rapid prototyping of computational models is essential.

In light of your feedback, we recognize the necessity to clarify these points in our manuscript. We will include a more detailed discussion of the advantages of our method compared to traditional approaches, backed by quantitative analyses where applicable.