Mila Ai4science

Quebec Mila

Alex Hernandez-
email: alex.hernandez-garcia@mila.quebec

Garcia Mila

Alexandre Duval
email: alexandre.duval@mila.quebec

Alexandra Volokhova
email: alexandra.volokhova@mila.quebec

Yoshua Bengio
email: yoshua.bengio@mila.quebec

Divya Sharma
email: dsharm23@jh.edu

Pierre Luc Carrier
email: pierre.luc.carrier@mila.quebec

Michał Koziarski
email: michal.koziarski@mila.quebec

Victor Schmidt
email: schmidtv@mila.quebec

Crystal-GFN: sampling crystals with desirable properties and constraints

Accelerating material discovery holds the potential to greatly help mitigate the climate crisis. Discovering new solid-state materials such as electrocatalysts, superionic conductors or photovoltaic materials can have a crucial impact, for instance, in improving the efficiency of renewable energy production and storage. In this paper, we introduce Crystal-GFN, a generative model of crystal structures that sequentially samples structural properties of crystalline materials, namely the space group, composition and lattice parameters. This domain-inspired approach enables the flexible incorporation of physical and structural constraints, as well as the use of any available predictive model of a desired physico-chemical property as an objective function. To design stable materials, one must target the candidates with the lowest formation energy, which is used as an objective to evaluate the capabilities of Crystal-GFN. The formation energy of a crystal structure is predicted here by a new proxy model trained on MatBench. The results demonstrate that Crystal-GFN is able to sample diverse crystals with low formation energy.

Introduction

Materials discovery plays a vital role in transforming numerous industries that are currently responsible for a significant fraction of the global greenhouse gas emissions. From developing new-generation photovoltaic panels or solid-state batteries to tackling carbon capture or catalysis, the quest for innovative materials with targeted properties has the power to reshape the technological landscape.

However, the discovery of new materials is an extremely complicated task at various levels of innovation: crystal structure generation, property prediction and developing successful synthesis protocols. In practice, before scientists can think of designing a material with a desirable property, they must search through the vast space of possible materials and use their domain expertise to investigate a few/ handful of promising structures. This traditional trial-and-error process is extremely costly both in terms of time and efforts, providing machine learning (ML) with an opportunity to greatly accelerate the generation and evaluation of promising candidates.

While ML has already made significant progress with respect to property prediction and generation of small molecules [START_REF] Long | Constrained crystals deep convolutional generative adversarial network for the inverse design of crystal structures[END_REF][START_REF] Ren | An invertible crystallographic representation for general inverse design of inorganic crystals with targeted properties[END_REF][START_REF] Pakornchote | Diffusion probabilistic models enhance variational autoencoder for crystal structure generative modeling[END_REF][START_REF] Zheng | Towards predicting equilibrium distributions for molecular systems with deep learning[END_REF] 2 , its impact in the field of solid-state materials lags behind. At the root of this discrepancy lies the particular structure of materials, whose periodicity in 3D space imposes periodic boundary conditions that make it difficult to use finite graphs for crystal representation. As a result, the sole task of generating a stable candidate is complex. It involves exploring the combinatorial space of chemical compositions, considering all possible periodic arrangements of atoms as well as respecting complex physics that make them viable, stable and metastable structures.

In this work, we introduce Crystal-GFN (Fig. 1), a generative model designed to sequentially sample crystal structures in a space inspired by theoretical crystallography. Instead of directly generating atomic positions, Crystal-GFN samples the space group, the composition and the lattice parameters of a structure (Section 3). This domain-inspired approach facilitates the incorporation of chemical and structural constraints, as well as the possibility to extend the sampling space, for example to obtain atomic positions. The GFlowNet framework offers the ability to sample proportionally to a reward function, which enables diverse sampling. Through an empirical evaluation (Section 4) of a Crystal-GFN trained with a predictive model of the formation energy, we demonstrate that our proposed method is able to sample diverse candidates with low formation energy.

Background

In this section, we briefly review the necessary background on GFlowNets and crystallography.

GFlowNets

Generative Flow Networks [START_REF] Bengio | Flow network based generative models for non-iterative diverse candidate generation[END_REF] (GFlowNets or GFN for short), were introduced as an amortised inference method to sample from high-dimensional distributions, where both traditional methods such as MC-MC and reinforcement learning are inefficient in terms of mode mixing. In essence, 2 We refer the reader to Appendix A for a more detailed related work overview.

Crystal-Lattice System triclinic, monoclinic, ... [START_REF] Gebauer | Symmetry-adapted generation of 3d point sets for the targeted discovery of molecules[END_REF] Point Symmetry centrosymmetric, polar, ... [START_REF] Callum | 3-d inorganic crystal structure generation and property prediction via representation learning[END_REF] Space Group GFlowNets are designed to learn a sampling policy π to generate objects x ∈ X proportionally to a non-negative reward function R(x), that is π(x) ∝ R(x). This facilitates the discovery of multiple modes of a given function, which is a desirable objective in scientific discovery [START_REF] Jain | GFlowNets for AI-driven scientific discovery[END_REF]. For instance, GFlowNets have been successfully applied for biological sequence design and molecular property prediction as part of active learning loops [START_REF] Jain | Biological sequence design with GFlowNets[END_REF][START_REF] Hernandez-Garcia | Multi-fidelity active learning with gflownets[END_REF].

A key property of GFlowNets is that objects are generated sequentially. Starting from a special state s 0 , transitions s t →s t+1 ∈ A are applied between states s ∈ S, forming trajectories τ = (s 0 → s 1 → . . . → x), where A is a predefined action space and S the state space. This sequential construction of objects, together with the fact that the transition policy is parameterised by a neural network (with parameters θ), provides amortisation and the potential of systematic generalisation. Here, we make use of the Trajectory Balance objective [START_REF] Malkin | Trajectory balance: Improved credit assignment in GFlowNets[END_REF], which has proven effective in various tasks [START_REF] Jain | Biological sequence design with GFlowNets[END_REF][START_REF] Hernandez-Garcia | Multi-fidelity active learning with gflownets[END_REF][START_REF] Malkin | Trajectory balance: Improved credit assignment in GFlowNets[END_REF].

Crystals

Crystals are highly structured solid materials defined by a unit cell made of atom types and coordinates that repeat in the three spatial directions, as specified by the crystal's periodic lattice. This ordered structure gives them their unique shapes and properties and their understanding is crucial in the development of new materials. While there exists an infinite number of crystal structures, the field of theoretical crystallography has developed ways to systematically parameterise and classify crystals.

Lattice An n-dimensional lattice Λ can be defined as the set of integral combinations of the linearly independent lattice basis vectors

a i ∈ R n : Λ . = { n i m i a i | m i ∈ Z} .
In 3D, the 3 lattice basis vectors can be converted into 6 parameters: a, b, c determine the lengths of each dimension and α, β, γ the angles between dimensions.

Composition We call composition the set of all atom types present in the crystal's unit cell and their respective quantities. Together with the lattice Λ, the atomic positions provide a complete description of the crystal structure.

