Efficient Approximations of Complete Interatomic Potentials for Crystal Property Prediction
Authors: Yuchao Lin
Keqiang Yan
Youzhi Luo
Yi Liu
Xiaoning Qian
Shuiwang Ji
Venue: Submitted to ICLR 2023 (2023)
Keywords: graph neural network
material property prediction
crystal property prediction
crystal structure modeling
interatomic potential

Abstract: We study the problem of crystal material property prediction. A crystal structure consists of a minimal unit cell that is repeated infinitely in 3D space. How to accurately represent such repetitive structures in machine learning models remains unresolved. Current methods construct graphs by establishing edges only between nearby nodes, thereby failing to faithfully capture infinite repeating patterns and distant interatomic interactions. In this work, we propose several innovations to overcome these limitations. First, we propose to model physics-principled interatomic potentials directly instead of only using distances as in existing methods. These potentials include the Coulomb potential, London dispersion potential, and Pauli repulsion potential. Second, we propose to model the complete set of potentials among all atoms, instead of only between nearby atoms as in prior methods. This is enabled by our approximations of infinite potential summations with provable error bounds. We further develop efficient algorithms to compute the approximations. Finally, we propose to incorporate our computations of complete interatomic potentials into message passing neural networks for representation learning. We perform experiments on the JARVIS and Materials Project benchmarks for evaluation. Results show that the use of complete interatomic potentials leads to consistent performance improvements with reasonable computational costs.

Under review as a conference paper at ICLR 2023
EFFICIENT APPROXIMATIONS OF COMPLETE INTER- ATOMIC POTENTIALS FOR CRYSTAL PROPERTY PRE- DICTION
Anonymous authors Paper under double-bling review
ABSTRACT
We study the problem of crystal material property prediction. A crystal structure consists minimal unit cell that is repeated infinitely in 3D space. How to ac- curately represent such repetitive structures in machine learning models remains unresolved. Current methods construct graphs by establishing edges only between nearby nodes, thereby failing to faithfully capture infinite repeating patterns and distant interatomic interactions. In this work, we propose several innovations to overcome these limitations. First, we propose to model physics- ppriniiepr inter- atomic potentials directly instead of only using distances as in existing methods. These potentials include the Coulomb potential, London dispersion potential, and Pauli repulsion potential. Second, we propose to model the complete set of poten- tials among all atoms, instead of only between nearby atoms as in prior methods. This is enabled by our approximations of infinite potential summations with prov- able error bounds. We further develop efficient algorithms to compute the approx- imations. Finally, we propose to incorporate our computations of complete inter- atomic potentials into message passing neural networks for representation learn- ing. We perform experiments for evaluation. Results show that the use of complete interatomic potentials leads consistent performance improvements with reasonable computational costs.
the JARVIS and Materials Project benchmarks
INTRODUCTION
The past decade has witnessed surge of interests and rapid developments in machine learning for molecular analysis (Duvenaud et al. 2015). These initial studies mainly focus on the prediction and generation problems of small molecules. To enable computational analyses, molecules need to be featurized in appropriate mathematical representation form. Recently, with the advances of graph neural networks (GNNs) (Gilmer et al.] 2017 Battaglia lia et al., 2018), molecules are more commonly represented graphs in which each node corresponds a chemical bond (Stokes et al., 2020 Wang et al., 2022ch. A variety molecular graph predic- tion (Stokes et al.) 2020 Wang al. 2022c) and generation (Shi et al. 2019: [Jin et al.) 2018 Luo let al., 2021) methods have been developed based on 2D molecular graph representations. A key lim- itation of the 2D graph representations is that the 3D geometries of molecules are not captured, but such information may be critical in many molecular property prediction problems (Hu et al. 202 To enable the encoding of 3D molecular geometries in GNNs, series of 3D GNN methods have been developed for prediction (Schütt et al. 2017 Gasteiger et al., 2019 Liu et al. 2022b; Wang al. (022b) and generation (Liu 20.22 Luo & Ji 2022 Hoogeboom al. (022) problems. In these 3D graph representations, each node is associated with the corresponding atom's coordinate in 3D space. Geometric information, such as distances between nodes and angles between edges, used during message passing in GNNs. Recently, these methods have been extended to learn representations for proteins (Jing et al. 2020 Wang et al. 2022a).
atom, and each edge corresponds to
Inspired by the success of GNNs on small molecules, Xie & Grossman 2018) developed the crystal graph convolutional neural network (CGCNN) for crystal material property prediction. Different from small molecules and proteins, crystal materials are typically modeled by minimal unit cell (similar to a small molecule) that is repeated in 3D space with certain directions and step sizes. In theory, the unit cell is repeated infinitely in 3D space, but any real-world material has finite size.

---

Under review as a conference paper at ICLR 2023
However, given that our modeling is at the atomic level, modeling crystal materials as infinite repeti- tions of unit cells is approximately accurate. Therefore, a key challenge in crystal material modeling is how to accurately capture the infinit-raange interatomic interactions resulted from the repetitions of unit cells in 3D space. Current GNN-based crystal property prediction methods construct graphs by creating edges only between atoms within a pre-specified distance threshold (Xie & Grossman 2018 Chen ct al. 2019; et al., 2020) Schmidd et al. [202]] Choudhary & 1021), Thus, they fail to capture interactions between distant atoms explicitly.
In this work, we propose a new graph deep learning method, PotNet, with several innovations to significantly advance the field of crystal material modeling. First, we propose to model interatomic potentials directly as edge features in PotNet, instead of using distance as in prior methods. These potentials include the Coulomb potential (West 1988). London dispersion potential (Wagner & Schreiner 2015), and Pauli repulsion potential (Krane 1991). Second, a distinguishing feature of PotNet is to model the complete set of potentials among all atoms, instead of only between nearby atoms in prior methods. This is enabled by our approximations of infinite potential summations with provable error bounds. We further develop efficient algorithms to compute the approximations. Finally, we propose incorporate our computations of complete interatomic potentials into message passing neural networks for representation learning. We performed comprehensive experiments on the JARVIS and Materials Project benchmarks to evaluate our methods. Results show that the use of complete interatomic potentials in our methods leads to consistent performance improvements with reasonable computational costs.
BACKGROUND AND RELATED WORK
2.1 CRYSTAL REPRESENTATION AND PROPERTY PREDICTION
A crystal structure can be represented as periodic repetitions of unit cells in the three- dimensional (3D) Euclidean space, where the unit cell contains the smallest repeatable structure of a given crystal. Specifically, let n be the number of atoms in the unit cell, crystal can be rep- resented as M = (A, L). Here,
= {a:}"= = {( p;) describes one of the unit cell structures of M. where æ, E Rº and Pt E R denote the b-dimensional feature vector and the 3D Cartesian coordinates of the i-th atom in the unit cell, respectively. L [1, 12, 13] E R333 is the lattice matrix describing how a unit cell repeats itself in the 3D space. In the complete crystal structure, every atom in unit cell repeats itself periodically in the 3D space. Specifically, from an arbitrary integer vector & E 2 and the unit cell structure A, we can always obtain another repeated
unit cell structure A {af Jam = {<ff pt ) where x = 2j p Pi + Lk. Hence, the complete crystal structure À of M with all unit cells can be described as
Ā =
A
kEZ
In this work, we study the problem of crystal property prediction. Our objective is to learn a property prediction model f: M + y E R that can predict the property y of the given crystal structure M. We will focus on predicting the total energy, or other energy-related properties of crystals.
2.2 CRYSTAL PROPERTY PREDICTION WITH INTERATOMIC POTENTIALS
Most of the classical crystal energy prediction methods are based on interatomic potentials. Ac- cording to the studies in physics (West) 1988 Daw et al.) 1993; Brown 2016), the total energy of crystal structure can be approximated by the summation of interatomic potentials in the crystal. Particularly, the three following categories of interatomic potentials are widely used in erystals, and they considered sufficient for accurate energy approximation.
Coulomb potential is caused by the electrostatic interaction of two atoms with charges. Coulomb potentials are closely related to ionic bonding and metallic bonding in crystals (West) 1988). For any two atoms a and b, let zn and on denote the number of charges in the atom and b. and let d(a, b) be the Euclidean distance between the atom a and b. The Coulomb potential V Coutomb is
defined as V Coulomb (a, b) In)]d(a, Here is the elementary charge constant, and En is the permittivity constant of free space.
in =n=p"
2

---

Under review as a conference paper at ICLR 2023
London dispersion potential describes the Van der Waals interaction between atoms. London dispersion potential is often considered in energy estimation because its contribution is cumula- tive over the volume of crystals (Wagner & Schreiner, 2015) and sometimes very strong in bulk crystals, such as the crystals containing sulfur and phosphorus. The mathematical form of this potential can be described as V London b) e//ddaa, b). where is a hyperparameter.
Pauli repulsion potential is resulted from the Pauli exclusion principle generally exists in all crystal structures. The Pauli exclusion principle forces any two atoms to be sufficiently far away from each other so that the electron orbits of them do not overlap. Such exclusion interactions lead Pauli repulsion potential with the form of V Pauli (a, b) cdd(a.b) where or is hyperpa- rameter (Buckingham
1938 Slater 1928).
2.3 CRYSTAL PROPERTY PREDICTION WITH DEEP LEARNING
While physics- based methods have been used for predicting the crystal energy for long time, these methods are usually crystal-specific, i.e., one method can only achieve accurate approximation for one specific type of crystals. Recently, thanks to the advances of deep learning, many studies have been done develop a general crystal property predictor for a variety of different crystals with pow- erful deep neural network models. Some studies (Wang et al. 20211] Jha et al., 2018 2019] Goodall & Lee]|222 2020) represent crystals as chemical formulas, and adopt sequence models to predict proper- ties from these string representations. However, more recent studies consider crystals as 3D graphs and employ expressive 3D GNNs (Schütt et al. 2017 Klicpera et al., 2020b] Liu et al. 2022b). a family of deep neural networks specifically designed for 3D graph-structured data, to crystal repre- sentation learning. CGCNN (Xie & Grossman) 2018) is the first method that proposes to represent crystals with radius graphs and adopts a graph convolutional network to predict the property from the graph. Based on the pioneering exploration of CGCNN, a lot of subsequent studies (Schmidt al. 2021 Louis et al. 2020 Chen al., 2019; Choudhary & DeCost [021] Batzner et al. 2022) propose various 3D GNN architectures to achieve more effective crystal representation learning. Particularly, by enhancing the input features with angle information, ALIGNN (Choudhary De- Cost 2021) develops the currently most powerful 3D GNN architecture for crystals and achieves the
best crystal property prediction performance.
3 METHOD
Although existing GNN-based methods have achieved impressive performance in crystal property prediction, they struggle in further boosting the performance due the approximation of interatomic interactions using functional expansions based on distances and failing in capturing complete inter- atomic interactions. In this section, we present PotNet, a novel crystal representation model that can overcome these limitations of prior methods. Based on the physical modeling of crystal energy, Pot- Net explicitly uses infinite potential summations as input features to capture complete interatomic interactions. The infinite potential summations are incorporated into the message passing mecha- nism of graph neural networks and efficiently approximated by a fast algorithm. To the best of our knowledge, PotNet is the first work that bridges the classical crystal energy computation methods based on potentials and the data-driven methods based on deep neural networks.
3.1 APPROXIMATING CRYSTAL ENERGY WITH COMPLETE INTERATOMIC POTENTIALS
According to the density functional theory (DFT) in physics, for any crystal M = (A, L) with the complete structure À defined in Eqn. d1 its total energy E(M) can be accurately approximated by the embedded atom method (Daw & Baskes 1984 Daw et al. 1993 Baskes 1987 [ee et al. 2016
Riffe et al., 2018) in the form of
E(M) = Σ Σ ac D
Σ
V(a. b) + Σ F(pa),
(2)
2
acA
where V(a, b) denotes the interatomic potentials between the atoms and b, capturing the magni- tude of interactions; Pa is the local electron density of the atom a. determined by the coordinate and number of charges of the atom according to the HlhenbeeggKKohh theorem; F(.) is a parametrized
3

---

Under review as a conference paper at ICLR 2023
function to embed the electron density Pa. Actually, existing studies (Jakkanen & Müser 2013 show that Pa can be considered as function of 2 bbya, (a, b) mathematically. Hence, Eqn.
can
be rewritten in the following form:
B(M) Σ
Σ ()
Σ byAaEE
v(a,
aEA bya.be
where G(.) is a parametrized function. Eqn. 3 can be considered as a way to compute the energy from the complete interatomic potential summation Ebza, (a, b) of every atom in the unit cell A. However, in practice, the function G is computationally expensive if not infeasible. Hence, more and more studies have turned to the powerful learning capability of modern deep neural net- work models to approximate it effectively.
3.2 LIMITATIONS OF EXISTING DEEP LEARNING METHODS
Currently, most of the existing graph deep learning methods for crystals (Xie & Grossman 2018
Chen al.) 2019: Louis al., 2020: Choudhary & DeCost 202T) use radius graph representations and distance- based features as inputs to predict the crystal energy in Eqn. ( Specifically, for crystal M (A, L), the radius graph is constructed by adding edges between any atom a in the unit cell A and any other atom b in the complete crystal structure A whose distances are smaller than pre-specified distance threshold r. In addition, some functional expansions of distances, e.g., radial basis functions (RBF), are used to model interatomic interactions and form the input edge features to 3D GNN models. Hence, let a aFa Pa), (æb. Pt), the crystal energy prediction É(M) ) of these methods can be generally described as
E(M) = Σ Σ a DEN
H (IPpa - Pbll2)),
(4)
where (a(a)= (b: b a, b E A, Ilpa poll2 < <r), o(-) denotes the functional expansions, and H(+) is non-linear function based on 3D GNN models.
However, we argue that predicting or approximating the energy with Eqn. 44 is a suboptimal so- lution. Actually, compared with Eqn. (3). which is physics-principled, there exist non-negligible approximation eтrors in Eqn. 44). First, Eqn. (4) captures the interatomic interactions based on interatomic distances, while the energy can more accurately approximated by function of in- teratomic potentials as in Eqn. (3. Though according to Sec. 2e2. interatomic potentials themselves are also functions of distances, we argue that directly using functional expansions of distances not the best solution to crystal energy prediction. The commonly used functional expansions in existing methods, such as RBF (6.), have different mathematical forms from potentials defined in Sec. 2.2 Intuitively, this poses more challenges to 3D GNN models since they need to learn mapping from o(.) to the energy E, while the energy E is not a direct function of oc.). Hence, we argue that di- rectly employing the physics-principled potential functions instead of ((.) as input features is more suitable for crystal energy prediction.
Second, different from Eqn. (3), Eqn. (4) does not capture the complete set of interatomic inter- actions because the summation set (,(() (a) of atoms b is constrained to be the atoms whose dis- tances the atom ignoring long-range interatomic potentials, i.e., interatomic potentials between distant atoms. Dif- ferent from molecules with finite structures, long-range interatomic potentials cannot be ignored for crystals with infinite structures. By the first principles in physics, interatomic potentials decay algebraically when pairwise interatomic distances become larger. Hence, for a finite structure like molecules, the potentials from atoms that are far away from a given atom is limited and can be ig- nored. However, long-range interatomic potentials have a significant influence on a given atom in the infinite crystal structure. Taking Coulomb potentials as example, assuming that we are given a ID crystal structure where there is only one atom repeating itself with Euclidean distance 1, each atom has only unit of charge, and the total energy is simply the sum of all interatomic Coulomb potentials. As defined in Sec. 2.2 Coulomb potential energy V(a b) between atoms a and b satisfies V(a, b) x 1/d, where d is the distance between atoms and b. Considering the
smaller than This can lead to a significant approximation error due to
4

---

Under review as a conference paper at ICLR 2023
Coulomb potentials between given atom and all other atoms, the total potential V of them sat- isfies V x - 1/n. If only the pairwise atoms within the distance threshold considered, an infinite error of energy calculation is introduced. Specifically, the smallest possible prediction error AV satisfies AV x
[n+1r+1j 1/n, which is infinite. In other words, ignoring interatomic Coulomb potentials between atoms with distances larger than will cause significant prediction error of the total energy. We can observe from this example that the failure to capture complete interatomic potentials due the use of radius graphs is a key factor that prevents accurate energy prediction in existing GNN-based methods.
3.3 MESSAGE PASSING WITH COMPLETE INTERATOMIC POTENTIALS
follows from the analysis in Sec. 3.2 that major limitations existing deep learning methods for crystal representation learning lie in (I) not making predictions from physics-principled interatomic potentials, and (2) not considering complete interatomic interactions. To overcome these limitations, propose explicitly
complete interatomic potential summations in GNN models. Since proposed method is tightly related to potentials, we name it PotNet.
By reformulating Eqn. (3), PotNet incorporates the crystal energy computation with complete interatomic potentials into the message passing scheme of GNN models. For any material structure
M = L), we can rewrite the definition of its complete structure A in Eqn. L as
Ú- U
= U
U (b*)
U Ao.
(5)
KEZ
KEZ BE A
BE KEZS
b
where As UKEzN {b*} denotes the of atoms containing the atom b from the unit cell A and all its periodically repeated duplicates in the complete crystal structure. With Eqn. (5) we can reformulate Eqn. (3) as
Σ Σ v(a.c) BEA bEA AD..
B(M) Σ
Σ c(a,c)+G
a
bEA Ao
(6)
Σ acA
((a,b(()))))) +G
Σ Es(aab))
bEA
bEA
where the infinite potential summation S(a, b) -D)==Effaacc V(a, c) denotes the sum of the inter- atomic potentials from the atom together with its all periodic duplicates the atom a. Eqn. (6 can be integrated into the message passing scheme of GNN models. Specifically, G for M = (A, L). where each atom the unit cell A corresponds to a node in the graph. For any two nodes u, in the graph, there is edge connecting them, and every node in the graph is also connected to itself by a self-loop edge. we consider the infinite potential summation S(a, b) the feature of the edge from node b to a, we can use the message passing based non-linear neural network model in GNN to fit the function 2 Ebe A S(a, b) + G (E S(a, b)).
create graph
Based on this design of directly using interatomic potentials as edge features, a GNN model with multiple message passing layers on the graph G' to predict the crystal energy of M. The computational process of the E-th message passing layer for the node a can be described as
PotNet employs
h = U2
h
(f-1)
Σ So h
(f-1)
(8-1) hb
S(a.b))
bEA
DEA
where (1) denotes the embedding vector of node a generated from the f-th message passing layer, (0) is initialized to the atom feature vector of the atom a, and (,) fo(.) are both neural net- work models with trainable parameters ş and 6, respectively. Here, the model fo plays the role of capturing information from both atomic features and complete interatomic potentials. Detailed information about model architectures of fe and g is provided in Appendix D.I Note that our Pot- Net is actually 3D GNN model even though 3D geometric information is not explicitly involved in Eqn. (7). This is because the edge feature S(a, b) is related to potential functions, and by Sec. 2.2 know that they are computed from interatomic distances. In other words, PotNet can be consid- ered to encode 3D geometric information with potential functions, though our direct motivation of using potential functions comes from the physical modeling of crystal energy.
5