Space group A space group is the symmetry group of the crystal's repeating pattern in space, whose elements are the rigid transformations of the pattern that leave it unchanged. A space group element can be described as a tuple (W, t), where W is the linear part of the transformation and t a translation. An element maps a vector x ∈ R n to Wx + t. In three dimensions, space groups are classified into 230 types.

Crystal-GFNs

In this section, we describe the details of Crystal-GFN. Drawing inspiration from crystallography, we propose to represent crystals as the concatenation of three distinct components, or subspaces: space group (SG), composition (C) and lattice parameters (LP) of the unit cell of a crystal. The union of the spaces of these three components would make the state space S = S SG × S C × S LP and sample space X = X SG × X C × X LP of a naive, unrestricted implementation with GFlowNets. An important advantage of the GFlowNet framework and of the domain-inspired data representation that we propose to use is that it allows us to flexibly introduce domain knowledge such as well-studied geometrical constraints from crystallography as well as chemical constraints. These constraints are mostly introduced by restricting the action space A of the Crystal-GFN, both inter-and intra-subspace. We design the Crystal-GFN such that trajectories first select the space group, then the composition and finally the lattice parameters. Below we describe the three subspaces and the constraints we have introduced.

Space group

The sample space of the space group subspace of our Crystal-GFN is X SG = {1, 2, . . . 230}, corresponding to the space of the 230 symmetry groups that exist in three dimensions. As discussed in Section 2.1, GFlowNets rely on the decomposition of objects into multiple steps to facilitate generalisation. Therefore, we draw inspiration from theoretical crystallography to incorporate additional structure into the space group subspace. Crystals are classified into seven crystal systems: triclinic, monoclinic, orthorhombic, tetragonal, trigonal, hexagonal and cubic. Each crystal system comprises a set of point groups. Another related but slightly different category is the lattice system, with also seven systems (see Fig. 3), each comprising a set of Bravais lattices. The lattice system category is convenient for our purposes because it imposes specific constraints on the lattice parameters. Here, we use a derived categorisation resulting from the combination of a crystal system and a lattice system. We refer to it as crystal-lattice system and it has eight categories: triclinic, monoclinic, orthorhombic, tetragonal, trigonal-rhombohedral, trigonal-hexagonal, hexagonal-hexagonal and cubic.

Besides the crystal-lattice system, we introduce another category in the Crystal-GFN space group subspace: the point symmetry (also known as site symmetry), which defines the type of symmetry of a point group. We use the following five categories of point symmetries: centrosymmetric, noncentrosymmetric, enantiomorphic, polar and enantiomorphic-polar. We choose these categories because the combination of point symmetry and crystal system gives rise to one of the 32 crystal classes or crystallographic point groups.

In sum, our space group subspace S SG is a three-dimensional set where the entries correspond to the crystal-lattice system, the point symmetry and the space group. We define the action space such that, starting from the source state s 0 , the transitions can set any of the three dimensions. However, the selection of an option along any dimension restricts the valid options for the remaining dimensions. For example, if point symmetry "non-centrosymmetric" is selected from the source state, only three crystal-lattice systems remain valid and the number of valid space groups is reduced to 25. The "stop" action is valid if and only if the space group has been selected.

Composition

We represent a composition as the number of atoms of each element present in the unit cell. Specifically, for a vocabulary of D elements, where each element can have up to K atoms, we construct D-dimensional vectors where each entry indicates the number of atoms of element d. This yields the space

S C = {(k 1 , . . . k D )|k d ∈ {0, 1, . . . K}, d = 1, . . . D}.
The action space of the composition subspace consists of the choice of element and number of atoms, that is A C = {1, . . . D}×{1, . . . K}, plus the special "stop" action to finish a trajectory.

As per the law of electroneutrality, the sum of positive and negative charge in a solid must be equal. Since we have control over the action space of the GFlowNet, we incorporate a hard constraint to ensure that all the generated compositions can have a neutral charge. We argue that this is an advantage over methods that aim at learning such properties implicitly from the data [START_REF] Xie | Crystal diffusion variational autoencoder for periodic material generation[END_REF].

Finally, we also have an opportunity to incorporate inter-subspace constraints between the space group and the composition. Because of the symmetry imposed by each space group, not all compositions are possible given a space group, and given a composition not all space groups are valid. By way of illustration, the most restrictive space group-international number 230-cannot accommodate compositions with fewer than 16 atoms per element. Since here the trajectories of a Crystal-GFN first sample a space group, we restrict the number of atoms per element in the composition given the space group.

Lattice parameters

The lattice parameters (a, b, c, α, β, γ), in combination with the lattice system (or crystal system) determine the shape of the unit cell in three dimensions. While both the composition and the space group subspaces are discrete, the lattice parameters are real-valued and thus yield a continuous sample space X LP = {(ℓ min , ℓ max ) 3 × (θ min , θ max ) 3 } , where ℓ min and ℓ max are the minimum and maximum lengths and θ min and θ max are the minimum and maximum angles, respectively. In order to sample sets of continuous lattice parameters, the GFlowNet policy model outputs the parameters of a six-dimensional mixture of Beta distributions.

Here, we incorporate an inter-subspace constraint from the space group. Every space group is classified into a lattice system which, in turn, constrains the values of the lattice parameters. For instance, for the cubic lattice system, a = b = c and α = β = γ = 90 • . The remaining constraints are illustrated in Fig. 3. This drastically reduces the state space and removes the need of having to learn such relationships and constraints from the data.

Reward function

As explained in Section 2, our Crystal-GFN is trained to sample proportionally to a reward function. This offers the flexibility of using as reward any quantity of interest, such as a desired property of materials. As an initial proof of concept for the method, we choose the reward to be based on the formation energy (FE) of the sampled crystal structure, in order to generate materials with higher likelihood of being thermodynamically stable. Because the sampled crystal structure may be unknown or not characterised in existing data bases, its true FE may be unknown too. This is why we train a proxy machine learning model to predict the FE given a crystal x ∈ X , which is parameterised as the output of Crystal-GFN. In particular, we train a physics-informed Multi-Layer Perceptron (MLP) on the MatBench [START_REF] Dunn | Benchmarking materials property prediction methods: The matbench test set and automatminer reference algorithm[END_REF] data set. Its architecture is detailed in Appendix D.1.

Empirical evaluation

In this section, we present the results of experiments designed to serve as proof of concept for Crystal-GFN. The motivation of using GFlowNet to generate crystal structures is two-fold: One, we are interested in discovering structures with high scores of a property of interest, in this case low formation energy. At the same time, we want to discover not just one but multiple and diverse structures. The main reason to seek diverse samples is that often the target function is underspecified and the true objective is multifaceted or unknown. A natural way of dealing with underspecification is to try multiple candidates to increase the likelihood of finding successful structures for the downstream applications [START_REF] Jain | GFlowNets for AI-driven scientific discovery[END_REF][START_REF] Hernandez-Garcia | Multi-fidelity active learning with gflownets[END_REF]. Thus, we here aim at analysing whether Crystal-GFN is able to discover diverse crystal structures with low formation energy.

Experimental setup

Exploring the infinitely large space of crystals to find the ones with certain properties is a daunting "needle-in-a-haystack" task. In order to make the haystack a little smaller-though still infinitely large-we restrict the search task to a subset of space groups and compositions, and to a limited range of the lattice parameters, as described in Appendix C. We train the Crystal-GFN for 50,000 iterations, which amounts to 500,000 queries to the proxy model and about 12 hours on a CPU-only machine. Further details about the experimental setup are provided in Appendix C.

Results

In this section, we present results related to our proxy model, to the distribution of the formation energy in Crystal-GFN samples and to the diversity of these samples.

Proxy model of the formation energy As presented in Section 3.4, our reward function is based on a proxy MLP trained on MatBench. Since part of our analysis of the performance of the proposed Crystal-GFN depends on the prediction of the formation energy by this proxy model, it is important to first verify its accuracy. In contrast to other methods, it does not rely on atomic positions to predict the FE of a crystal, but rather on a higher-level description of the crystal, namely its space group, composition and lattice parameters, akin to the outputs of Crystal-GFN. It achieves a mean absolute error (MAE) of 0.10 ± 0.005 eV on the validation set. We further describe its hyper-parameters and detailed performance in Appendices D.2 and D.3.

Predicted formation energy of samples Evaluating generative models is well-known to be a hard, task-dependent problem. Here, in order to gain insights about the sampling policy of the trained Crystal-GFN, we sample 10,000 crystals and compare the distribution with 1) the validation set from MatBench and 2) 10,000 samples from a randomly initialised and untrained Crystal-GFN. Section 4.2 shows the kernel density estimation of the formation energies predicted by our proxy model in the aforementioned sets of samples. As main conclusion, we observe that training the Crystal-GFN shifts the formation energy of the sampled crystals towards lower values, with a median of -3.1 eV, even lower than the median energy in the validation set. Given the vast search space, it is remarkable that with a relatively short training, Crystal-GFN learns to sample from a distribution where 95 % of the structures have a formation energy lower than -2 eV (as predicted by the proxy model). Diversity of samples As discussed before, an important goal in certain materials discovery applications is to find diverse candidates. In other words, sampling crystals with low formation energy would be almost useless for most applications if all the crystals were the same or very similar. To gain insights about the diversity, we analyse the crystal structures sampled by the trained Crystal-GFN. The main conclusion is that Crystal-GFN samples candidates with very high diversity. In particular, regarding the composition, we find that all 12 elements are found in the samples from the Crystal-GFN, and 10 of them appear in the 100 samples with lowest predicted energy. Regarding crystal categories and space groups, we find that 5 out of the 8 crystal-lattice systems and all 5 point symmetries are found in just the top-100 samples. While not all 113 space groups are found in the 10,000 samples, we find 80 (70 %) of them. Finally, in terms of lattice parameters, we also observe relatively similar distributions of lengths and angles, compared to the MatBench data set. Further details and visualisations of these results are provided in Appendix E.

Conclusions and future work

In this paper, we have introduced a new generative model, Crystal-GFN, to sample inorganic crystal structures proportionally to a property of interest. A key feature of Crystal-GFN is its flexibility to be trained with any available reward function and to incorporate domain constraints. The latter is due to the fact that Crystal-GFN constructs crystals sequentially in the space of space groups, compositions and lattice parameters of the unit cell. This parameterization, inspired by theoretical crystallography, has allowed to flexibly incorporate constraints regarding the neutral charge of the composition, the compatibility between composition and space group and between space group and lattice parameters.

By training a Crystal-GFN with a reward function based on the formation energy predicted by a proxy model trained on the MatBench database, we have shown it learns to sample candidates with low predicted formation energy (lower mean and median that the distribution of the MatBench data set). Importantly, it does so by sampling highly diverse crystals, in terms of coverage of space groups, compositions and lattice parameters.

Interesting directions for future work include extending Crystal-GFN with more domain-inspired constraints as well as with additional subspaces to sample the atom positions in the unit cell. Alternatively, crystal structure prediction from the samples could be performed using pre-trained machine learning models that generate atomic coordinates. Furthermore, it would be interesting to explore other properties of interest beyond the formation energy of the sampled crystals.

A Related Work

Material and Molecule generation. There exists a myriad of generative methods for molecular and crystal generation, generally falling into one of these catagories: Variational Auto Encoders (VAE) [START_REF] Ren | An invertible crystallographic representation for general inverse design of inorganic crystals with targeted properties[END_REF][START_REF] Pakornchote | Diffusion probabilistic models enhance variational autoencoder for crystal structure generative modeling[END_REF], Generative Adversarial Networks (GANs) [START_REF] De | Molgan: An implicit generative model for small molecular graphs[END_REF][START_REF] Asma Nouira | Crystalgan: Learning to discover crystallographic structures with generative adversarial networks[END_REF][START_REF] Long | Constrained crystals deep convolutional generative adversarial network for the inverse design of crystal structures[END_REF], Normalizing flows [START_REF] Garcia Satorras | E(n) equivariant normalizing flows[END_REF][START_REF] Ahmad | Free energy calculation of crystalline solids using normalizing flows[END_REF] and Diffusion models [START_REF] Zheng | Towards predicting equilibrium distributions for molecular systems with deep learning[END_REF][START_REF] Xu | Geodiff: A geometric diffusion model for molecular conformation generation[END_REF]. In general, these methods are trained to reconstruct the training data from a latent distribution by maximising its likelihood or minimising the discrepancy between generated sampled and the real distribution. Once trained, they can generate new instances by sampling from the learned latent space. [START_REF] Hoffmann | Data-driven approach to encoding and decoding 3-d crystal structures[END_REF][START_REF] Long | Constrained crystals deep convolutional generative adversarial network for the inverse design of crystal structures[END_REF][START_REF] Callum | 3-d inorganic crystal structure generation and property prediction via representation learning[END_REF] treat materials as 3D voxel images, but the process of decoding images back to atom types and coordinates often results in low validity, and the models are not rotationally invariant. Several methods [START_REF] Ren | Inverse design of crystals using generalized invertible crystallographic representation[END_REF][START_REF] Kim | Generative adversarial networks for crystal structure prediction[END_REF][START_REF] Zhao | High-throughput discovery of novel cubic crystal materials using deep generative neural networks[END_REF] encode unit cell related information to generate crystals. GANCSP [START_REF] Kim | Generative adversarial networks for crystal structure prediction[END_REF] and CubicGAN [START_REF] Zhao | High-throughput discovery of novel cubic crystal materials using deep generative neural networks[END_REF] use as input fractional coordinates, element properties, and lattice parameters, to build models that generate crystals conditioned on composition or both composition and space group. Such models are generally not invariant to any Euclidean transformations. One close baseline, CDVAE [START_REF] Xie | Crystal diffusion variational autoencoder for periodic material generation[END_REF], employs a diffusion model in conjunction with a VAE to sample new materials, leveraging lattice parameters. However, it is does not enable one to impose domain-expertise based constraints into the generation process. [START_REF] Cheng | Crystal structure prediction by combining graph network and optimization algorithm[END_REF] leveraged physics and symmetry-based constraints to solve an optimization problem over crystal graphs to obtain atomic coordinates, but requires a composition and lattice information as input.