---

Under review as a conference paper at ICLR 2023
Σ
V(
S(w. v?)
Σ MHA.NOU
V(w.u)
((w.v))
(((,a)
u
S(v. v)
((v.a))
V(w.a)
(s)
(b)
(c)
Figure 1: Schematic illustrations of how complete interatomic itneractions captured in PotNet. Note that PotNet models 3D crystals while we have 2D illustration for simplicity. (a) An example crystal in which each unit cell contains of atoms are captured. For simplicity, we only show the potentials from all v atoms to a u atom. (h) The complete set of potentials (a) can be grouped into four categories, including + v. + u, + and v + u. (c) We propose to compute an approximate summation for each category of potentials.
atoms u and v. In PotNet, the potentials between all pairs
Intuitively, the message passing process Eqn. 77 over the graph G can be considered general case of employing a radius graph where the distance threshold r goes to infinity, i.e., r + +oo. In this case, as shown Fig. D a), for any atom in the crystal, all the other atoms in the complete crystal structure have been included to interact with it. If we follow the radius graph construction process in the previous methods (Xie & Grossman 2018 Chen et al.) 2019 Louis et al. 2020 Choudhary & DeCost 2021 , we obtain graph G in which there exist an infinite number of edges between every pair of nodes. However, PotNet simplifies this complicated graph & to the graph G' in which only one edge exists between every node pair. Specifically, PotNet directly models interatomic interactions as potentials and for any two nodes G, PotNet aggregates all edges between them to single edge by the use of infinite potential summation S(a, b) (sce Fig. [b)). In other words, PotNet provides an effective solution that enables GNN models to capture complete interatomic interactions through the use of infinite potential summations.
3.4 EFFICIENT COMPUTATION OF INFINITE POTENTIAL SUMMATION
Although we have effectively incorporated infinite potential summations into the message passing based GNN models, the computation of these infinite potential summations is not easy. Basically, there are two challenges to achieve accurate and efficient computation of the infinite potential sum- mations. For accuracy, the computation algorithm needs to have provable error bounds. For effi- ciency, it needs be fast for scalable GNN training and fast crystal property prediction. To tackle these two challenges, based on the Ewald summation method (Ewald) 1921). Specifically, we unify the summations of three infinite potentials between atom and all duplicated positions of atom b into integral form such that the Ewald summation method can be applied for efficient implementation in Pot- Net (Fig. [()). The key idea of the Ewald summation is that a slowly converging summation in the real space is guaranteed to be converted into a quickly converging summation in the Fourier space (Woodward 2014). Based on that, the Ewald summation method divides summation into two parts. One part has quicker converging rate in the real space than the original summation. The other filower-to-converee
derive fast approximation algorithm for infinite potential summations
part is transformed into the Fourier space and becomes quickly conver- gent. In our method, the Ewald summation method is used with the infinite summations by dividing the integral into two parts, including one part that converges fast in the Fourier space and another part that converges fast in the real space, to obtain fast approximations with provable error bounds.
To apply the fast approximation algorithm of infinite summations proposed by Ewald (921). a unified integral view of infinite potential summations is needed. Following notations Sec. and 3.3 we denote the positions of the atoms in the set Ab as P. {pÊ p = Ps + Lk, & z The Euclidean distances between the atom and all atoms in Ab can be represented as
{d = |p, Lk Pall: k 2 .. As described in Sec. 3.3] we simplify the Coulomb potential function as V Coutomb (a, b) -e1/ e/d((a b), where is hyperparameter scaling the potential. Based

---

Under review as a conference paper at ICLR 2023
Table I: Comparison between PotNet and other baselines in terms of test MAE on the Materials Project dataset. To make the comparison clear and fair, we retrain baseline methods using the same dataset settings. We also show results from original papers in parentheses and mark missing results as The best results are shown in bold and the second best results are shown with underlines.
Formation Energy
Band Gap
Bulk Moduli
Shear Moduli
Method
eVaatom
eV
log(GPa)
log(GPa)
CGCNN
0.031 (0.039)
0.292 (0.388)
0.047 (0.054)
0.077 (0.087)
SchNet
0.033 (0.035)
0.345 (-)
0.066 (-)
0.099 (-)
MEGNET
0.030 (0.028)
0.307 (0.33)
0.051 (0.050)
0.088 (0.079)
GATGNN
0.033 (0.039)
0.280 (0.31)
0.045
0.075
ALIGNN
0.0221
0.0218
0.051 (-)
0.078 (-)
PotNet
0.0196
0.0204
0.042
0.069
on that, we can represent Coulomb potentials from all atoms in As to the atom a as { [|p6 + Lk Pall.k k ) Similarly, London dispersion potentials from all atoms in As to the atom can be represented as { fo d | 1 pall, A Z²}. ). It is worth noting that Coulomb potentials and London dispersion potentials can be represented in unified view { connttant d
Vaab + Lk]|, Vab Ph Pa: k- }²), where p is a positive number. And we represent Pauli potentials from all positions of atom to atom i as {e-ad d = Pb + Lk Pall, k E Z²}. We provide detailed proofs in Appendix C.I that the summations of these three potentials can be unified in an integral form as
S(a, b) -
Σ kE?
C- A Pdt,
where A, B, C, D are constants derived from the corresponding specific potential forms. We then apply the Ewald summation method (Ewald 1921) to Equ. (8) and split it into two parts as
S(a, b) -D Σ [ (C----ArLAAttssttt
Σ 2E23
CG-le-AnlLktwad 2t- dt
KE23
=S (a, b) + Sdiner(a, b),
where Sairect denotes the part that converges fast in real space, and S denotes the other part that quickly converges in Fourier space when the total summation converges as shown by Ewald (192T). Based on this, we further show that Sarce and Spourier can be expressed as summations of incomplete Bessel functions (z, y) in Appendix C.2 and the approximation error can be analyzed and shown to be bounded. Note that for London dispersion potentials and Pauli potentials, the transformed summations of incomplete Bessel functions can be directly approximated. However, for Coulomb potentials, the direct summation diverges. Concretely, when p > for potentials of Pa, k E Z3), the corresponding potential summation
d= Vab + Lk]], Pab ps Pb-Paa the corresponding potential summation diverges. To tackle this problem, we follow previous mathematical derivations (Harris) 2008 Slevinsky & Safouhi] 2010 Jones 2007)
form constant
di
converges, and when
and use the analytically continued incomplete Bessel functions to approximate infinite summations
of Coulomb potentials as shown in Appendix C We then provide detailed mathematical proofss that the summations of incomplete Bessel functions are convergent and can approximated with an error bounded by the Gaussian Lattice Sum in Appendix B.2 Detailed implementation of the proposed summation algorithm can be found in Appendix C.4 It is worth noting that PotNet is the first method to use the incomplete Bessel function to compute the potential summations and is able to compute the Pauli potential summation, while previous methods (Crandall 1998 Lee & Cail 2009 Nestler et al. 2015) cannot achieve this. Also, are able to compute other interatomic potential summations including Lennard-Jone Jones potential, Morse potential, and screened Coulomb potential by using method as shown in Appendix