Finally, auto-regressive methods like G-SchNet [START_REF] Gebauer | Symmetry-adapted generation of 3d point sets for the targeted discovery of molecules[END_REF] or G-SphereNet generates 3D molecules sequentially, placing atom-by-atom. They are unaware of periodicity and cannot generate the lattice, making them ill-suited for crystals. This also applies to flow models [START_REF] Shi | Graphaf: a flow-based autoregressive model for molecular graph generation[END_REF].

As these models directly sample 3D positions, they struggle to embed essential physical constraints 3 , often leading to physically unrealistic and unstable structures. Unlike them, by selecting the composition, lattice, and space group and imposing validity constraints, our approach inherently incorporates domain knowledge about materials science and crystallography ensuring that the generated crystals adhere to fundamental principles and constraints, resulting in more physically valid structures.

GFlowNets: GFlowNets are a new [START_REF] Bengio | Flow network based generative models for non-iterative diverse candidate generation[END_REF] class of trained generative models that can be used to sample objects through a series of constructive steps, allowing for high reward samples and also diversity in sampling. They have the potential of being employed in a variety of different materials discovery avenues [START_REF] Jain | GFlowNets for AI-driven scientific discovery[END_REF], starting from sample design to experiment design. In particular, they have already shown promise in the area of drug-discovery, having been used to sequentially construct small molecules [START_REF] Bengio | Flow network based generative models for non-iterative diverse candidate generation[END_REF] and biological sequences such as peptides and DNA/protein sequences [START_REF] Jain | Biological sequence design with GFlowNets[END_REF].

B Crystals

(a) (b) (c) (d) 2 1 0 (f)
Figure 3: The seven lattice systems and the constraints they impose on the lattice parameters. In order a-g: triclinic, monoclinic, orthorhombic, tetragonal, rhombohedral, hexagonal and cubic. Source of the figures: Wikimedia Commons, licensed under the Creative Commons Attribution-Share Alike 3.0 Unported license.

Materials encompass the substances that make up objects. Their distinctive characteristics arise from a combination of factors, including their composition, structure, and manufacturing processes.

While materials come in various forms, our primary focus will be on solid-state materials, specifically crystals, owing to their significance in many machine learning applications we have encountered.

Solid state materials consist of a multitude of smaller building blocks, such as atoms or molecules, securely positioned in fixed locations. They exhibit high-density packing and strong mutual attraction, resulting in a stable structure with a well-defined volume. Materials that possess a well-ordered, infinitely repetitive structure are categorized as crystals, while those lacking such long-range positional order are termed amorphous.

Crystals, also known as crystalline materials, represent a category of solid substances characterised by their periodic structure. This entails a consistent predictable repetition of atom or molecule placement within the crystal across spatial dimensions. It is this inherent periodicity that sets crystals apart from smaller molecules and confers unique properties upon crystals. They manifest in diverse forms and showcase exceptional attributes like transparency and elevated melting points. The study of crystals offers a pathway to precise characterization and controlled manipulation, opening up avenues for engineering materials with targeted properties to fulfill a wide range of needs. The subsequent section delves into how machine learning models address and harness this periodicity.

A crystal structure includes the type and position of every atom as well as the translational axes which allow the structure to repeat. Since crystal structures extend infinitely, we need to define some small section that can be repeated. We call this small section a unit cell. Unit cells are made of a lattice and basis. Lattices tell you how the crystal is repeated, the basis tells you what is repeated. There are six parameters of a unit cell: 3 edges a, b, c and the angles between the edges α, β, γ. The edges of a unit cell may be or may not be perpendicular to each other. A crystal structure combines the lattice (how atoms are translated) with a basis (which atoms are translated) which theoretically describes every atom in the crystal. The lattice is just a mathematical idea, but the basis is the "real" component which is repeated according to the mathematical lattice.

Lattice The lattice is fully specified by its basis vectors a i and is associated with a group of translations T Λ . The lattice Λ is the underlying structure of the unit cells, which have the property of tilling the space R n when translated by lattice vectors. Of particular prominence is the primitive cell U for the basis a i :

U . = { n i x i a i | 0 ≤ x i < 1} .
In a material, the unit cell comprises a set of atomic positions S = {(Z i , x i ) | x i ∈ U }, where x i is the atom position (called the basis vector) and Z i its atomic number. S can contain an arbitrary number of atoms and is not constrained to a particular structure.

Periodic boundary conditions (pbc) are commonly used to represent the infinite repeating nature of crystals and other periodic systems. Rather than confining the analysis to a single unit cell housing all atoms, pbc envision an infinite array of identical unit cells extending indefinitely in all directions. This means that if an atom or molecule exits the cell on one side, it effectively re-enters from the opposite side, as if traversing a continuous lattice. This approach enables scientists to simulate the behavior of an infinitely extending crystal using a finite computational domain.

C Experimental setup

In this section, we provide additional details about the hyper-parameters and experimental setup used to train the Crystal-GFN to obtain the results presented in Section 4.

In order to reduce the search space in our experiments, we apply the following restrictions:

• Compositions consist of up to 5 different elements from the subset of these 12 elements: H, Li, C, N, O, F, Mg, Si, P, S, Cl and Fe. These are the most common elements in the MatBench data set used to train the proxy model.

• Compositions can contain up to 50 atoms in total and up to 16 atoms per element. This number is obtained by finding the lowest Wyckoff position multiplicity for each space group, and then computing the maximum of these values across all space groups. This means that 16 is the lowest possible value that still makes it possible for the Crystal-GFN to generate samples from all space groups while respecting their symmetry constraints.

• Space groups are the intersection of train and validation space groups, from the proxy model's data set. There are 113 of them.

• The minimum and maximum lengths of the unit cell are 0.9 and 100 angstroms, respectively; the minimum and maximum angles are 50 • and 150 • , respectively.

We train the Crystal-GFN by sampling 10 trajectories per iteration from the current forward policy.

In order to encourage further exploration, 10 % of the steps in the trajectories are sampled at random from a uniform distribution. In total, we train for 50,000 iterations, which amounts to 500,000 queries to the proxy model. This took about 12 hours on a CPU-only machine.

The architecture of both the forward and the backward GFlowNet policy models is a 3-layer multilayer perceptron with 256 units per layer. We trained with the Adam optimiser and a learning rate of 0.0001. As is common with the Trajectory Balance objective [START_REF] Malkin | Trajectory balance: Improved credit assignment in GFlowNets[END_REF], we set a higher learning rate (0.01) for the 16 learning weights used to parameterised partition function. In order to sample structures with lower (negative) formation energy, we set a temperature of 8 in the reward function defined in Section 3.4.

The distribution to sample the increments of the lattice parameters subspace is a mixture of 5 Beta distributions. In order to ensure numerical stability during training, we restrict the values of the coefficients of the Beta distributions to the range [0.1, 100]. One of the conditions that must be satisfied by generalised GFlowNets is that trajectories must have finite length. To this end we set a minimum increment of 10 % the range of each dimension.

D Proxy MLP D.1 Architecture

For a given crystal structure, the input to the proxy MLP is the concatenation of: a physical embedding of the crystal's elements using PhAST [START_REF] Duval | Phast: Physics-aware, scalable, and task-specific gnns for accelerated catalyst design[END_REF], a learned embedding for its space group and the standardized lattice parameters. We stratify the MatBench data set into train, validation and test sets, controlling for the distribution in target FE. The final reward function is a Boltzmann version of the proxy: R(x) = exp -MLP(x)

T
where T is a temperature hyper-parameter. This ensures that lower FE yields higher positive reward, and we can control our preference for lower energies with the temperature T .

Overall, the architecture consists of concatenated representations for the composition (h C ), the space group h SG and the lattice parameters (h LP ) that are used as input to a prediction MLP ŷ = MLP out ([h C , h SG , h LP ]). In the following we describe how each h is obtained from the data:

• h C : Each element Z i in the composition is embedded with PhAST [START_REF] Duval | Phast: Physics-aware, scalable, and task-specific gnns for accelerated catalyst design[END_REF]: it is a concatenation of 1/ static physical properties P Z projected into some latent space with a small MLP 2/ a learned embedding for the atomic number Z, Period P and the Group G of the element, noted E i :

h C = MLP C MLP P (P Z ), E Z , E P , E G (1)
• h SG : Each space group is embedded as a lookup in a learnable embedding matrix.

• h LP : Lattice parameters are first standardize element-wise using the training set statistics. This 6-dimensional vector is then passed through a small MLP:

h LP = MLP LP LP(x) -µ train (LP) σ train (LP) (2)

D.2 Performance

Our proxy MLP model achieves an overall Mean Absolute Error of 0.10 ± 0.005 eV4 . For reference, best performing methods in Matbench leaderboard5 can achieve up to 10× better performance.

Figure 1 :

 1 Figure 1: A schematic of the crystal generation process of Crystal-GFN.

Figure 2 :

 2 Figure2: A comparison of the distributions of the formation energy predicted by our proxy model in three relevant distributions of samples: in blue, samples from Crystal-GFN after training; in orange, the validation set, representative of the MatBench database; in pink, samples from an untrained Crystal-GFN. As a main conclusion, we observe that Crystal-GFN, after training, manages to sample crystals with predicted formation energies in the range of the validation set.

Figure 5 :

 5 Figure 5: Distribution of FE values in the validation set and associated Proxy MLP MAE, with 25%, 50% and 75% MAE quantiles.

Figure 6 :

 6 Figure 6: Distribution of total element prevalence per data set (10k Crystal-GFN samples, MatBench val, MatBench train). Each element is counted as per its stoichiometry in the

Figure 7 :

 7 Figure 7: Distribution of binary element prevalence per data set (10k Crystal-GFN samples, MatBench val, MatBench train). Each element is counted only once per crystal.

Figure 10 :

 10 Figure 10: Distribution of lattice parameters

Table 1 :

 1 Proxy MLP hyper parameters

	Hyper parameter	Value	Description
	properties_proj_size	64	Projection size of atomic properties h = Wp Z
	group_emb_size	16	Embedding size of element's group
	period_emb_size	256	Embedding size of element's period
	z_emb_size	128	Embedding size of element's atomic number
	sg_emb_size	128	Crystal space group embedding size
	lat_hidden_channels	284	Hidden channels for layers processing the lattice parameters
	lat_num_layers	1	Number of hidden layers for the lattice parameters MLP
			Number of layers for final MLP processing composition,
	num_layers	5
			space group and lattice parameters hidden representations
	hidden_channels	576	Size of final MLP hidden layers.
	optimizer	Adam
	lr	0.0017	Learning rate
	batch_size	448
			Divide learning rate by 2 when validation MAE
	scheduler	ReduceLROnPlateau
			does not improve for 4 epochs
	es_patience	11	Early stopping patience (in epochs)
	epochs	100	Total epochs (effective after early stopping: 97)

like energy minimization, dihedral angles geometric constraints or bonding preferences between atom types.

95% confidence interval, computed modeling the MAE as a log-normal distribution. This is further justified by Fig.

[START_REF] Callum | 3-d inorganic crystal structure generation and property prediction via representation learning[END_REF] https://matbench.materialsproject.org/Leaderboards%20Per-Task/matbench_v0.1_ matbench_mp_e_form/

Acknowledgments and Disclosure of Funding

We thank Yasmine Benabed for her contribution to the conception and first steps of this project. We also thank Bruno Rousseau, Simon Blackburn and Mickael Dollé for their valuable insights from the materials science perspective. Sékou-Oumar Kaba and Michael Kilgour have also participated in helpful discussions. We thank Mila's IDT team for their support. Finally, we thank Quebec's Ministry of Economy, Innovation and Energy for their financial support. This research was enabled in part by compute resources provided by Mila (mila.quebec).

However, they all use 3D positions of the crystals, and remarkably, our approach outperforms some of them, indicating that leveraging composition, space-group and lattice parametrs serves as a viable representation in the formation energy estimation task.

In Figures 4 and5 we detail further the performance of our Proxy MLP on the MatBench validation set. In particular we observe that it maintains good performance for structures with FE up to ∼ 0.3ev. The MAE beyond that FE level becomes high and unstable, but this is expected as those are the regions of the target space with least data available. Additionally we can verify in Figure 4 that our stratification algorithm works as expected, with similar FE distributions in the validation and train splits.

D.3 Hyper parameters

We present the hyper parameters of our Proxy MLP in Table 1 along with a description of their role in the architecture.

E Additional results

In this section, we provide additional results demonstrating the diversity of the samples generated by the Crystal-GFN. They are displayed in Figure 6, Figure 8, Figure 9 and Figure 10.

Space group

Crystal-GFlownet MatBench Val MatBench Train


Alvaro Carbonero
email: alvaro.carbonero@mila.quebec

Alexandre Duval

Victor Schmidt

Santiago Miret

Alex Hernandez-

Garcia Mila

Yoshua Bengio

David Rolnick

On the importance of catalyst-adsorbate 3D interactions for relaxed energy predictions

The use of machine learning for material property prediction and discovery has traditionally centered on graph neural networks that incorporate the geometric configuration of all atoms. However, in practice not all this information may be readily available, e.g. when evaluating the potentially unknown binding of adsorbates to catalyst. In this paper, we investigate whether it is possible to predict a system's relaxed energy in the OC20 dataset while ignoring the relative position of the adsorbate with respect to the electro-catalyst. We consider SchNet, DimeNet++ and FAENet as base architectures and measure the impact of four modifications on model performance: removing edges in the input graph, pooling independent representations, not sharing the backbone weights and using an attention mechanism to propagate non-geometric relative information. We find that while removing binding site information impairs accuracy as expected, modified models are able to predict relaxed energies with remarkably decent MAE. Our work suggests future research directions in accelerated materials discovery where information on reactant configurations can be reduced or altogether omitted.

Introduction