---

Under review as a conference paper at ICLR 2023
Table Comparison between PotNet and other baselines in terms of test MAE on JARVIS dataset. The best results are shown in bold and the second best results are shown with underlines.
Formation Energy
Bandgap(OPPT)
Total energy
Bandgap(MBBJ)
Ehull
Method
eVaatom
eV
eVaatom
eV
eV
CFID
0.14
0.30
0.24
0.53
0.22
CGCNN
0.063
0.20
0.078
0.41
0.17
SchNet
0.045
0.19
0.047
0.43
0.14
MEGNET
0.047
0.145
0.058
0.34
0.084
GATGNN
0.047
0.17
0.056
0.51
0.12
ALIGNN
0,0331
0.142
0.037
0.31
0.076
PotNet
0.0308
0.136
0.034
0.28
0.050
EXPERIMENTAL STUDIES
4.1 EXPERIMENTAL SETUP
We conduct experiments on two material benchmark datasets, including The Materials Project and JARVIS. Baseline methods include CFID (Choudhary et al. 2018). SchNet (Schütt et al. 2017
CGCNN (Xie Grossman 2018), MEGNET (Chen et al. 2019). GATGNN (Louis et al. and ALIGNN (Choudhary & DeCost 2021J. Unless otherwise specified, we report the results ported by referred papers or provided by original authors. All PotNet models are trained using the Adam (Kingma & Ba 2014) optimizer with weight decay (Loshchilov & Hutter 2017) and one cycle learning rate scheduler (Smith & Toppn 2019) with learning rate of 0.001, training epoch of 500, and batch size of 64. We use Pytorch to implement
2020
models. For all tasks two benchmark datasets, we use one NVIDIA RTX A6000 48G GPU for computing, Other detailed configurations of PotNet for different tasks provided in Appendix D.I
To capture the global infinite-range interactions without losing details of local interactions, PotNet uses both local and infinite crystal graphs. Specifically, for the local crystal graph, we use the radius crystal graph proposed by CGCNN but replace Euclidean distances with interatomic potentials for edge features. be concrete, because the influences of London dispersion potentials and Pauli potentials are limited and can be ignored when only considering nearby regions, we only use the Coulomb potentials the radius crystal graph. The infinite crystal graph is constructed as described in Sec. 3.3 where Coulomb potentials, London dispersion potentials and Pauli potentials are used.
4.2 EXPERIMENTAL RESULTS
The Materials Project. We first evaluate PotNet on The Materials Project-2018.6.1, which is widely used large-scale material benchmark with 69239 crystals. We follow previous works (Xie & Grossman 2018 Chen al.) 2019 Choudhary & DeCost 2021 Louis al. 2020) and use four crystal properties including formation energy, band gap, bulk moduli, and shear moduli. We notice that previous works (Xie & Grossman 2018 Chen et al. 2019 Choudhary & DeCost. 2021 Louis ett al. 2020 [Scüüt et al. (017) compare with each other using different versions of splitting train- ing, evaluation, and testing datasets with different random seeds. For instance, the original CGCNN paper only uses 28046 training samples for formation energy prediction, resulting in the original result of 0.039 as shown Table 7 To make the comparisons fair, we follow the settings of the pre- vious satte-of-the-ass energy and band gap prediction, and retrain all other baselines using the same dataset setting for these two tasks. For bulk and shear moduli, we follow the dataset setting of GATGNN (Louis al.) 2020) which has the best prediction performances for these two tasks and retrained other base- lines. Detailed configurations are also shown in Appendix D.2 We present our results in Table where PotNet consistently outperforms other SOTA methods on all four tasks by large margins.
(SOTA) ALIGNN (Choudhary & DeCost 2021) for tasks including formation
JARVIS Dataset. We then evaluate PotNet on JARVIS, newly released benchmark dataset pro- posed by Choudhary et al. (2020) with 55722 crystals. We evaluate PotNet on five crystal prop- erty prediction tasks including formation energy, bandgap (OPT), bandgap (MBJ), total energy, and Ehull. We follow ALIGNN and use the same training. validation, and test sets for these tasks. Since there missing results for baseline methods, we retrain corresponding baseline methods follow-
8

---

Under review as a conference paper at ICLR 2023
ing the same settings with ALIGNN discussed in Appendix D.2 As shown in Table achieves the best performances on all five tasks consistently by significant margins.
2 PotNet
Efficiency of PotNet. Beyond the superior modeling capacity for crystals, our PotNet faster and more efficient than ALIGNN. To demonstrate the efficiency of PotNet, compare PotNet with ALIGNN in terms of training time per epoch, total training time, and inference time for the task of JARVIS formation energy prediction. From Table PotNet four times faster in terms of total training time and inference time compared with ALIGNN.
Table 3: Model complexity and runtime compared with ALIGNN on JARVIS formation energy.
Method
4Params
Time/Epoch
Total Training Time
Total Testing Time
ALIGNN
15.4 MB
327
27
156x
PotNet
7,9MB
48
6.7h
41x
We also analyze the time cost of infi- nite summations as shown in Table We calculate the infinite summations during preprocessing. Unlike previ- ous methods including ALIGNN, we need to calculate the complete poten- tial set besides constructing graphs. Therefore, we do require more com-
Table 4: Preprocessing time compared with ALIGNN on the JARVIS dataset. The first and second columns show the preprocessing time on the whole JARVIS dataset with 55722 crystals. The third column denotes the mean inference time considering both preprocessing and model time.
Method
Preprocessig Timo/Crystal
Total Preprocessing Time
Inference/Crystall
AI KINN PoiNet
2.7 ms
152
34 ms
10.8
603
puting time for preprocessing. However, for a single material, the preprocessing time of our method is at the level of milliseconds. Even if additionally compute the infinite potential summations, the preprocessing time of our method and ALIGNN is still within the same order of magnitude. If considering both preprocessing time and model inference time for single material screening, we have a faster inference speed than ALIGNN as illustrated in the table. Overall, the computational cost of our method reasonable.
4.3 ABLATION STUDIES
In this section, we demonstrate the importance of two core components of PotNet, including in- teraction modeling using potentials and infinite summation of potentials for crystal prediction. We conduct experiments on the JARVIS formation energy task, and use test MAE as the evaluation metric.
Table 5: Ablation studies for the effects of adding
Interaction Modeling using Potentials. We demonstrate the importance of interaction mod- eling using potentials by replacing potentials with Euclidean distances used by previous works in our PotNet with exactly the same model architecture. Specifically, we denote PotNet with only local crystal graph as the base model. We use 'Base + Euclidean" to represent the base model with Euclidean distances as edge features and 'Base + Potential' to represent the base model using Coulomb potentials as edge features. As shown in Table 5 by replacing Eu- clidean distances with Coulomb potentials, PotNet without considering infinite potential summation already obtains significant performance gain from 0.0363 to 0.0318, revealing the importance of interaction modeling using potentials in PotNet.
Coulomb potentials and infinite summation.
Method
JARVIS Formation Energy
Base Euclidean
0.0363
Base Potential
0.0318
Base Potential + Infinite
0.0308
Infinite Summation of Potentials. The importance of infinite summation of potentials is demon- strated by comparing the previous hase models with 'Base + Potential + Infinite', denoting the full PotNet model with infinite summation in infinite crystal graph. It can be seen from Table (5) that by using infinite crystal graphs introduced in Sec. 3.3 the global repeating patterns of crystal structures are captured, resulting in performance gain from 0.0318 to 0.0308 for formation energy prediction.
5 CONCLUSION
We study the problem of how to capture infinite-range interatomic potentials in crystal property prediction directly. As radical departure from prior methods that only consider nearby atoms, we develop a new GNN, PotNet, with the message passing scheme that takes efficient approximations to capture the complete set of potentials among all atoms. Experiments show that the use of com- plete potentials leads to consistent performance improvements. Altogether, our work provides a theoretically principled and practically effective framework for crystal modeling.

---

Under review as a conference paper at ICLR 2023
REFERENCES
MI Baskes, Application of the embedded-atom method to covalent materials: a semiempirical potential for silicon. Physical review letters, 59(23):666,, 1987.
Harry Bateman. Tables of integral transforms [volumes I & II1, volume 1. McGraw-Hill Book Company, 1954.
Peter W Battaglia, Jessica B Hamrick, Victor Bapst, Alvaro Sanzaez- Vinicius Zambaldi, Mateusz Malinowski, Andrea Tacchetti, David Raposo, Adam Santoro, Ryan Faulkner, et al. Relational inductive biases, deep learning, and graph networks. arXiv preprint arXiv: 1806.00261, 2018.
Simon Batzner, Albert Musaelian, Lixin Sun, Mario Geiger, Jonathan P Mailoa, Mordechai Ko- rbbluhh, Nicola Molinari, Tess E Smidt, and Boris Kozinsky. E (3)-equivariant graph neural networks for data-efficient and accurate interatomic potentials. Nature communications, 13(1):
1-11, 2022.
Laurent Bétermin, Markus Faulhuber, and Stefan Steinerberger. A variational principle for gaussian lattice sums. arXiv preprint (rXiv:2110.06000 2021.
Max Born. Über elektrostatische. gittcrpotentiale. Zeitschrift für Physik, 7(1):124-140, 1921.
David Borwein, Jonathan M Borwein, and Keith F Taylor. Convergence of lattice sums and madelung's constant. Journal of mathematical physics, 26(11):999-3300) 1985.
JM Borwein and O- Yeat Chan. Uniform bounds for the incomplete complementary gamma function. Mathematical Inequalities and Applications, 12:115-121, 2009.
I David Brown. The chemical bond in inorganic chemistry: the bond valence model, volume 27. Oxford University Press, 2016.
Richard A Buckingham. The classical equation of state of gaseous helium, neon and argon. Pro- ceedings the Royal Society of London. Series A. Mathematical and Physical Sciences, 168 (933):264-283, 1938.
Chi Chen, Weike Ye, Yunxing Zuo, Chen Zheng, and Shyue Ping Ong. Graph networks universal machine learning framework for molecules and crystals. Chemistry of Materials, 31(9):3564- 3572, 2019.
Kamal Choudhary and Brian DeCost. Atomistic line graph neural network for improved materials property predictions. npj Computational Materials, 7(1):1-8, 2021.
Kamal Choudhary, Brian DeCost, and Francesca Tavazza. Machine learning with force-field- inspired descriptors for materials: Fast screening and mapping energy landscape. Physical review materials, 288)088801, 2018.
Kamal Choudhary, Kevin F Garrity, Andrew CE Reid, Brian DeCost, Adam Biacchi, Angela R Hight Walker, Zachary Trautt, Jason Hattrick Hattric--Simees, The joint automated repository for various integrated simulations (jarvis) for data-driven materials design. npj Computational Materials, 6(1) 6(1):1-13, 2020.
Gilad Kusne, Andrea Centrone, et al.
Richard E Crandall. Fast evaluation of epstein zeta functions. manuscript, at http:== www. perfsci. com, 1998.
Richard E Crandall and Joseph P Buhler. Elementary function expansions for madelung constants.
Journal of Physics A: Mathematical and General. 20(16):5497, 1987.
Murray S Daw and Michael I Baskes. Embedded-atom method: Derivation and application to impurities, surfaces, and other defects metals. Physical Review B, 29(12):6443, 1984.
Murray S Daw, Stephen Foiles, and Michael Baskes. The embedded atom method: a review of theory and applications. Materials Science Reports, 9(7-8):251-310, 1993.
10

---

Under review as a conference paper at ICLR 2023
Bernard Deconinck, Matthias Heil, Alexander Bobenko, Mark Van Hoeij, and Marcus Schmies. Computing riemann theta functions. Mathematicc of Computation, 73(247):1417-1442, 2004.
David Duvenaud, Dougal Maclaurin, Jorge eguilera-lparragen
Tim- othy Hirzel, Alán Sapuru-Guzise and Ryan Adams. Convolutional networks on graphs for learning molecular fingerprints. In Proceedings of the 28th International Conference on Neural Information Processing Systems, pp. 2224-2232, 2015.
Rafael
Paul P Ewald. Die berechnung optischer und elektrostatischer gitterpotentiale. Annalen der physik, 369(3):253-287, 1921.
Matthias Fey and Jan E. Lenssen. Fast graph representation learning with PyTorch Geometric. In ICLR Workshop Representation Learning Graphs and Manifolds, 2019.
Mark Galassi, Jim Davies, James Theiler, Brian Gough, Gerard Jungman, Patrick Alken, Michael Booth, Fabrice Rossi, and Rhys Ulerich. GNU scientific library. Network Theory Limited Go- dalming. 2002.
Johannes Gasteiger, Janek GroB, and Stephan Günnemann. Directional message passing for molec- ular graphs. In International Conference Learning Representations, 2019.
Johannes Gasteiger, Florian Becker, and Stephan Günneemann Gemnet: Universal directional graph neural networks for molecules. Advances in Neural Information Processing Systems, 34:6790- 6802, 2021.
Philippe Gaudreau, Richard M Slevinsky, and Hassan Safouhi. Computation of tail probabilities via extrapolation methods and connection with rational and padé approximants, SIAM Journal on Scientific Computing, 34(1):B65-B85, 2012.
Justin Gilmer, Samuel Schoenholz, Patrick F Riley, Oriol Vinyals, and George E Dahl. Neural message passing for quantum chemistry. In International conference on machine learning, pp- 1263-1272. PMLR, 2017.
Rhys EA Goodall and Alpha Lee. Predicting materials properties without crystal structure: Deep representation learning from stoichiometry. Nature communications, 1 11(1):1-9, 2020.
Frank E Harris. Incomplete bessel, generalized incomplete gamma, leaky aquifer functions. Journal of Computational and Applied Mathematics, 215(1):260-26) 2008.
Emiel Hoogeboom, Victor Garcia Satorras, Clément Vignac, and Max Welling, Equivariant dif- fusion for molecule generation in 3D. International Conference on Machine Learning, pp- 8867-8887. PMLR. 2022.
Weihua Hu, Matthias Fey, Hongyu Ren, Maho Nakata, Yuxiao Dong, and Jure Leskovec. OGB- LSC: A large-scale challenge for machine learning on graphs. In Thirty-fifth Conference on Neural Information Processing Systems Datasets and Benchmarks Track, 2021.
Jari Jalkanen and Martin H Müser. Systematic analysis and modification embedded-ane -tto poten- tials: case study of copper. Modelling and simulation in materials science and engineering, 23 (7):074001, 2015.
Dipendra Jha, Logan Ward, Arindam Paul, Wei-keng Liao, Alok Choudhary, Chris Wolverton, and Ankit Agrawal. Elemnet: Deep learning the chemistry of materials from only elemental compo- sition. Scientific reports, 8(1):1-13, 2018.
Dipendra Jha, Logan Ward, Zijiang Yang, Christopher Wolverton, Ian Foster, Wei-keng Liao, Alok Choudhary, and Ankit Agrawal. IRNet: A general purpose deep residual regression framework for materials discovery. In Proceedings of the 25th ACM SIGKDD International Conference on Knowledge Discovery & Data Mining, pp. 2385-2393, 2019.
Wengong Jin, Regina Barzilay, and Tommi Jaakkola. Junction tree variational autoencoder for molecular graph generation. International conference on machine learning, pp. 2323-222 2-322. PMLR, 2018.
11

---

Under review as a conference paper at ICLR 2023
Bowen Jing, Stephan Eismann, Patricia Suriana, Raphael John Lamarre Townshend, and Ron Dror. Learning from protein structure with geometric vector perceptrons. International Conference on Learning Representations, 2020.
Jones. Incomplete bessel functions. Proceedings of the Edinburgh Mathematical Society, 50 (1):173-183, 2007.
Diederik Kingma and Jimmy Ba. Adam: method for stochastic optimization. arXiv preprint 0rXiv:1412.6980 2014.
Klaus Kirsten. Generalized multidimensional epstein zeta functions. Journal of Mathematical Physics, 35(1):459-470, 1994.
Johannes Klicpera, Shankari Giri, Johannes T Margraf, and Stephan Günnemann. Fast and uncertainty-awar directional message passing for non-equilibriun molecules. arXiv preprint arXiv:2011.14115, 2020a.
Johannes Klicpera, Janek GroB, and Stephan Günnemann. Directional message passing for molec- ular graphs. arXiv preprint 2rXiv:2003.03123, 2020b.
Kenneth S Krane. Introductory nuclear physics. John Wiley & Sons, 1991.
Joseph PS Kung and Chung-Chun Yang. Complex analysis. Encyclopedia of Physical Science and Technology, 2003.
Eunkoo Lec, Kwang-Rye Byeng Lee, MI Baskes, and Bycong-Joo Lec. A modified embedded-atonn method interatomic potential for ionic systems: 2 nnmeam+ qeq. Physical Review B, 93(14): 144110, 2016.
Hark Lee and Wei Cai. Ewald summation for coulomb interactions in a periodic supercell. Lecture Notes, Stanford University, 3(1):1-12, 2009.
JE Lennard-Jones and Beryl M Dent. Cohesion crystal surface. Transactions of the Faraday Society, 24:92-108, 1928.
David Levin and Avram Sidi. Two new classes of nonlinear transformations for accelerating the convergence of infinite integrals and series. Applied Mathematics and Computation, 9(3):175- 215, 1981.
Meng Liu, Youzhi Luo, Kanji Uchino, Koji Maruhashi, and Shuiwang Ji. Generating 3D molecules for target protein binding. In Proceedings of The 39th International Conference on Machine Learning, pp. 13912-13924, 2022a.
Yi Liu, Limei Wang, Meng Liu, Yuchao Lin, Xuan Zhang, Bora Oztekin, and Shuiwang Ji. Spherical message passing for 3d molecular graphs. In International Conference Learning Representa- tions, 2022b.
Ilya Loshchilov and Frank Hutter. Decoupled weight decay regularization. arXiv preprint a1Xiv:1711.05101,22 2017.
Yves Louis, Yong Zhao, Alireza Nasiri, Xiran Wang, Yuqi Song, Fei Liu, and Jianjun Hu. Graph convolutional neural networks with global attention for improved materials property pre- 2020.
diction. Physical Chemistry Chemical Physics,
Youzhi Luo and Shuiwang Ji. An autoregressive flow model for 3D molecular geometry generation from scratch. In International Conference on Learning Representations, 2022.
Youzhi Luo, Keqiang Yan, and Shuiwang Ji. GraphDF: A discrete flow model for molecular graph generation. In Proceedings The 38th International Conference on Machine Learning, pp. 7192- 7203, 2021.
Philip M Morse. Diatomic molecules according to the wave mechanics. ii vibrational levels. Phys- ical review, 34(1):57, 1929.
12

---

Under review as a conference paper at ICLR 2023
Franziska Nestler, Michael Pippig, and Daniel Potts. Fast ewald summation based on nfft with mixed periodicity. Journal of Computational Physics, 285:280-313. 2015.
DM Riffe, Jake D Christensen, and RB Wilson. Vibrational dynamics within the embedded-atom- method formalism and the relationship to born-von-kármán force constants. Journal of Physics: Condensed Matter, 30(45):45570, 2018.
Jonathan Schmidt, Love Pettersson, Claudio Verdozzi, Silvana Botti, and Miguel A. L. Marques. Crystal graph attention networks for the prediction stable materials. Science Advances, 7(49): eabi7948, 2021. doi: 10.1 12b/cciadv.
Kristof Schütt, Pieter-Jan Kindermans, Huziel Enoc Sauceda Felix, Stefan Chmiela, Alexandre Tkatchenko, and Klaus-Robert Müller. Schnet: A continuous-filter convolutional neural network for modeling quantum interactions. Advances in neural information processing systems, 30, 2017.
Kristof Schütt, Oliver Unke, and Michael Gastegger. Equivariant message passing for the prediction of tensorial properties and molecular spectra. In International Conference on Machine Learning, pp. 9377-9388, PMLR, 2021.
Atle Selberg and Sarvadaman Chowla. On epstein's 's zeta-function. 1967.
Chence Shi, Minkai Xu, Zhaocheng Zhu, Weinan Zhang, Ming Zhang, and Jian Tang. GraphAF: a flow-based autoregressive model for molecular graph generation. In International Conference on Learning Representations, 2019.
John C Slater. The normal state of helium. Physical Review, 32(3):349, 1928.
Richard M Slevinsky and Hassan Safouhi. A recursive algorithm for the g transformation and accurate computation of incomplete bessel functions. Applied numerical mathematics, 60(12): 1411-1417, 2010.
Leslie N Smith and Nicholay Topin. Super-convergence: Very fast training of neural networks using large learning rates. In Artificial intelligence and machine learning for multi-domatt operations applications, volume 11006, pp. 369-386. SPIE, 2019.
Jonathan M Stokes, Kevin Yang, Kyle Swanson, Wengong Jin, Andres Cubillos-Ruiz, Nina M Donghia, Craig R MacNair, Shawn French, Lindsey A Carfrae, Zohar Bloom-Ackermann, et al. A deep learning approach to antibiotic discovery. Cell, 180(4):688-700 2020.
Audrey A Terras. Bessel series expansions of the epstein zeta function and the functional equation. Transactions of the American Mathematical Society, 183:477-486. 1973.
J Philipp Wagner and Peter R Schreiner. London dispersion in molecular chemssst-ryr--- steric effects. Angewandte Chemie International Edition, 54(42): 2015.
Anthony Yu-Tung Wang, Steven K Kauwe, Ryan J Murdock, and Taylor D Sparks. Composition- ally restricted attention-based network for materials property predictions. Npj Computational Materials, 7(1):1-10, 2021.
Limei Wang. Haoran Liu, Yi Liu, Jerry Kurtin, and Shuiwang Ji. Learning protein representations via complete 3D graph networks. arXiv preprint arXiv:2207.12600 2022a.
Limei Wang. Yi Liu, Yuchao Lin, Haoran Liu, and Shuiwang Ji. ComENet: Towards complete and efficient message passing for 3D molecular graphs. In The 36th Annual Conference on Neural Information Processing Systems, 2022b.
Zhengyang Wang, Meng Liu, Youzhi Luo, Zhao Xu, Yaochen Xie, Limei Wang. Lei Cai, Qi Qi, Zhuoning Yuan, Tianbao Yang, and Shuiwang Ji. Advanced graph and sequence neural networks for molecular property prediction and drug discovery. Bioinformatics, 38(9):2579-2555 25886, 2022c.
Anthony West. Basic Solid State Chemistry. Wiley, 1988.
Philip M Woodward, Probability and information theory, with applications to radar: international series of monographs electronics and instrumentation, volume 3. Elsevier, 2014.
13

---

Under review as a conference paper at ICLR 2023
Tian Xie and Jeffrey C Grossman. Crystal graph convolutional neural networks for accurate and interpretable prediction of material properties. Physical review letters, 120(14): 2018.
Hideki Yukawa. On the interaction of elementary particles. i. Proceedings of the Physico- Mathematical Society of Japan. 3rd Series, 17:48-57, 1935.
A GAUSSIAN LATTICE SUM
The Gaussian Lattice Sum (Bétermin et al. 2021) computes the summation of Gaussian functions centered the points given by shifted lattice, formally defined as
Go Gn(L, Σ
<Lk+o|"
(10)
KEE4
where [Rd the lattice matrix, v E Fd is vector inside a unit cell and c R+ is prefixed
constant. One the characteristics of the Gaussian Lattice Sum is that the term e ecc] rapidly decays as l becomes large, resulting fast convergence of Gq 2, c). As shown by Deconinck al. 22004), for R we have
d
d
)
Σ Vel Lkto
cf
r
(11)
2
x is the incomplete Gamma function and min(Ve]LA, |k zd. + 0}. We can obtain the upper bound of the Gaussian Lattice Sum by setting R such that
where 122,x) Ieee
Σ Eelkktuu
d
(
d
r
(12)
2
kE2d
B INCOMPLETE BESSEL FUNCTION
The incomplete Bessel function Harris 2008 is defined as
K2 Kz(x,y) v)
eut-le-ttu/ll
(13)
R. x 2 0. y 2 0, and v > O. In existing studies (Harris 2008 Slevinsky & Safouhi 2010 Jones 2007), K2 (z. y) can be analytically continued to R. In this work, we follow these studies and also use the analytic continuation for our calculations, i.e., we also consider IR.
where z,
B. FAST APPROXIMATION OF THE INCOMPLETE BESSEL FUNCTION
Computing the incomplete Bessel function extremely challenging as there does not exist any explicit closed-form solution. In this work, we investigate fast approximation of the incomplete Bessel function. Specifically, we adopt the algorithm in the work by [Slevinsky & Safouhi (2010), where the incomplete Bessel function K, (z y) is approximated by the G transformation with linear time complexity of (n)), and n is the number of iterations. As shown by Levin & Sidi 1981), the approximation (m(m) to f f(t) is given as the solution of
d'
((m) 1-" f(t)dt
Σ k=0
(())
3
0,
(14)
dx
r
where BL and (m) are unknowns, dan d (m)) 0, VI 0. σg min(sk, 1), and sy. is the
largest of E such that lim x ((k) (x) for l 0. m 1. Slevinsky & Safouhi
(210)) proved that the incomplete Bessel function satisfies Eqn.
14 with m I, through which
we can obtain
G(D) //" f(t)dt 20 f(x)
Σ 1==
(15)
14

---

Under review as a conference paper at ICLR 2023
To eliminate all the unknowns A Slevinsky & Safouhi (2010) applied the (x² dz operator times such that
(x² d dx
G
So n(I)d
ym
0. J f()ddt
(16)
By doing this, obtain
(x² d 10
(2.(z)
G
dx
r f(x))
(17)
(x.2 d )
(n(())
de
in which can have
d
N and Dn(x) = (x² d
Nx(x)) = (x² de
1n
(18)
dx
dx
with
Sod f(t)dt No(x) = x f(x)
and Do xo f(x)
(19)
This leads to recursive algorithm to approximate G To compute the incomplete Bessel function (x, v). Slevinsky Safouhi (2010) investigated the following property
Ku(r,y) + x² {" t--le - dt x² [o- t 1r--e e 5y dt,
(20)
in which the term Jodh 1-v-ee "p can approximated by G( Therefore, to approximate (x, y) we have
G = /- 1- 1 - ddt)
Nn((x) Dr(x)) So カ-ーーー -
(1(
(21)
(-) n- (x) (x² d drr-l (x (x- 1ov-lee---0
(x)
(z) (x)
K,(x,)) v) by recursively solving Nx and (x: (Gaudreau et al. 2012). The detailed expressions of Nº, (x) and (x) given in Slevinsky & Safouhi (2010). In addition, we follow Nestler et al. (2015) to further optimize the approximation of the incomplete Bessel function when the remaining part of the Taylor expansion of Ko (z, y) is small and we can approximate Ko(a, y) by the first terms of Taylor series such that
As result, we obtain the approximation (O)
and æ, are both small, e.g. x² +y² < 1. In this case,
(-1) n!
Ko( v)
="v"T(-n,z))
(22)
The detailed error bound by this expansion is shown in Nestler et al. 2015))
B.2 CONVERGENCE OF INCOMPLETE BESSEL FUNCTION SUMMATION
We define the summation of incomplete Bessel functions on a lattice as
Σ (o Lk + v + 77,
(23)
kE
where 1/, a. B,y E R are constants, a > 0, 3 0, 2 0, E Rd a vector inside a unit cell and L E ddd is the lattice matrix. We aim to prove the summation of incomplete Bessel functions is convergent and can be approximated with an error bounded by the Gaussian Lattice Sum introduced in Appendix A
15

---

Under review as a conference paper at ICLR 2023
Proof.
The incomplete Bessel function has upper bound such that
(2 v) ,
e mtt at- w/t di
- I 1- e mt de 1c-
(24)
where æ > and I is the incomplete Gamma function described in Appendix obtain
Based on this, we
[a/ka + of + 7,33)| < ['(-V, a Lk w/ + +2) (ax Lk + 112 + 7)
(25)
where 0: Lk + w/2 + > 0. As shown by Borwein & Chan such that
(2009).
[(
has upper bound
T(2, x) x²a
(1 + s)=-1ds
(26)
2>1
2
where 2 E R, x E R and x > 0. For a prefixed value R E R, R² > -V, R² > and R² > we have
(R)
Σ
Ko w/2 +7,8)
kEZz o kEZ Latoo Σ
KA(a Lk + u + + 8)1
kEZ t
(27)
Σ
caalLktol2
Ge
where Gn(L, a:) is the Gaussian Lattice Sum as described in Appendix A Therefore, the incom- plete Bessel function summation can be divided into two parts such that
Σ (a Lk + 2 +7,9)= kEZd
Σ
Ko Lk + w1 + 7,8) + c()).
(28)
kEZ L
where the first part is finite part inside an ellipsoid with a size of 1aP/a-//aa 7/aa /R/aa-/aa, part c(R) is bounded by the Gaussian Lattice Sum G a) which is convergent. Thus, the incomplete Bessel function summation convergent. Consequently, to approximate the incomplete Bessel function summation, we choose to evaluate the summation inside an ellipsoid with the size of 11//a-7/a for prefixed R R, such that R² > -v, R² > 1 and R² 7. Then the error <(Rt bounded by Gaussian Lattice Sum Gp L, v, v.a). We can further bound the error by the inequality [I] of Gaussian Lattice Sum as introduced in Appendix A such that
and the second
Go(L,
c(R) 1
Σ
e ak
- d dT
2 Q 2)9) -
(29)
2
KEZ .a LA+
where p = minmmin(aaa & E Zd K
+0}. 0}. This completes the proof.
C FAST ALGORITHM OF POTENTIAL SUMMATION
C.I INTEGRAL TRANSFORMATION
We denote G(L,v) v) as the potential summation and ULL, v) as the potential function with a lattice matrix L E Rddd and a vector v E Rd between two atoms inside a unit cell. For a potential summation S(a, b) of crystal with lattice matrix L, have S(a, b) G(L, vab). Based on these notations, we prove that the summation of the three introduced potentials can be transformed into integral form as
(
C(L = Σ
(30)
kEZd
where A, B, C, D constants derived from the corresponding specific potential forms.
16

---

Under review as a conference paper at ICLR 2023
Proof. We first prove that both the potential forms ULL 1/| + w and UL e e Lkttl
Lkttl can be written in the following integral form
UL, 0) D (-/² As Lk+ul dt.
(31)
1). For the potentials in the form of ULL v) 1/ Lk + =01² apply the Mellin transform such that
(L = / 10 pp-te-1laae -| dt
(32)
r(p)
Ln 01-p"
Thus, we obtain
ULL v) = LLk 2p [(p) J
t -t m| dt.
(33)
Apparently, we can obtain
1/T, 0, C = p and 1pT for the integral form of
U(L, v) = 1/ 1/ Lk + u 20
2). For the potentials in the form ( = e ol we consider the inverse Laplace trans-
form on e -ays as shown by Bateman 1954). such that
L (-Iee-vv
(34)
2
Therefore, apply the Laplace transform in Eqn.
34 to derive the integral form of e Lk
ULL v) e [b
a
1 1
a2
2
(35)
-ē L x²
2π J
(t + Tt).
Apparently, we can obtain A = 1, B U(L v) a|(boo
o²
for the integral form of C = - and D = a 2
4
Finally, we conduct a summation for these two types potentials U (L, v) in the space 2d as
G(L = Σ U(L
(36)
REZ4
This completes the proof.
In fact, the summation (LL v) = Lk is a special case of multidimensional zeta func-
tion (Crandall & Buhler 1987 Terras 1973 Crandall] 1998), which is generalization of Riemann
zeta function. The multidimensional zeta function (Crandall & Buhler 1987) is defined as
2niu Lk
(s
(37)
Lk 0/87
kezg
where E C, L E Rdxd and 0 E Rdd is not hard to show summation of U(L, v) Lk + v) 22 can be expressed as (2p: -v), and therefore functional properties of the mul- tidimensional zeta function also hold for summation of U(L v) = 1/| + v/²". For instance, Zs u. v) has an analytic continuation the entire complex plane, except for simple poles at s = and = (Crandall & Buhler 1987. As a result, summation U(L,v) = 1/| + v 2p are all convergent for p C except 0 and d/2 by analytic continuation. Moreover, the multidimensional zeta function can be written the form of an integral summation by Eqn. [3 as
2 Lk
((s/2) /
Zz (s u,
Σ kazza
1 20iu Lk ( dt.
(38)
Lk - w/²
17

---

Under review as a conference paper at ICLR 2023
Based on this, also split the integral and apply Poisson summation Equ. obtain two summations of incomplete Bessel functions to evaluate this series. For more details on the multi- dimensional zeta function, refer readers to Crandall & Buhler (9877: Teraas 1973): Crandall
1998
Kirsten 1994); Selberg & Chowla 1967
C.2
CALCULATING INTEGRAL SUMMATION
As shown in Sec, 3.4 G v) can be written the summations in the Euclidean space directly and then Fourier space:
G(
= D Σ
Ar-le-AA/ ddt
REZS
D Σ
1 A P d
Σ ( C dt
(39)
KEZ -0
REZd
GFourier (L, L, + GdircaLL,
dt denotes the summation Fourier denotes the summation in direct space. Apparently, Gaireaa,, (L, v) is already the form of the incomplete Bessel function summation. We apply the analytic continuation to Gdirec [L, v) to expand the domain of C in Gaireca [L, v) as detailed in Appendix C.3 we have
where Goourie L.v) = DE So tc-
space, and Gairecl L) D E S.
Gäiree v) =
Σ kEzz
(40)
for all constant C R. Below, we prove that Grearier L, v) be deduced into the incomplete Bessel function summation.
Proof. Inspired by the Ewald summation (Ewald 1921 Crandall 1998
we consider Grourier (L, v)
on the reciprocal lattice using the Poisson summation (Crandall 1998):
2niw Lk-u-LL-xt xt -- 2riiww
Σ
2mi k.u F =
(41)
det L
KEE
kEZd
where E Rdd is vector and 0 in our case, and L' = L(L"L) is the lattice matrix for the reciprocal lattice. As a result, We obtain
Grouter () D Σ
-- - Anlk-o
KETA
C-
Lkto
Anl di(t + A
1 A
AC
kEZ
J0
1 D
c- 4- n
n- =/ AB dt (Eqn. 41
det
L AC
KEZA J0
(42)
1
Σ
c- 1 2 k. L'kk' dt(t + At)
det L A
KEZd
J0
1 D
t² -1
2 1- -Bt dt 1
det A
1
D
2 k.u K Ke_
2 B,O).
det
A
REZZ
We also apply the analytic continuation to Grourier (L, v) in the last step. Apparently, Grounier (L, v) deduced into the incomplete Bessel function summation, and this completes the proof. 1
18

---

Under review as a conference paper at ICLR 2023
Therefore, both Gruurier v) and Gairea (L, v) can be expressed by the incomplete Bessel function summation forms:
G(L, 0) Gdireal L, 1) + Grourier L, v))
=D
Σ K- K-c(Ar Lk + v
kE?d
(43)
D
2 k. K Kc- T |
+ B,0).
det L A
kEZ
As shown in Appendix B.2 the incomplete Bessel function summation ΣnEZA (a LLk w/2 7/ 7) is convergent and can be approximated. Therefore, G(L, v) can also be approximated.
C.3 ANALYTIC CONTINUATION OF POTENTIAL. SUMMATIONS
To represent the series that is not convergent, including the inverse summation ΣKEZ [Lk+ol as shown in Sec. 3.4 we need to investigate the analytic continuation of potential summation. As shown in Appendix C.2 the potential summation can be written
G(L = Σ
cC-le Aon-AALLLffffee Edtt
kEZ4
(44)
E ( 2 D Σ C-- (- di.
kEZ
kEZd
These two parts can be deduced into incomplete Bessel function summations by analytic con- tinuation. Analytic continuation is technique to extend the domain P of a given analytic function f(x). Consequently, we denote far as an analytic continuation of f(x) to Q, where we denote the domain Q containing P. and denote function f(z) that is analytic Q, and f(x) fx) holds for all x in P. As shown by |Kung & Yang (2003), the analytic continuation is unique and satisfies the permanence of functional relationships, i.e., the equations holding for f(x) will also
hold for f(x).. In our case, we can expand the domain of constant C' in G(L, v) to C' R such that R. For example, for summation of the potentials in the form of EkEZe L/|Lk + 0/2p, we have C = p as shown Appendix C.1 Analytic continuation enables us to compute the summation when
G(L, v) is well-defined for any C E
0.5, which is initially divergent. Formally, we assume that C originally defined in the domain P C IR. As derived in Appendix
C.1 for C' E P, we have
G(L,v) Σ
-
D
Σ kEZ
C--e Am|Lk+ $ dt
KEZS
=D
Σ
K c-c(AA + o)2
+
(45)
KE24
D
L'ke² A
Σ
2 Ko_
+ B,0).
det L A
kE24
Furthermore, have -C and C 2 in the incomplete Bessel function, which is analyt-
ically continued to E R (Jones 2007). That is, P is contained by v's domain R. Therefore, the incomplete Bessel function summations the analytic continuation of the corresponding potential summations.
Our approach ultimately to capture the total contribution of the potential in the crystal system. Based on this, to explain the advantage of analytic continuation, we consider the total infinite poten- tial summation inside unit cell such that
S= >>> Σ
S(a, b)
Σ Σ V(a, b).
(46)
aEAM BEA
a bea.bEĀ
where A and À are described in Sec. mate convergent summation S(a, b), we can directly approximate by
2.1 If we use analytically continued function approxi-
S
ΣΣ
S(a, b).
(47)
READEA
19

---

Under review as a conference paper at ICLR 2023
As for non-convergent summation S(a, b), due to the permanence of functional relationships of analytic continuation, we still obtain
S= Σ Σ S(a. b).
(48)
aEA BEA
This implics that we can use analytically continued summations to approximate the total contribution of potentials where individual potential summations are initially divergent. In addition, for a non- convergent summation S(a, b), analytic continuation will result in unusual value such negative energy instead of infinity. This is useful in practice since the infinity will cause a numerical explosion in training. A famous example of using analytic continuation for crystal energy prediction is the Madelung constant of NaCl, which is derived by the summation of Coulomb potentials among Na and Cl ions and calculated by analytic continuation of an absolutely convergent series (Borwein et al. 1985). We further show energy calculation of NaCl by analytic continuation in Appendix
C.4
IMPLEMENTATION AND NUMERICAL EXAMPLES OF APPROXIMATION
Here we describe the implementation of our algorithm as in Eqn. (43) Considering the fact that solving the inverse of the incomplete Gamma function is complicated, instead, provide proper value R and then calculate its corresponding error bound e based on Eqn. (29). Given a lattice
matrix L Raxd vector E Rd inside unit cell, and the constants A, B, C, D derived from specific potential functions as described in Appendix C.T we aim to evaluate these two parts
D. ExEza BriikkuKK Kc_g( D Exeza K K-o Lk + 0/2, B) and det A7
L'k A
B,0).
To evaluate Gdirec [, ExEza K- K-c(AT] Lk + = , ,B, we derive the following steps.
Step 1: Determine the value R such that 12 C. R2 1, and calculate the error bound
2 p ((2))T( (R 2)),, where p min (v Ar Lk] 1k gd 03.
Step 2: Select points inside an ellipsoid such that P (kl Lk wl S R/ Ar).
Step 3: Evaluate the incomplete Bessel function summation Gditeal (L, v) by calculating the every term DK Lk + 212 BB) based on Appendix B.1
To evaluate Groutier (L, v) = det L A ΣKEZA e 2 k Ko- lowing steps.
L'k|² A + B,0), we derive the fol-
Step I: Determine the value R such that 12 > C, 12 > 1, R2 > B, and calculate the error
bound ((3))T( (VR2 - 2) where p min / L'k
k 2d k + 0}.
Step 2: Select points inside an ellipsoid such that P (k) Lk wl < / A(R -B)T).
Step 3: Evaluate the incomplete Bessel function summation G Fourier L, v) by calculating every term
D 2 ko Ferl'kkKKC- Kc_ |''|" + B) based Appendix B.1
det L A
Our implementation is based on Cython and GNU Scientific Library (Galassi et al. 2002). ). in which the native incomplete Gamma function and Bessel function are used to implement the incomplete Bessel function. We conduct numerical experiments on Intel Xeon Gold 6258R CPU. We show the evaluation examples in Table the corresponding error bound and evaluation time. The running time is at the scale of milliseconds.
C.5 POTENTIAL SUMMATION EXTENSIONS
To highlight the generality our potential summation method, in this section, we introduce addi- tional potentials that can be converted to our general integral form in Equ. 31 including Lennard- Jones potential, Morse potential, and screened Coulomb potential. These potentials are only used for some specific types of materials, such as gas and fluid materials.
Lennard-Jones Potential (Lennard- Jones & Dent 1928) is an intermolecular pair potential that is usually used for gas organic materials. The commonly used expression for the Lennard-Jones
20

---

Under review as a conference paper at ICLR 2023
Table 6: Numerical examples of our algorithm. Here, ((x) E 1/n". can be found in Appendix C.4
Implementation details
Ground Truth
Evaluation
Estimated Error
Real Error
Time
2* ((2) 3.28986813
3.28068288
8e-1
9e-3
0.002
2* ((2) 3.28986813
3,28984070
6e-2
3e-5
0.003
2+( 3.28986813
3.28986812
7e-4
le-8
0.003
(+{(2) 3.28986813
3,28986813
le-6
le-8
0.003
2* = 2.40411381
2,404
le-6
< 1e-8
0.003
(+{(4) 2.16464647
2.16464647
le-6
1e-8
0.003
E 2 =0 | 6.02681204
5.99068949
3
4e-2
0.003
mEZ? n20 In] 6.02681204
6.02670959
2e-1
le44
0.003
n In 6.02681204
6.02681199
3e-3
5e-8
0.003
m n=0 |n| 6.02681204
6,02681204
5e-6
1e-8
0.003
mEZ e
2.16395341
2.16395326
4e-1
2c-7
0.002
L mEZ In = 2.16395341
2.16395341
3e-4
< le-8
0.004
ME23 E -mmi 25.39268269
25.39268214
2.5
5e-7
0,003
MEZZ E |nj 25.39268269
25.39268269
le-2
< le-8
0.003
potential is
(L, = 4e (
12
0
(49)
Lk + vl
Lk + 2
where e and a are hyperparameters. And the summation of UrJ(L, v) can converted to two potential summations of type 1/| + 2 with p = 3 and 6, such that
Lk + 2/12 t/+2 )
v) = UL,(L) = 4c 12
Σ kEz
1
0
Σ KEE
Lk
2/
KEZz
(50)
where we show the calculation of the potential summation of type 1/| u in Appendix C.1
Morse Potential (Morse 1929) an interatomic potential of diatomic molecule and can be used for simple molecular materials. The Morse potential has mathematical form of
U Morse (L, v) D. (e ( (2l[((Lkto)-rr) 2e - -a((Lk+to--r))
(51)
where De and re are hyperparameters. Similarly, the summation UMoore e [L, v) can be converted to two potential summations of type e a-| with a = a and a = 2a, such that
GMoree v)
Σ KEEZ
UMorse L. De
2are
Lk+ol
are
Σ kEZ4
-a
Lk+tvl
kEzZ
(52)
where we show the calculation of the potential summation of type e e- Lb in Appendix C.1
Screened Coulomb Potential or Yukawa potential (Yukawa 1935) represents the Coulomb in- teractions with damping of electric fields. an important potential reflecting the behaviors of charge-carrying fluids or particles in semiconductors. The screened Coulomb potential has an ana- lytic form of V(a, b) d(a.b) cxp(-ad(a, b)) where d(a, b) is the distance between atom a and b, Za. 2b are charges atom a and b, e elementary charge constant and a a scaling hyperpa- rameter. Since 2a, 2b. constants and can be extracted outside the summation, we can obtain the simplified screened Coulomb potential
=ge2
U Usereened L, v) = e a
(53)
Lk
21

---

Under review as a conference paper at ICLR 2023
Zi
e, s s
s
X1
Xj
iip Sij siy. st
92
256
Linear
Embedding
llo
256
256
768
Interaction
MLP
MLP
Interaction Block
Interaction
256
Interaction
256
AvgPooling
256
BN
236
255
Linear
ReLU
Output
X;
+: Plus
d:: Sigmoid
a Concatenate along Gra. dimension
O: Hadamard Product
Σj Aggregate over j
Figure 2: The developed network architecture for PotNet.
Consider the inverse Laplace transform on e-ava/V5 3, by using Bromwich contour with branch points, obtain
L- (e [33
(54)
Therefore, we can apply the Laplace transform in Eqn. 53 such that
eaalLk+ol
/ t - 1o-1kkooltttt of dt.
U scrreened (L, v)
(55)
Lk +
n
Then we can obtain A = 1/T, B B= a², and D 1/ in Eqn. 30 to fit screened Coulomb potential into our potential summation method.
D MODEL IMPLEMENTATION
D.1 POTNET IMPLEMENTATION
The employed network architecture shown in Fig. 2 Since our major contribution is to consider infinite interatomic potentials, we simply design our network architecture following the commonly used settings. Specifically, existing methods for 3D graphs (Xie Grossman, 2018 Schütt et al., 2017 Klicpera et al.] 2020b a Liu et al.) 2022b Gasteiger et al.) [221]] Schütt et al. (021) share a similar architecture, which usually contains an input block, an interaction block, and an output block. Without loss of generality, take the updating process for node example illustrate
the network.
The Inputs contain atomic features and potentials. z; the 92-dimensional atomic feature for any atom following CGCNN (Xie & Grossman 2018). Below we denote d as the interatomic distances and show our potential features. For our implementation of infinite potential summations S(a, b) in Eqn. (7). we add up the infinite summations of three potentials, including the Coulomb potentials,
London dispersion potentials, and Pauli repulsion potentials described in Sec. 2.2 We simplify
the mathematical form of Coulomb potentials V Coulomb (a, b) 4xredd hj to V VCoulomb b)
-e/d((a, b) where is hyperparameter, because e, म, to all known constants, and Ea; 2b can
22

---

Under review as a conference paper at ICLR 2023
be learned from atomic features. As explained in Sec. 4.1 eº -c/ /d denotes the Coulomb
potentials for the local crystal graph. s = - Eg c = e/d/d s" = =E e -ad denote the summations of Coulomb potentials, London dispersion potentials, and Pauli potentials for the infinite crystal graph, and we set < = 4.0 for ef, 6 = - 1.0 for s", 6 4.0 for s' and a = 3.0 for s", respectively. For simplicity, we apply them to the RBF embeddings with the same cutoff. Since ef and s' are negative values and have values exceeding the RBF cutoff, we use an exponential function to make their absolute values smaller.
The input block contains Linear layer and an Embedding layer. For each node i, the Linear layer employed to generate a 256-dimensional vector as the input node features to the first interac- tion layer. For each edge, the Embedding layer is employed to map the Coulomb potentials and summations of Coulomb potentials, London dispersion potentials, and Pauli repulsion potentials to 256-dimensional embeddings by using 256 RBF kernels with centers from -4.0 to 4.0.
The Interaction block contains several Interaction layers. Each layer updates the feature vector of node i based on features of the neighboring nodes and potential embeddings of the connected edges. Particularly, for any of the neighboring node j of node i, the corresponding potential embeddings e. Siy: s. and fij for edge ij are all produced by the Embedding layer.
The Readout block contains an ArgPooling layer and another Linear layer. We first use the Arg- Pooling layer to aggregate features from all nodes in graph and then use the Linear layer that maps the hidden dimension of 256 to the final output which is a scalar.
D.2 CONFIGURATIONS OF RETRAINED MODELS
In this section, we show detailed configurations of retrained models of CGCNN (iie & Grossman
2018). SchNet (Schütt et al.. 2017), MEGNET (Chen al. 2019), GATGNN Louis et al.j (0220)
and ALIGNN (Choudhary& & DeCost (02I). If not specified, models are trained with a radius cutoff of 8.0 using the Adam (Kingma & Bal (014)) optimizer with weight decay (Loshchilov & Hutter
2017] and one cycle learning rate scheduler (Smith & Topin 2019).
CGCNN (Xie & Grossman 2018). We directly use the publicly available code from Xie & Gross- man (2018) to build and train the CGCNN model. We build the model with 128 hidden dimensions and three message-passing layers and train the model for 1000 epochs with a batch size of 256 and an initial learning rate of le-2.
SchNet (Schitt al. 2017). We directly adopt the SchNet model from PyTorch Geometric Fey & Lenssen 2019) with 128 hidden dimensions and six message-passinng passing layers following the original paper. We train SchNet with an initial learning rate of 5e-4 and batch size of 64 for 500 epochs.
GATGNN (Louis et al. 2020), We directly use the publicly available code from Louis et al. (2020) to build and train the GATGNN model. We build the model with 128 hidden dimensions and three bassage-passing layers following the original default settings. We train the model with an initial learning rate of 5e-3 and batch size of 64 for 500 epochs.
MEGNET (Chen et al. 2019). We directly adopt the MEGNET model from the publicly available code from et al. (2019), Following the original paper, we use three message-passing layers with the same feature dimensions as mentioned in the original paper and use Sct2Set readout function. We train MEGNET with an initial learning rate of 1e-3 and batch size of 128 for 1000 epochs following the configuration settings mentioned in the original paper.
ALIGNN (Choudhary & DeCost 202T.. We directly use the publicly available code from Choud- hary & DeCost (2021). We use the official best model configurations of ALIGNN to train ALIGNN models with an intiill learning rate of Ic-3 and batch size of 64 for 500 epochs.
E LINEAR ENERGY MODELING USING INFINITE POTENTIAL SUMMATION
In this section, we provide examples of calculations where the true energies of these materials can be directly approximated by linear combinations of our infinite potential sums. These are special cases of Eqn. 6 with a linear embedded function G. Specifically, we evaluate the total energy per atom of NaCl and two other materials (MgO, LiF) whose crystal structures are similar to NaCl.
23

---

Under review as a conference paper at ICLR 2023
Since they pure ionic crystals and Coulomb interactions dominate the system, we first con- sider the total electrostatic energy and Coulomb potential summations. For a neutral system, the electrostatic energy is convergent because the Coulomb potentials cancel each other. To cal- culate the energy of crystal, we cannot derive it directly due to its complexity. Instead, we divide it into many individual infinite potential summations. Although these individual sum- mations might be originally divergent, they will all be well-defined under analytic continuation and the total energy will also be convergent.
Nat
cl-
Figure 3: Crystal structure of NaCl.
As shown in Appendix C.3 analytic continuation allows compute potential summations that
are initially divergent, and also allows us to compute the convergent summations by analytically continued potential summations. For example, we can calculate ΣnEZ, [-1)" by analytic con-
tinuation such that
Σ n
(-1) n
Σ n
20 - n2 n
20
-(11)))((11))) -
)
(56)
where ((s, a) is Hurwitz zeta function and ((s.a) = (ntaja when 1, + 0.11, and its analytic continuation elsewhere. That is, we can use the analytically continued zeta function C(s, a) to precisely approximate a convergent series.
Inspired by analytic continuation, we can approximate the total energy per atom of crystals of NaCI by analytically continued infinite potential summations. Due to the symmetry of the NaCl cell, we only involve atoms a, c, d in our calculation. Specifically, as shown Fig. atom represents the body center Nat, atom b represents the face center CI atom represents the edge center Na and atom d represents the corner Cl- Based on this, assuming that the side length of the unit cell 1, the total energy of Nat is approximated by the total Coulomb interactions with atom such that
ENa = -N
1 Naa
Σ
1
+
6
2Nab
Σ
4TeD
d(a, u)
4TED
d(a, u
Uc
u
ueAr
+ - 12.
NNa2Nel
+
-8-
Σ
4TEO
d(a, u)
4TEE
dia u)
ueA
ueAd
--NA
1
+
2C12
Σ
4TEOrO
d(a,u)
4TeOrO
d(a,u)
uEA
ueAD
ZNabNael
Σ
ZNre
Σ
+3
d(a u)]
4TEOTO
Ma.u)
4TTEro
ME
ueAs
(57)
NalzNellzaule
-3- +3 $(a,c) $(a,
4mcoro
loceeellzcule2
141864874 - 3 (-0.04796615) +3 (-0.29126077)
4T
(-0.40096799)
NalzNol|zcule
(-1.7475646)
4
8.81eV.
where NA is Avogadro constant, d(a, u)) is the distance between atom a and u,, is the minimum distance between Na and Cl, d d= d/ro is the normalized distance, Aa denotes the set of atoms containing atom and all its repetitions, $ is the infinite potential summation, approximated by
24

---

Under review as a conference paper at ICLR 2023
infinite potential summation method in Sec. 3.4 Na 2Ct charges of Nat and Cl is the elementary charge constant, 6g is the permittivity constant of free space, and the coefficients 1, 6,, $12, $8 denote the fraction of atoms in a unit cell. We finally obtain constant 1.744766 from our infinite potential summations. In fact, this constant is exactly the famous Madelung con- stant M (Borwein al. 1985). Also, by considering additional repulsion term, we can derive the calculation result in Equ. 57 to the famous Born- equation (Bom 1921
Nal= 1| e
E
(58)
where 2 =* are the charges cation and anion, M is the Madelung constant computed from infinite Coulomb potential summations, and n. is the Born exponent mea- suring the effect of repulsion. By choosing n 9, we can derive an approximation for the total energy of NaCl of 7.84 eV.. Similarly, we can apply Eqn. 58 to MgQ and LiF. We further show these approximations in Table. and it can be noticed that these approxima-
Table 7: Total energy per atom approxi-
mation of NaCI, M0 and LiF.
Formala
To
"
Ground
Egn
NoCI
282 pm
9
8.15
c
MgO
210 pm
6
39.33 eV
29.45
LAF
201 pm
7
10.67 eV
10.60 cV
tions already give rough results compared the ground truth energy. This implies that our features can serve as good starting point for machine learning models to learn the ground truth energy. Apparently, previous methods cannot achieve this due the lack of such informative features. worth noting that the Madelung constant is typically unknown because those coefficients for the infinite potential summations depend on the charge distribution in the system, which we do not know the beginning- Also, these crystals special cases of Eqn. 6 with linear embedded function G, while G is typically a nonlinear function (Daw & Baskes 1984. Therefore, the network serves the purpose of learning those coefficients learn the Madelung constant and providing nonlinearity.
25


Reviews:
[{"review_id": "x2HG5GEWIoR", "reviewer": "Reviewer_q5TT", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "The paper proposes a graph neural network (GNN) for property prediction for crystals. A key issue in modeling crystals is to handle long range interactions. This is handled by using interatomic potentials as edge features in the GNN and approximating the infinite summation over all pairwise interactions using Ewald summation.\n", "summary_of_the_review": "Good paper with a simple but interesting and novel idea. Reasonable empirical evaluation and comparison. The point raised by FZse should be addressed.", "strength_and_weaknesses": "\nStrengths:\nGreat example of how incorporation of physical knowledge into a GNN structure leads to a more efficient architecture.\n\nThe methods is demonstrated on several examples, and ablation reveals the importance of the infinite summation to handle long range interactions.\n\nWeaknesses:\n\nThe description of Ewald summation in the main paper is very concise - it would be nice with a bit more intuition about the methods for readers that are not familiar with this. Also, it is not clear to what extent this is just application of a known result/method, or if this is a novel mathematical result. In case there is something substantially new in the methods, this should be in the main paper in my opinion. I did not carefully check the derivations in the appendix.\n\nI expected to see comparison with other GNN approaches where the cutoff was varied, to be able to more directly understand the importance of the infinite summation. It is not clear from the paper (but from the references) what the cutoff was set to in e.g. the SchNet comparison (5Å i think)\n\nEDIT: (after reading FZse's review) Since the electrostatic (Coulomb) potential is summed infinitely for each (charged) atom, does that not mean that these individual sums in general will diverge, even though the system is electrically neutral? \n", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "XWwsIsBawN", "reviewer": "Reviewer_FEhN", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "This paper proposes a new graph deep learning method, PotNet, with several innovations in the field of crystal material modeling. The proposed method models interatomic potentials directly as edge features instead of only using distances and models the complete set of potentials among all atoms with approximations of infinite potential summations instead of only between nearby atoms. With the computations of complete interatomic potentials, the representation learning is performed through message passing neural networks. Evaluation experiments are also performed on the JARVIS and Materials Project benchmarks, showing the proposed method leads to consistent performance improvements with reasonable computational costs.", "summary_of_the_review": "This paper proposed a novel graph deep learning method, PotNet, which can improve the performance of the graph neural networks potentials in the field of crystal material modeling with certain computational costs. However, compared with the existing method, the gain of the proposed method in accuracy is not remarkable with the time costs of about 4 times. It is better to illuminate the specific application situations where the proposed method has more pronounced effect or further improve the algorithms to achieve faster implementation.", "strength_and_weaknesses": "Strength:\n\nThe proposed method improves the accuracy of the GNN-based interatomic potentials for crystal property prediction.\n\nWeaknesses:\n\n1. Compared with ALIGNN method, PotNet method consumes about four times as much time when infinite potential summations are considered, but the gain in accuracy is not that obvious.\n\n2. The comparisons in efficiency between the proposed method and direct DFT calculations was not made. Similarly, the advantages and disadvantages of the proposed method compared with the machine learning potentials using the constructed vectors as crystal structure input.", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "LsGQR6oJbv", "reviewer": "Reviewer_FZse", "summary": "", "questions": "", "limitations": "", "rating": 3, "confidence": 5, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "This paper proposed to incorporate Coulomb potential, London dispersion potential and Pauli repulsion potential in the input features of deep learning models, in order to deal with the long range interactions in PBC systems like crystal materials. It is proved that the proposed model outperforms SOTA in two material datasets: The Materials Project and JARVIS. ", "summary_of_the_review": "This motivation for using closed-form formula to describe long-range interactions in neural network is a great idea. However, a fatal mathematical error (see above) make the computation of Coulomb term completely wrong, which diminish the contribution of this paper.", "strength_and_weaknesses": "Strength: the motivation and idea of the paper is awesome. Traditional neural networks for molecular property prediction usually use local information as atom features, which limits the description for long-range interactions. However, unlike short range interactions which involves complicated electronic structure, the formula of long-range interactions are usually easy. So incorporate fix-formula is a promising solution for long-range interactions.\n\nWeakness: The paper has a fatal mathematical error, which invalidate the computation of Coulomb part. The error occurs in Appendix C.1. In formula (31), it is ... e^{-|Lk+v|t}, but in formula (32), it is miswrote as ... e^{-t|Lk+v|^2}. In appendix C.2, C.3 and C.4, it is e^{-t|Lk+v|^2} everywhere. This error is fatal because if we change e^{-t|Lk+v|^2} back to the correct e^{-t|Lk+v|}, the Poisson summation (39) does not hold, so the summation of Fourier term (40) also does not hold. In fact, if the authors did not make the mistake, the analytic continuation of infinite summation of Coulomb potentials should be diverged, similarly to the fact that Riemann zeta function diverges at s=1.", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "1: The main claims of the paper are incorrect or not at all supported by theory or empirical results.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "Not applicable"}, "extra_text": {}}]

Author rebuttal:

Decision: Reject