Materials discovery is a strong driver of innovation, unlocking new materials with tailored properties that benefit various domains including energy efficiency, transportation systems and electronics [START_REF] Keith | Machine learning for molecular and materials science[END_REF]]. Yet, it faces significant hurdles. Characterizing new materials is computationally expensive [START_REF] Artem R Oganov | Structure prediction drives materials discovery[END_REF], even when replacing lab experiments by quantum mechanical simulations like the Density Functional Theory (DFT) [START_REF] Kohn | Density functional theory of electronic structure[END_REF]. Besides, exploration of materials is hindered by the vastness of the search space, which encompasses myriad compositions, atomic arrangements and properties [START_REF] Edward O Pyzer-Knapp | Accelerating materials discovery using artificial intelligence, high performance computing and robotics[END_REF].

To overcome these challenges, researchers have turned to Machine Learning (ML) for accelerated materials discovery [START_REF] Zhang | A survey on graph diffusion models: Generative ai in science for molecule, protein and material[END_REF], Miret et al., 2022[START_REF] Long | Matsciml: A broad, multi-task benchmark for solid-state materials modeling[END_REF][START_REF] Song | MatSci-NLP: Evaluating scientific language models on materials science language tasks using text-to-schema modeling[END_REF] for two primary reasons: First, ML holds the power to quickly model materials' properties (i.e., to evaluate candidates) using Geometric Graph Neural Networks (GNNs). Second, generative ML can automatically propose new and consistent material candidates.

In recent years, ML has emerged as a crucial tool for the discovery of electro-catalysts, which play a key role in promoting renewable energy processes and sustainable chemical production, including the production of ammonia for fertilizers and hydrogen [START_REF] Lawrence Zitnick | An introduction to electrocatalyst design using machine learning for renewable energy storage[END_REF]. The Open Catalyst Project (OCP) [START_REF] Chanussot | Open catalyst 2020 (oc20) dataset and community challenges[END_REF] has significantly contributed to this field by releasing an extensive dataset of 1,281,040 DFT relaxations of catalyst-adsorbate pairs, selected from a pool of meanigful candidates1 . This dataset was specifically designed to train ML models to predict the relaxed energy of 3D adsorbate-catalyst (adslab) systems, a critical property that influences a catalyst's activity and selectivity [START_REF] Tran | Active learning across intermetallics to guide discovery of electrocatalysts for co2 reduction and h2 evolution[END_REF], or its effectiveness for a specific chemical reaction. 2 . OCP has facilitated significant advancements in catalysis discovery, with ML models increasingly bridging the performance gap with DFT simulations, while offering a speed advantage of several orders of magnitude [START_REF] Musaelian | Learning local equivariant representations for large-scale atomistic dynamics[END_REF][START_REF] Gasteiger | Gemnet-oc: developing graph neural networks for large and diverse molecular simulation datasets[END_REF][START_REF] Passaro | Reducing so (3) convolutions to so (2) for efficient equivariant gnns[END_REF].

While research has primarily focused on predicting material properties [START_REF] Wieder | A compact review of molecular property prediction with graph neural networks[END_REF]], it's equally crucial to efficiently explore the vast search space of potential catalyst candidates [START_REF] Wei | Machine learning in materials science[END_REF][START_REF] Zhang | Machine learning: a new paradigm in computational electrocatalysis[END_REF]. This involves generating consistent candidates automatically, a challenge due to the intricate process of creating adslab samples. This process, further detailed in Appendix A, includes positioning an adsorbate molecule (e.g., H2O) with a catalyst in 3D space, meaning cutting a surface through the catalyst bulk, selecting the adsorbate's spatial orientation, and sampling a plausible binding site [START_REF] Chanussot | Open catalyst 2020 (oc20) dataset and community challenges[END_REF]. These steps, which are time-consuming and challenging to model, determine the input configuration of the adslab sample and significantly influence the relaxed energy prediction [START_REF] Lan | Adsorbml: Accelerating adsorption energy calculations with machine learning[END_REF]. An additional challenge is that this process and its relevance to the actual in-lab efficiency of a real material is not yet fully understood by chemists [START_REF] Deshpande | Graph theory approach to determine configurations of multidentate and high coverage adsorbates for heterogeneous catalysis[END_REF].

In light of these challenges, this paper explores the possibility of predicting the relaxed energy of an adslab without co-locating the adsorbate and the catalyst in the same 3D space. This direction of study is valuable for several reasons: (1) to better understand the role of adsorbate-catalyst geometry in determining the relaxed energy; (2) to reduce reliance on the exact input configuration, which often correlates with local energy minima, (3) to avoid the complexity and computational cost of determining a good input configuration (e.g., finding a binding site and adsorbate orientation). Our primary objective is quantify the loss in accuracy that occurs when the geometric relationship between inputs is unavailable.

To achieve this, we propose four modifications of existing GNN architectures, collectively referred to as Disconnected GNNs. All four modifications enable the base architecture to make relaxed energy predictions without explicitly modeling geometric interactions between the adsorbate and the catalyst. We then evaluate the trade-off of omitting these interactions through experiments on the OC20 dataset and suggest future directions to overcome the limits associated with the input adslab configuration.

Methods

Our methods test the assumption that meaningful results can be obtained in predicting the relaxed energy even when the geometric interactions between the adsorbate and the catalyst are unknown. Specifically, we propose four Disconnected GNN models. These models can leverage as their backbone any underlying GNN model that predict properties of 3D atomic systems using the general pipeline detailed under "Refresher" below, where we make targeted changes to the pipeline. In this paper, we investigate the effect of the proposed Disconnected GNN architecture on three backbone models: SchNet [START_REF] Schütt | Schnet: A continuous-filter convolutional neural network for modeling quantum interactions[END_REF], DimeNet++ [START_REF] Klicpera | Fast and uncertaintyaware directional message passing for non-equilibrium molecules[END_REF] and FAENet [START_REF] Duval | Faenet: Frame averaging equivariant gnn for materials modeling[END_REF].

Refresher. Recall that GNNs applied to 3D atomic systems mostly use the following pipeline. (1) Graph creation: construct a graph representation of the input point cloud systems (e.g. adslab), represented with atomic numbers and 3D atom positions. (2) Embedding: derive node/edge embeddings based on atomic numbers and geometric information (e.g. relative atom positions). (3) Interaction blocks: apply a fixed number of message passing layers [START_REF] Gilmer | Neural message passing for quantum chemistry[END_REF] to update the node and/or the edge embeddings, using geometric information and preserving data symmetries. (4) Output block: project final atom representations into scalar values (e.g. atom-properties). For graph-level prediction, a global pooling (e.g. sum) is performed to aggregate atom predictions. Figure 1: Architectures are presented in order of increasing complexity. 1a illustrates the standard property prediction pipeline for Geometric GNNs, with no modifications. In 1b, edges between parts are removed. In 1c, the pooling step is performed separately on atom representations for the adsorbate and the catalyst. Resulting embeddings are concatenated and passed through an MLP. In 1d, each part has its own model with distinct weights. In 1e, an interaction block with artificial edges and weights allows node embeddings between both parts to communicate.

Disconnected baseline

The baseline Disconnected GNN model does not create any edges between the catalyst and the adsorbate in the graph creation step. By removing all such edges, no relative geometric information will pass from the adsorbate to the catalyst when modelling the system, since GNNs only propagate information through graph edges. Since this modification only relates to the input data representation, this change can be applied to any backbone GNN architecture, like FAENet, SchNet or DimeNet++ for instance. It is illustrated in Figure 1b.

Note that we still use the original OC20 dataset where the adsorbate/catalyst atom positions have been determined after selecting a binding site and an adsorbate orientation. But since GNNs use atom relative positions and not absolute atom positions to preserve translation equivariance, disconnecting the adsorbate-catalyst graph is enough to discard totally the input configuration. Processing such a disconnected graph is equivalent to considering both components independently. Associated results thus inform about the importance of the relative adslab geometric information to each model.

Independent pooling

This method builds upon the disconnected baseline defined above, also removing all edges between the adsorbate and the catalyst. As illustrated in Figure 1c, there are three additional changes. First, the projection block now outputs a hidden representation h ∈ R ⌊H/2⌋ for each node instead of a scalar quantity, where H is the dimension of the embedding. Second, the resulting node representations are pooled separately for each component, leading to distinct hidden representations for the adsorbate and the catalyst. Lastly, we concatenate these two representations and pass them to an additional dense layer to compute the energy of the system. This final MLP layer gives more expressive power to the model than the disconnected baseline because we explicitly model non-geometric interactions between the adsorbate and the catalyst.

Independent backbones

This approach, illustrated in Figure 1d, builds upon independent pooling. The main difference is that catalysts and adsorbates are assigned distinct GNN models (hence "independent") instead of sharing backbone weights. The motivation is that the adsorbate and the catalyst have very different sizes and roles. Hence, it might be beneficial to have independent GNNs, where each one is trains for its particular input. Similarly to independent pooling, both models are modified to produce a hidden representation of their respective parts. They are then concatenated and passed through an MLP.

Attention

This model (see Figure 1e) builds upon the independent model by allowing the nodes of both graphs to communicate during the interaction layers.

We modify the graph creation process to produce a heterogeneous graph where we add the following weighted edges to previous models' edges. For every catalyst node i and adsorbate node j, we create the new edges (i, j) and (j, i) with weights z and -z respectively, where z is the z-axis coordinate of x i . Besides, in the embedding step, we pass the edge weights through an RBF layer and a linear layer to increase their dimension. Notice that this model can run without locating the adsorbate and catalyst in the same 3D plane since no information regarding the adsorbate's location (e.g. distance, relative position) is given. The motivation behind these weighted edges is to make the model aware of the closeness of a node to the catalyst's surface, while respecting the restrictions of the setup. As noted in [START_REF] Duval | Phast: Physics-aware, scalable, and task-specific gnns for accelerated catalyst design[END_REF], fixed atoms not near the surface of the catalyst (tag 0 nodes) can be removed without hurting efficiency, indicating that closeness to surface plays a factor.

We now describe how this graph is used. Let F denote a "normal" intra-interaction layer, let h l ads and h l cat denote the node embeddings of the adsorbate and the catalyst at interaction layer l, and let e ads and e cat stand for the edges of the adsorbate and catalyst, respectively. To produce h l+1 ads and h l+1 cat , we first calculate h ′ ads = F (h l ads , e ads ) and h ′ cat = F (h cat , e cat ). Subsequently, we calculate

h ′ = GATConv(Concat(h ′ ads , h ′ cat ), e ′
) where e ′ is the set of new weighted edges and GATConv is a convolution layer as described in [START_REF] Veličković | Graph attention networks[END_REF]. Finally, Metrics: We measure accuracy on each validation split through the energy MAE. Running time is measured with the throughput at inference time, i.e. the average number of samples per second that a model can process in its forward pass, on similar GPU types. We can see the the disconnected baseline consistently achieves the best performance. In terms of throughput, it is on par with other modifications for DimeNet++ and FAENet, while it is slower than all other modifications for SchNet (except for the Attention modification).

h l+1 ads , h l+1 cat = Norm(h ′ ads + h l ads ), Norm(h ′ cat + h l cat ).
Baselines: We compare base results of the original models (considering 3D adsorbate-catalyst interactions) with the four disconnected architectures detailed in Section 2. Note that we apply PhAST components [START_REF] Duval | Phast: Physics-aware, scalable, and task-specific gnns for accelerated catalyst design[END_REF] on all models to improve performance and inference time.

Results. From Figure 2, we see that the independent model architecture performs best among all proposed disconnected GNNs. However, it presents a performance drop in MAE of 20%, 31%, and 28% with respect to baseline architectures, for Schnet, DimeNet++ and FAENet respectively. A complete numerical table of results is reported in Table 1.

Performance is significantly worse for disconnected GNNs than for GNNs with access to the adslab configuration; however, these results are promising in several ways.

Firstly, a drop in performance was expected since we neglect important information by omitting to modelise the adsorbate-catalayst 3D interaction. The important take-away is that our disconnected GNN manages to produce decent energy predictions while being agnostic to the adslab's input configuration. Future work will probably manage to reduce further this performance gap. Note that we discuss some axes of improvements in Appendix C.

Secondly, the OC20 dataset often contains multiple input configuration samples for the same catalystadsorbate pair, each with a different relaxed energy target. This cannot be captured by our disconnected GNN, leading to multiple targets for the same sample, which can impair training dynamics. Furthermore, when evaluating on a dataset that does not present this feature, e.g. the 10k split of OC20 IS2RE, disconnected GNNs match (or outperform) baseline models, which is very promising. We present and briefly discuss these results in Appendix D.

The OC20-Dense dataset [START_REF] Lan | Adsorbml: Accelerating adsorption energy calculations with machine learning[END_REF] was recently released, containing multiple input configurations (e.g. binding site, adsorbate orientation) for each adsorbate-catalyst pair along with the associated relaxed energy. Building upon this dataset, we will be able to train disconnected GNN models to predict the minimum relaxed energy of each adslab over all possible input configurations.

The obvious benefit of such approach is that it marginalises over adslab binding sites and adsorbate orientations, avoiding the need to explore all configurations to actually find the relaxed energy global minima of this adsorbate-catalyst pair, which is expected to be the main reaction driver in real life experiments [START_REF] Lan | Adsorbml: Accelerating adsorption energy calculations with machine learning[END_REF].

Conclusion

In this paper, we have explored predicting an adsorbate-catalyst system's relaxed energy without co-locating them in 3D space. We proposed four disconnected models built on top of existing Geometric GNNs architectures, and showed through experiments that this task was possible with a manageable performance loss. Looking forward, we see the potential of using such disconnected models to find the global energy minima of adslab systems while circumventing the need to consider all possible input adslab configurations. This sounds particularly promising for generative methods and should contribute to the acceleration of catalysis discovery. set of 82 molecules that are chosen for their utility to renewable energy applications. For the catalyst surface, the process is divided in three stages. First, we choose the number of distinct chemical elements composing it. It can be a unary material (5% chance), a binary material (65%) or a ternary material (30%). These elements are chosen from a set of 55 elements comprising reactive nonmetals, alkali metals, transition metals, etc. Next, a stable bulk material is randomly selected from the 11K samples of the Materials Project [START_REF] Jain | The materials project: A materials genome approach to accelerating materials innovation[END_REF] with the number of elements chosen in the first step. Lastly, all symmetrically distinct surfaces3 from the material with Miller indices [START_REF] Neil | Solid state physics[END_REF] less than or equal to 2 are enumerated, including possibilities for different absolute positions of surface plane. From this list of surfaces, one is randomly selected.

Input adslab configuration. The objective of the second step is to place the adsorbate and the catalyst surface in the same 3D plane. Using pymatgen's Voronoi tesselation algorithm [Ong et al., 2013] on surface atoms and adsorbate possible binding sites, Catkit [START_REF] Jacob R Boes | Graph theory approach to high-throughput surface adsorption structure generation[END_REF] enumerates a list of symmetrically distinct adsorption sites along with suggested per-site adsorbate orientations. From this list, an adsorption configuration is randomly selected, yielding the initial adslab structure. This initial structure is then relaxed using Density Functional Theory simulations, performed using the Vienna Ab Initio simulation Package (VASP) [START_REF] Kresse | Ab initio molecular-dynamics simulation of the liquid-metalamorphous-semiconductor transition in germanium[END_REF].

B Full results table C Axes of improvement

There are two main ways in which Disconnected GNN models can be further improved. The first one is to improve upon dense layers as a way to combine the final hidden representations of the adsorbate and the catalyst. The second one is to find a way to effectively allow the nodes of the adsorbate and catalyst to communicate, and to do this while ommitting the relative position and orientation between the adsorbate and the catalyst. This is attempted in the attention models through the inter-interaction layers. Given that they do not consistently outperform the independent backbone models, it is likely that the way in which attention models use GAT convolution layers hinders the model's prediction.

We still hope that this can be fixed since the attention models can outperform baseline models when trained on the 10k dataset, as we described in Appendix D. Table 2: MAE and inference time for various GNNs and their disconnected counterpart on OC20 10k IS2RE, averaged over 3 runs. Average MAE is computed over all validation splits. Disconnected models all show lower accuracy. In every model, every subsurface catalyst atom (tag 0 atoms) is removed as suggested in [START_REF] Duval | Phast: Physics-aware, scalable, and task-specific gnns for accelerated catalyst design[END_REF] to improve inference time.

D Results in the 10k dataset

Only 89% of the OC20 IS2RE dataset consists of adsorbate-catalyst systems which are unique up to the combination of adsorbate id, bulk id, and bulk cell. That is, to Disconnected GNNs, 11% of the dataset consists of inputs with multiple targets due to distinct binding sites, which impairs training dynamics. When compared to the OC20 10k IS2RE dataset, every single adsorbate-catalyst system has a unique combination of adsorbate id, bulk id, and bulk cell. Thus, evaluating disconnected GNNs on the 10k dataset is a useful way to assess their performance on a suitable data tailored to them. This is corroborated by the results in Table 2.

From Table 2, we see that independent models perform on par or better (for FAENet) than baseline models when trained on the 10k dataset. In particular, we see a performance improvement of 18%, -4%, and 5% for SchNet, DimeNet++, and FAENet respectively. Different than in the all dataset, Attention-FAENet outperforms Independent models-FAENet by 7.5%. These results are a sign that Disconnected GNNs can perform well when tested against the data for which they were designed for. The OC20-Dense dataset [START_REF] Lan | Adsorbml: Accelerating adsorption energy calculations with machine learning[END_REF] is thus an exciting development as it will allow for the creation of the datasets that these models need.

  We use the OC20 IS2RE[START_REF] Lawrence Zitnick | An introduction to electrocatalyst design using machine learning for renewable energy storage[END_REF] dataset, which involves the direct prediction of the relaxed adsorption energy from the initial atomic structure, i.e. a graph regression task requiring E(3)-invariance. OC20 comes with 450K training samples and a predefined train/val/test split. We evaluate models on the four distinct splits of the validation set (∼ 25K samples each): In Domain (ID), Out of Domain Adsorbates (OOD-ads), Out of Domain catalysts (OOD-cat), and Out of Domain Adsorbates and catalysts (OOD-both).

Figure 2 :

 2 Figure2: Comparison of MAE performance drops across our four proposed modifications to three base architectures (left), and their associated impact on model throughput (right). We can see the the disconnected baseline consistently achieves the best performance. In terms of throughput, it is on par with other modifications for DimeNet++ and FAENet, while it is slower than all other modifications for SchNet (except for the Attention modification).

Table 1 :

 1 MAE and inference time for various GNNs and their disconnected counterpart on OC20 IS2RE, averaged over 3 runs. Average MAE is computed over all validation splits. Best results are shown in bold, second best underlined. Overall, disconnected models show a lower accuracy than baselines but they are able to predict relaxed energy relatively well despite omitting geometric interaction between parts in the adslab.

	Baseline / MAE	ID	OOD-ad OOD-cat OOD-both Average Throughput (samples/s)
	SchNet	0.631	0.687	0.625	0.626	0.642 0.642 0.642	16001 ± 2264
	Base disconnected-SchNet	0.940	1.037	0.932	0.947	0.964	18863 ± 3034 18863 ± 3034 18863 ± 3034
	Independent pooling-SchNet	0.766	0.847	0.762	0.789	0.791	18590 ± 2230
	Independent models-SchNet	0.761	0.805	0.764	0.764	0.773	12220 ± 1593
	Attention-SchNet	0.779	0.852	0.762	0.797	0.797	5652 ± 876
	DimeNet++	0.577	0.693	0.568	0.621	0.615	631 ± 125
	Base disconnected-DimeNet++	0.922	1.059	0.914	0.989	0.971	667 ± 143 667 ± 143 667 ± 143
	Independent pooling-DimeNet++ 0.781 0.0.992	0.767	0.901	0.860	650 ± 139
	Independent models-DimeNet++ 0.767	0.884	0.754	0.818	0.806	374 ± 68
	Attention-DimeNet++	0.777	0.952	0.763	0.870	0.841	329 ± 57
	FAENet	0.554	0.623	0.546	0.578	0.575	1263 ± 766
	Base disconnected-FAENet	0.884	0.937	0.878	0.870	0.892	1549 ± 1027
	Independent pooling-FAENet	0.698	0.813	0.697	0.766	0.743	1553 ± 1049
	Independent models-FAENet	0.690	0.782	0.689	0.760	0.731	803 ± 327 803 ± 327 803 ± 327
	Attention-FAENet	0.732	0.814	0.715	0.756	0.754	902 ± 565

37th Conference on Neural Information Processing Systems (NeurIPS 2023).

The adsorbate refers to the molecule involved in the electrochemical reaction that is accelerated through the introduction of a catalyst material, represented by a semi-infinite periodic substructure.

In other words, it quantifies the extent to which the catalyst reduces the energy required for a chemical reaction to occur.

a surface is a cut through the bulk operated using Miller indices

A Adslab generation

This section summarises the adslab creation process of the OC20 dataset. Refer to the original paper [START_REF] Chanussot | Open catalyst 2020 (oc20) dataset and community challenges[END_REF] for more details.

Adsorbate and catalyst surface selection. The first step is to select the adsorbate and the catalyst that will compose the adslab. For the adsorbate, it extremely simple: it is sampled randomly from a