CLOUD: A Scalable and Physics-Informed Foundation Model for Crystal Representation Learning

Domain: physics
Published: 2026-03-17
Editorial outcome: minor_revision

Abstract
Predicting crystal properties is essential for understanding structure-property relationships and accelerating material discovery. However, conventional approaches like experimental measurements or density functional theory calculations are resource-intensive, limiting their scalability. While machine learning offers a promising alternative by learning complex structure-property relationships from data, existing models often rely on labeled data, adopt representations insufficiently capturing essential structural characteristics, and lack integration of physics, limiting their generalizability and interpretability. Here, we introduce CLOUD (Crystal Language mOdel for Unified and Differentiable materials modeling), a transformer-based framework trained on a Symmetry-Consistent Ordered Parameter Encoding (SCOPE) that encodes crystal symmetry, Wyckoff positions, and composition in a compact, coordinate-free string representation. Pre-trained on over six million crystals, CLOUD is fine-tuned on downstream tasks and achieves competitive performance across diverse material properties, demonstrating strong scaling with respect to data and model size. Furthermore, as a proof-of-concept of differentiable materials modeling, CLOUD is applied to predict the phonon-related properties by integrating with the Debye model. This approach enforces thermodynamic consistency and enables temperature-dependent property prediction without requiring additional data. These results demonstrate CLOUD’s potential as a scalable and physics-informed foundation model for crystalline materials, unifying symmetry-consistent representations with physics-grounded learning for property prediction and materials discovery.

Introduction
Predicting the properties of crystalline materials is of scientific and technological significance, as it enables insight into structure-property relationships and facilitates the discovery of functional materials [1–5]. Accurate knowledge of material properties underpins the design of functional materials across a broad range of applications, including electronics [6], energy storage [7, 8], and catalysis [8, 9]. One illustrative example is the heat capacity at constant volume (C_v) which serves as a representative example of a thermodynamic property that plays a key role in vibrational analysis and phase stability [10–12]. As such, the ability to predict these properties from a material’s structure is essential for accelerating the development of new materials and deepening our understanding of structure–property relationships.

Conventionally, crystal properties are evaluated through experimental characterization or theoretical simulations. In the lab, properties such as crystal structure, conductivity, and stability are measured using techniques like X-ray diffraction, spectroscopy, or electrochemical analysis. On the theoretical side, Density Functional Theory (DFT) is widely used to compute formation energies, electronic structures, and elastic constants from first principles. However, these approaches are often time-consuming and resource-intensive, limiting their scalability for high-throughput screening and materials discovery [13, 14]. Consequently, materials development has been constrained to a limited portion of the vast chemical space of synthesizable crystalline materials due to the high costs associated with evaluating material properties within specific systems [15, 16]. To overcome these limitations, the materials science community has increasingly embraced data-driven methods for predicting material properties with high throughput, reducing reliance on costly simulations and resource-intensive experiments. Among these, machine learning surrogates have proven especially effective, as they can learn complex structure-property relationships from data and enable rapid, scalable property predictions. Leveraging these strengths, ML models have the potential to substantially broaden the scope of material performance evaluation in the vast chemical space of crystals [17]. Pioneering work has significantly improved the prediction accuracy for crystals [18–23]. Rational design of representations that map crystals to continuous vector space has enabled the development of more accurate and generalizable machine learning models for crystal property prediction. Early models focused on composition-based representations [20, 24, 25], which are simple and broadly accessible but suffer from limited expressiveness due to their inability to distinguish between structural polymorphs. Structure-based models overcome this limitation by leveraging atomic coordinates, often through graph neural networks, to capture geometric and spatial information [18, 19, 21, 26, 27]. However, these models struggle with long-range interactions [28] and require costly structural data. To bridge the gap, coordinate-free models have emerged, offering a balance between structural awareness and computational efficiency by encoding symmetry and connectivity without explicit coordinates [29–31]. Nevertheless, current coordinate-free representations still face challenges in efficiently encoding critical structural information for accurate property predictions, motivating the need for improved representations that balance structural richness and computational simplicity.

Despite the advances in material property prediction with machine learning, state-of-the-art (SOTA) models typically rely on labeled training data obtained through costly wet-lab experiments or DFT calculations. This reliance limits their applicability and scalability due to the scarcity, high cost, and variability of labeled materials datasets [32]. Moreover, these models often exhibit poor out-of-distribution (OOD) generalization, which is critical for real-world applications where crystals with novel compositions or structures are common [17]. Recently, foundation models (FMs) [33], which are trained on broad data at scale and can be adapted to a wide range of downstream tasks, have emerged as a promising solution. Most FMs adopt a transformer-based architecture [34], enabling them to be pre-trained at scale using vast amounts of unlabeled data [35–37], which are orders of magnitude more abundant than labeled data. The key motivation for pre-training on large-scale unlabeled datasets is to learn universal representations of crystals through self-supervised learning. Aside from language-model-style pre-training, alternative self-supervised strategies have been explored in the crystal domain, such as Crystal Denoising Self-Supervised Learning (CDSSL) [38], which pre-trains graph neural networks by perturbing atomic positions and learning to recover the original interatomic distances. These representations, whether learned from sequence-based or geometric-denoising objectives, can be effectively transferred to a wide range of downstream tasks, improving model performance even in low-data regimes. Furthermore, the implicit parallelism of attention allows for massively scalable transformer architectures, whose performance can improve predictably with the size of the model and training data, following neural scaling laws [39–41]. If such a scaling behavior holds, FMs are expected to achieve SOTA performance simply by increasing the data and model capacity. FMs leverage self-supervised pretraining to capture high-level structural and chemical patterns from unlabeled data [33, 42], which has shown promise in predicting crystal properties [43–45]. In the molecular domain, models like ChemBERTa [46] and MoLFormer [47] have demonstrated that performance improves with larger pre-training datasets across both regression and classification tasks. Systematic studies on neural scaling laws have been conducted for molecular property prediction [42, 48], yet such analyses remain largely unexplored for materials-focused foundation models.

Another key question of developing foundation models for crystalline materials lies in the incorporation of physics. Even though high predictive accuracy could be achieved with machine-learning models, the prediction consistency in physics is not necessarily guaranteed. For example, previous work [28, 49] that learns the mapping from crystal structures to heat capacity (C_v) in a pure data-driven manner fails to learn the temperature dependency of C_v, hence these models are only capable of predicting C_v values under exactly the same temperature as the training data, which limits their applications in real-world material science tasks. Gurunathan et al. [50] predicted the phonon density of states (pDOS) from which C_v could be calculated with the standard Bose-Einstein distribution formalism, accomplishing enhanced predictive accuracy for temperature-dependent phonon-based properties. However, the model itself is still predicing the target property in one-shot given the input structure without injecting physics prior. The machine learning model can be further integrated with physical laws by predicting the variables within these laws, which, in turn, determine the target property. This framework can be trained end-to-end, provided that the physical law is differentiable. Recently, DiffMix [51] is developed as a differentiable geometric deep learning framework for mixtures that learns physical coefficients in predefined mixture physics laws, showing promising performances than purely data-driven methods, especially in small datasets. DiffMix exhibits strong extrapolation to higher temperatures, enhancing predictive reliability beyond training conditions. By integrating physics-based modeling with geometric deep learning, it enables more accurate and generalizable predictions, supporting efficient material design and optimization. Recent advances in differentiable materials modeling [52–54] have unified physical simulation and learning by enabling gradient-based optimization through physical laws, such as thermodynamic potentials, finite element solvers, or lattice mechanics. These frameworks allow models to predict intermediate physical quantities and propagate gradients through differentiable physical pipelines, supporting data-efficient and generalizable material design.

Herein, we design Symmetry-Consistent Ordered Parameter Encoding (SCOPE) to represent crystal structures as sequences and consequently propose Crystal Language mOdel for Unified and Differentiable materials modeling (CLOUD), a transformer-based foundation model for accurate and generalizable crystal property predictions. The model consists of a BERT [35] encoder and a multi-layer perceptron (MLP) prediction head. SCOPE is a symmetry-consistent string representation for crystals that integrates symmetry, equivalent sites, and compositions in a coordinate-free encoding of the crystal structures, as shown in Figure 1b. Using this representation, CLOUD is pre-trained via masked language modeling (MLM) on ∼6.3M unique crystals collected from the OPTIMADE [55], representing one of the largest high-quality datasets of DFT-relaxed crystals to date, then fine-tuned on downstream datasets to learn task-specific structure-property relationships (shown in Figure 1c). CLOUD exhibits competitive model performance when evaluated for predicting various DFT-calculated material properties, indicating its potential as the surrogate model in the loop of material discovery. Additionally, we investigate the scaling performance of CLOUD, suggesting potential performance gains with further scaling of data and model size. We further integrate CLOUD in a differentiable physics framework for physics-consistent property predictions. To showcase this capability, we fuse the model with Debye model [56] (CLOUD-DEBYE, Figure 1d) and assess its ability to predict phonon internal energy and heat capacity, two properties that rely on long-range interactions of crystal structures. CLOUD-DEBYE exhibits improved predictive performance than CLOUD itself and other GNN or GNN-hybrid models, in the meantime its predictions for the heat capacity under different temperatures match the experimental data and ensure thermodynamic consistency. Such results suggest the potential of CLOUD to learn long-range interactions in crystal structures and observe physics consistency by integrating with physics-based models in materials applications. Overall, our model extends the machine learning models for crystal representation learning and property predictions by developing a string representation for effective encoding, building an accurate and scalable foundation model, and unifying physics models and data-driven techniques, thus possessing the potential to boost real-world material design and discovery.

Results
SCOPE crystal representation

Crystal representations in machine learning often fall into two broad categories: composition-based and structure-based. Composition-only representations, while easy to obtain, lack the specificity required to uniquely determine a structure due to the non-uniqueness of composition-structure relationships. On the other hand, structure-based graph representations typically rely on 3D atomic coordinates, which are expensive to compute and often unavailable. Moreover, both categories often overlook a fundamental feature of crystals: symmetry, which governs many physical properties such as electronic behavior, optical response, and mechanical stability [57]. To capture symmetry explicitly, many recent models incorporate equivariance constraints [58–60], although these introduce complexity into the architecture. Huang et al. [31] took a step forward by incorporating space group numbers and symbols into a sequence-based representation suitable for language models. Meanwhile, recent work has demonstrated the importance of explicitly incorporating crystallographic symmetry and Wyckoff information into generative and predictive models. For example, WyCryst [61] introduces a Wyckoff-based representation that preserves symmetry during generation and relaxation; WyFormer [62] encodes crystals as fused element–Wyckoff tokens to enable symmetry-conditioned autoregressive generation; SymmCD [63] decomposes crystals into asymmetric units and symmetry operations, ensuring symmetry-preserving diffusion even for rare space groups. These works underscore that symmetry-aware, Wyckoff-based encodings serve as effective and expressive representations for crystal structures, substantially improving the physical validity and diversity of generated crystals. However, limitations exist in existing symmetry-based represtations for crystals. For example, the symbolic identifiers used in Huang et al. [31]—such as P 2_1/c or #14 — are compact summaries rather than full encodings of symmetry operations. They do not convey the complete set of transformations that define the space group, such as mirror plane, rotation axis, or screw axes, nor do they capture the relationship between different space groups.

To address these limitations, we propose SCOPE, a coordinate-free string representation that integrates three essential components of a crystal structure: the space group, the Wyckoff positions, and the composition. By omitting explicit atomic coordinates and instead encoding symmetry constraints and atomic assignments, SCOPE captures the essential building blocks of a crystal structure in a compact and expressive form.

A crystal structure can be described mathematically by its lattice and basis. Following the notations in Jiao et al. [64], the lattice is a periodic grid in three-dimensional space defined by three linearly independent lattice vectors forming a matrix L = [l_1, l_2, l_3] \in R^{3\times3}. Atomic positions within a unit cell are specified by fractional coordinates F = [f_1, f_2, ..., f_N] \in [0,1)^{3\times N}, which relate to Cartesian coordinates by X = LF. A crystal is fully specified by its lattice L, the set of atomic positions X, and the chemical species occupying those positions. However, not all configurations are physically meaningful—most crystals exhibit symmetry, which constrains how atoms can be arranged.

The full symmetry of a crystal is captured by its space group, a finite set of symmetry operations g \in G that map the crystal onto itself. Each operation acts on the atomic coordinates X \in R^{3\times N} as:

g \cdot X := OX + t1^\top (1)

where O \in O(3) is an orthogonal matrix representing rotations, reflections, or rotoinversions, and t \in R^3 is a translation vector. The crystal structure M is said to be invariant under g if g \cdot M = M (here, “=” denotes equivariance under the group action).

Traditionally, space groups are denoted by Hermann–Mauguin symbols or numbers (1–230), which provide compact, human-readable summaries of the symmetry elements. However, these representations are not sufficiently expressive for machine learning—they do not enumerate the full set of symmetry operations or distinguish subtle relationships between different groups.

To overcome this, SCOPE uses generator strings [65], which provide a symbolic representation of the symmetry generators for each space group. These generators are a minimal set of operations from which the entire space group can be constructed, allowing for a more granular and interpretable encoding of symmetry.

Among all the possible g \in G in a space group, a generator is a fundamental symmetry operator from which all other operations of the space group can be derived through the combination of those basis operators. By selecting a minimal set of generators, we can fully describe the symmetry content of a space group. Therefore, each generator is represented with a symbolic notation to form the generator string for space groups. There are 14 generator matrices for all the symmetry operations besides translational symmetry (namely O in g), and these matrices are represented by letters ranging from a to n. Combining the generator matrices with the translation components (t in g) which are represented with 10 upper-case letters, we can obtain all the symmetry operators. With the generators, it is possible to compile all 230 space groups in a short ASCII file that is only 4104 bytes long. The generator strings form a compact representation of all the allowable basic symmetries of the crystal structure, hence making the generator strings fingerprints of symmetries. The generator strings are particularly well suited for language modeling, as it preserves the symmetry operations and provides a richer input format than traditional number-based space group identifiers. More information on the generator strings can be found in the Methods section and Section 10.3.6 in De Graef et al. [65].

Following the symmetry part, SCOPE encodes the atom positions under the constraints of symmetry. Atom positions are usually represented by atomic coordinates in the form of (x, y, z) in 3D space. However, using coordinate information is computationally intensive in terms of both data acquisition and model training. The acquisition of coordinate information necessitates DFT calculations for structure relaxation, while the model will grow extremely large as the crystal system grows larger. Alternatively, Wyckoff positions [29, 66] describe sets of symmetry-equivalent points in the unit cell given space group could be used for atomic position representation. Each Wyckoff position corresponds to a class of points whose site-symmetry group is conjugate under G. The site-symmetry group of a point f_i \in [0,1)^3 is defined as the stabilizer subgroup:

G_i = {g \in G | g \cdot f_i \sim f_i}, (2)

where g \cdot f_i = Rf_i + t, and \sim denotes equivalence under lattice translations. This group contains all symmetry operations that leave the point invariant up to translation.

Applying the full space group G to a representative point f_i generates a set of symmetry-equivalent positions:

W_i = {g \cdot f_i | g \in G}, (3)

whose cardinality is called the multiplicity of the Wyckoff position. Each Wyckoff position is denoted by a label such as 4d, where the number indicates the multiplicity |W_i|, and the letter encodes the relative site-symmetry: labels earlier in the alphabet correspond to higher-symmetry positions.

The set of available Wyckoff positions is determined by the space group, and each atom in a crystal can be assigned to a Wyckoff position according to its symmetry environment. In our SCOPE representation, we encode the occupied Wyckoff positions along with the corresponding chemical elements. This yields a compact and symmetry-consistent description of the crystal structure that is invariant under symmetry-preserving transformations and independent of explicit atomic coordinates.

We complete the representation by including the elements that make up the material along with their fractions in the material. As a result, the representation is informed of the material’s composition, ensuring that the model captures both the identity and proportion of elements.

Combining the three parts above, we obtain the symmetry-consistent string representation for crystal structures. Examples of representations can be found in Figure 2a. The representation forms a top-down description of the crystal structure – starting from symmetry which is global information of the structure, followed by the equivalent sites under the symmetry and the stoichiometry of the crystal – without explicitly encoding the structure itself, accomplishing efficient encoding of crystals. This provides a more granular and physically complete encoding than categorical site-symmetry tokens [62] or binary symmetry matrices [63]. Unlike WyCryst’s one-hot representation with degrees of freedom [61], SCOPE captures the complete set of symmetry generators in a compact symbolic form, enabling fine-grained distinction between closely related space groups and naturally aligning with MLM for large-scale transformer pretraining. With the SCOPE representation for crystal structures as the input, we tokenize the sequence so that each chemically meaningful unit is treated as a single token. Specifically, each symbol in the generator string is tokenized individually, each Wyckoff position symbol is treated as a single token, and both element symbols and fractional values are represented as distinct tokens. All the float numbers are rounded to 2 decimal places.

Transferable and scalable crystal representation learning with CLOUD

CLOUD adopts the BERT architecture [35] and follows the pretrain-finetune paradigm for training. The transformer model is based on the attention mechanism [34], enabling it to effectively models long-range interactions and contextual relationships within sequences. Specifically, it employs scaled dot-product attention, which operates through query, key, and value matrices. Further details on attention are provided in the Methods section. In our implementation, the transformer encoder consists of 12 hidden layers, each comprising 12 attention heads. The hyperparameters of CLOUD are initially set based on standard BERT configurations and subsequently tuned to optimize model performance. The hyperparameters that are experimented with are summarized in Table S1. The transformer encoder is pretrained via MLM [67–69] for effective representation learning from the unlabeled sequences, where 15% of the tokens in a sequence are randomly selected for potential replacement, and the pretraining objective is to reconstruct the original tokens by leveraging contextual information.

The pretrained model is subsequently fine-tuned in a supervised manner for crystal property prediction. Specifically, the final hidden states (feature vector) of the special token ‘[CLS]’, positioned at the beginning of the sequence, is passed through the prediction head to generate the prediction output.

We use the Open Databases Integration for Materials Design (OPTIMADE) [55] dataset for CLOUD pre-training. OPTIMADE comprehensively collects crystal data from various databases and provides an API for users with easy access to download the data. We download and de-duplicate the initial ∼13M CIF files from OPTIMADE, resulting in ∼6.3M data for pre-training. With the pre-trained model, we evaluate the predictive performance of CLOUD on MatBench datasets [70] which comprise diverse materials properties computed via DFT. We compare CLOUD to diverse baseline models in term of mean absolute error (MAE) normalized by mean absolute deviation (MAD), which is given by MAE/MAD = \sum |y_{i,true} - y_{i,pred}| / \sum |y_{i,true} - \bar{y}_{true}| so that the metric is invariant to scaling, enabling comparison across datasets of different units [21]. The baseline models are primarily divided into three groups: structure-based models which utilize coordinates of atoms as input: CGCNN [19] and ALIGNN [21], structure-agnostic models which do not require structural information: Roost [20], and coordinate-free models which rely on the structures to build the input while do not require the inclusion of atomic coordinates: Wrenformer [71], MatInFormer [31], and SLICES-BERT [30]. CLOUD achieves competitive results across nearly all eight regression tasks compared to structure-agnostic and coordinate-free models. As shown in Figure 2b, CLOUD outperforms previous coordinate-free models in 7 out of 8 tasks and overall performs the best in jdft2d and dielectric datasets Notably, the model attains strong performance even on data-limited tasks, suggesting that pre-training on large unlabeled datasets enhances generalization. Additionally, CLOUD surpasses MatInFormer [31], indicating that explicit symmetry encoding via generator strings contributes to improved predictive accuracy. Compared to SLICES [30] which encodes bond connectivity while neglecting overall symmetry, CLOUD reaches lower prediction error with much more compact representation (the mean sequence lengths for CLOUD and SLICES are 31 and 513 respectively for the same pre-train dataset), and consequently less computing budgets. SLICES-BERT only surpasses CLOUD on perovskites on which none of the structure-agnostic or coordinate-free models perform well while SLICES benefits from connectivity information which conveys more subtle differences between perovskite structures. We provide the MAE results for each model on the eight regression tasks in Table S3. Figure 2c presents the comparison between the results of CLOUD with/without pre-training. Model performance on downstream tasks is significantly enhanced if the model is pre-trained compared to training from scratch.

Aside from predicting the material properties from the relaxed crystal structures, we are also interested in assessing CLOUD’s capability for structure stability prediction without knowing the relaxed crystal structures from DFT calculations. Therefore, we fine-tune the pre-trained model on MatBench Discovery [71] which is designed to evaluate machine learning models for materials discovery, particularly for predicting the thermodynamic stability of materials. The dataset aims to reflect practical challenges in the discovery process by requiring predictions based on unrelaxed crystal structures, which avoid reliance on expensive DFT calculations. We fine-tune CLOUD on Materials Project [72] formation energy training data and test on WBM dataset [73] which contains ∼257K OOD crystal structures generated by systematically substituting elements in pre-existing structures from Materials Project. Model performance is evaluated in terms of the ROC-AUC score, which is calculated by varying the stability threshold, recording the false positive rate (FPR, the fraction of unstable structures being misclassified as stable) and true positive rate (TPR, the fraction of stable structures that are correctly identified), plotting the ROC (Receiver operating characteristic) curve, and calculating the area under curve (AUC) score. In contrast to directly predicting the formation energy for relaxed structures where CLOUD gives a larger prediction error than many structure-based GNNs, CLOUD manifests comparable classification performance despite not using atomic coordinates, given that the model achieves the AUC score of 0.81 which is close to that of ALIGNN [21] and higher than that of CGCNN [19] and MEGNet [74], as shown in Figure 2d. These results highlight CLOUD’s potential for accelerating materials discovery by effectively screening unrelaxed structures and directly predicting the stability of the relaxed ones, thereby minimizing reliance on costly DFT calculations. To further showcase the capability of CLOUD in material discovery, we conduct a prospective discovery-style evaluation. We use CLOUD’s predicted formation energies to derive E_hull, rank candidates by predicted stability, and quantify discovery efficiency using Precision@k, Enrichment Factor, and Cost-to-M-stable metrics. These metrics emulate a realistic screening workflow under finite resource constraints. As shown in Figure S2a–c, CLOUD substantially improves efficiency by reducing the number of DFT calculations required compared to random screening. In a case study on the Al–Fe system, CLOUD successfully identifies Al_2Fe as stable, despite its absence from the Materials Project convex hull, in agreement with results from the WBM dataset (Figure S2d). These results highlight the strong potential of CLOUD to accelerate the discovery of novel stable materials.

Noticeably, the crystal structures from MatBench or Matbench Discovery datasets are ‘regular’ crystals that predominantly feature well-ordered structures and lack the structural diversity, dimensional variation, or disorder commonly seen in more complex or ‘unconventional’ materials which are, in fact, ‘conventional’ in real-world under finite temperatures. Therefore, we fine-tune the model on crystal structures from UnconvBench [75], which includes crystals with defects, large unit cells, or low-dimensional structures. Different from the cases for MatBench in which the sophisticated structure-based GNN models like coGN [22], coNGN [22], and ALIGNN [21] show superior accuracy, CLOUD demonstrates comparable performance to the SOTA model CrysToGraph [75] and significantly outperforms the other structure-based models (Figure 2e). Such evidence reveals the promising potential of CLOUD to be applied for material discovery. In addition, the best model for each UnconvBench task is either CrysToGraph or CLOUD, both of which take into account global information about crystal structures. Distinct from CrysToGraph which learns the long-range interactions with a graph-wise transformer, CLOUD directly encodes symmetry information via generator strings and Wyckoff positions, providing an alternative for encoding global information for crystal structures. We provide the MAE results for each model on the six predictive tasks in Table S6.

To further understand how the SCOPE representation contributes to better modeling of crystal structures, we examine the attention scores between the [CLS] token and the other tokens. Attention scores indicate the correlation between tokens, thereby reflecting the model’s interpretation of the encoded features. We focus on the attention scores between [CLS] and the other remaining tokens in the sequence because the [CLS] embedding is fed to the prediction head for property prediction and high attention scores likely imply a significant contribution of the corresponding tokens to the material property prediction.

To illustrate, we compute attention scores across hidden layers and attention heads for crystal structures in the gvrh test set from MatBench. The attention scores between [CLS] and other tokens are aggregated across attention heads, from which we identify the top-k important tokens with the highest attention scores for each data point. We prioritize analyzing the attention scores from the first hidden layer as it is closest to the input and consequently learns more semantic information from the sequences [76]. To quantify the relationship between [CLS] and space group tokens, we introduce two probabilities: p_1, the probability that at least one space group token appears in the top-k list for a given sequence, and p_2, the proportion of space group tokens among the top-k tokens. We vary k and find that the probability of retrieving at least one space group token among the top-k tokens (p_1) increases from 0.93 at k = 1 to 0.99 at k = 3, while the average proportion of space group tokens among the top-k tokens (p_2) remains consistently high, ranging from 0.91 to 0.93. This confirms that space group tokens receive significant “attention” (contributions to the machine-learned feature) considering their average portion in a sequence in this dataset is ∼0.64. To assess the statistical significance of this prioritization, we perform a one-sample t-test comparing the observed values of p_1 and p_2 to a baseline of 0.5, which serves as the null hypothesis representing random or uniform attention. The test rejects the null hypothesis (p < 0.001), indicating that p_1 and p_2 are significantly greater than the expected 0.5. The 95% confidence intervals for p_1 ([0.9236, 0.9444] for k = 1, [0.9820, 0.9916] for k = 2, [0.9886, 0.9960] for k = 3) and p_2 ([0.9236, 0.9444], [0.9165, 0.9329], [0.9002, 0.9156]) consistently exceed 0.64, demonstrating that space group tokens are systematically prioritized. These results validate that the model effectively learns symmetry-consistent representations rather than distributing attention based solely on token frequency.

In addition to the evaluation of model prediction accuracy, it is crucial to understand how performance improves with increased model and data—an essential consideration for foundation models as they are often deployed at large scales to achieve generalization and transferability across tasks. Therefore, we further examine how the model performance scales with the number of pre-training data and model parameters. The model performance (cross-entropy loss for the pre-training task) empirically scales with the number of tokens D and non-embedding model parameters N (excluding the model parameters for vectorizing tokens) in a power law relationship, which is expressed by Hoffmann’s scaling law [40]:

L(N,D) = A/N^\alpha + B/D^\beta + E (4)

We vary the number of pre-training data and non-embedding parameters in the encoder of the model to fit the scaling law accordingly. We solve the optimization problem with the L-BFGS algorithm [77] to obtain the following result:

A=90.7464, B =21.9854, E =0.2296, \alpha=0.4339, \beta =0.3518

With the fitted scaling law, we estimate the optimal model size for a given compute budget. Following the derivation in literature [40], the optimal model size and data size can be written as:

N_{opt}(C) = G(C/6)^a, D_{opt}(C) = G^{-1}(C/6)^b (5)

where

G = (\alpha A / \beta B)^{1/(\alpha+\beta)}, a = \beta/(\alpha+\beta), b = \alpha/(\alpha+\beta) (6)

Given the \alpha and \beta we fitted, we can readily derive that a = 0.45 and b = 0.55 for CLOUD. Interestingly, the a and b from the empirical scaling law that we fit for CLOUD are close to the values in Hoffmann’s scaling law where a = 0.46 and b = 0.54 (derived from parametric modeling of the loss) [40]. Given a \approx b, the number of parameters and number of tokens should be scaled evenly, unlike the case for Kaplan’s scaling law which infers that the model size should scale faster than the number of tokens (a = 0.73, b = 0.27) [39].

The fitted scaling law is also visualized in Figure 3. The fitted entropy term E, which reflects the natural entropy of text data, is equal to 0.2296 for CLOUD, significantly lower than the E = 1.69 reported in [40], suggesting that the underlying data distribution for CLOUD is more structured or predictable, allowing for a smaller irreducible loss floor. The result aligns with our expectation, as the string representation that we design follows certain rules distilled from chemical knowledge.

Despite the strong scaling performance, the visualization in Figure 3b, however, also reveals the limitation of scaling. The entropy term E, which reflects the natural entropy of the text data, the benefits of scaling diminish because the loss asymptotically approaches E as the compute budget increases, highlighting that beyond a certain compute threshold, additional resources yield diminishing returns. This suggests that while scaling remains beneficial, achieving further gains may require innovations in architecture, optimization, or data efficiency rather than simply increasing model size and training tokens. Specifically, our symmetry-consistent string representation mainly focuses on efficient encoding while not ensuring completeness since this version does not encode free variables in Wyckoff positions and lattice parameters for unit cells. Relevant discussions are also provided in section S1 in Supplementary Information. Therefore, future work that includes richer information while maintaining efficiency is required for the design of foundation models for crystal structures.

CLOUD-DEBYE: Integration of CLOUD in a differentiable and physics-consistent materials thermodynamics framework

To advance beyond accurate but purely data-driven predictions, we propose integrating CLOUD with physical models to form a differentiable, physics-consistent framework. This enables the model not only to make predictions aligned with known physical laws but also to propagate gradients through those laws, allowing end-to-end training.

As a representative case study, we demonstrate this approach by combining CLOUD with the Debye model to predict temperature-dependent phonon properties such as phonon internal energy (U) and heat capacity (C_v). They arise from collective lattice vibrations that are governed by the global periodicity of the crystal and the vibrational modes spanning multiple unit cells. According to the Debye model [56], U and C_v can be expressed as:

U \approx 9Nk_B T (T/\Theta_D)^3 \int_0^{x_D} x^3 / (e^x-1) dx (7)

C_v \approx 9Nk_B (T/\Theta_D)^3 \int_0^{x_D} x^4 e^x / (e^x-1)^2 dx (8)

x_D = \Theta_D / T (9)

\Theta_D = \hbar \omega_D / k_B (10)

where \omega_D is the Debye frequency, N is the number of atoms, \hbar is the reduced Planck constant, k_B is the Boltzmann constant, and T is the temperature. The Debye frequency could be further written as \omega_D = (6\pi^2 N / V)^{1/3} v where V is the volume of the N atoms and v is the speed of sound in the material. The sound velocity v can be approximated using the first-order Hooke’s law as:

v \approx \sqrt{C/m} (11)

Here, C denotes the effective spring constant, m is the atomic mass within the primitive cell, and d is the effective interplanar spacing along the direction of vibration. Therefore, U and C_v depend on lattice volume, bonding strength, and sound velocity — all of which are inherently long-range characteristics of the crystal structure. As a result, accurate prediction of U and C_v requires models capable of capturing these long-range structural dependencies and periodic constraints. However, recent systematic analyzes have shown that convolution-based Graph Neural Networks (GNNs) face inherent limitations in this regard [28, 78]. Specifically, GNNs are prone to over-smoothing [79, 80] and over-squashing [81, 82], which restrict their receptive fields and impair their ability to learn long-range orders across the crystal. As a result, they often struggle with modeling properties that are sensitive to global structure. Transformer-based language models offer a promising alternative. Unlike GNNs, they capture interactions between all tokens without relying on sequential message passing, thereby enabling the modeling of long-range orders. Our proposed model, CLOUD, is built upon such a language model architecture and leverages the SCOPE representation, which directly encodes the space group—a global symmetry descriptor of the crystal. This makes CLOUD particularly well-suited for capturing the periodic and long-range features required for phonon property prediction. To evaluate this capability, we evaluate CLOUD on the internal energy (U) and heat capacity (C_v) datasets curated by Gong et al [28]. These datasets, constructed from Materials Project entries [72], include DFT-computed phonon properties derived from phonon density of states (DOS) at 300K. By comparing performance on these tasks, we aim to assess whether language models can surpass state-of-the-art GNNs in learning the global structural representations necessary for accurate phonon property prediction.

Meanwhile, the original benchmarking experiments in the literature [28] are conducted by predicting the U and C_v values at T = 300K directly from the models. Such an experimental setting comes with an obvious limitation that the model cannot learn the temperature dependency of C_v and U. In other words, learning the mapping between crystal structures and those properties is intrinsically ill-posed. Possible solutions include training different models for different temperatures or encoding temperature explicitly in the crystal embeddings, but both of them are computationally expensive as C_v labels under different temperatures are required, and they do not satisfy the thermodynamic consistency. According to Debye model, the phonon internal energy and heat capacity per unit cell is given by Equation 7 and 8, respectively [56]. Contrary to C_v and U, the Debye temperature \Theta_D is not dependent on the temperature but rather an intrinsic property of crystal structures. Gong et al. [28] claimed that the concatenation of descriptors improves the modeling of \Theta_D, whereas the thermodynamic consistency is still not achieved. Therefore, we propose CLOUD-DEBYE which is an integration of CLOUD and Debye model – CLOUD predicts \Theta_D which is fed to Debye model for predicting C_v and U. The differentiable implementation of Debye model enables training with C_v or U labels end-to-end without requiring the knowledge of Debye temperature.

The evaluation results for predicting U and C_v are presented in Table 1. The model performance is evaluated according to the metric MAE/MAD. Usually, a model will be considered as a good predictive model if MAE/MAD < 0.2 [21, 83]. While all the GNNs result in high MAE/MAD values, CLOUD, without the integration of Debye model, reaches 0.16 and 0.15 on C_v and U, respectively. CLOUD-DEBYE further reduces the MAE/MAD metric significantly, outperforming all the descriptor-hybridized GNNs proposed by Gong et al. [28] which aimed to compensate the lack of long-range information by concatenating global descriptors. The results manifest CLOUD’s capability of learning long-range interactions in crystal structures as well as the significance of introducing thermodynamic consistency into the material machine learning framework. Interestingly, CLOUD-DEBYE achieves improved prediction accuracy than the purely data-driven CLOUD model pre-trained on ∼6.3M data, even when trained from scratch. This underscores the importance of correctly incorporating physical principles: by enforcing thermodynamic consistency through the Debye model, CLOUD-DEBYE benefits from strong inductive bias that leads to accurate and physically valid predictions. These results highlight that proper integration of physics is as important as, if not more important than, advances in machine learning techniques for building scientific foundation models.

We further visualize the error in C_v prediction (C_{v,target} - C_{v,predict}) by different models in Figure 4b. Compared to CGCNN [19] and ALIGNN [21], two competitive GNNs for material property predictions, which exhibit wide distributions and high variability of \Delta C_v, the CLOUD and CLOUD-DEBYE error distributions exhibit much smaller spread. Noticeably, CLOUD-DEBYE results in a much narrower error distribution and substantially reduces the maximum absolute error by more than a factor of six, suggesting that integration with Debye model leads to more consistent C_v predictions. Two representative structures, Ga_3Os (mp-570844) and ZrO_2 (mp-1565), are highlighted to illustrate this effect: while CGCNN and ALIGNN exhibit errors exceeding 150–270 J/(mol-cell*K), CLOUD-DEBYE produces much more accurate predictions, reducing the error to 6.5% of the target value from the dataset, or even lower. Furthermore, CLOUD-only fails to lower the error and still results in unphysical predictions without the integration of the Debye model. Such a result indicates that even when using the same model architecture and training objective, physics integration constrains unphysical outliers and enforces more reliable predictions. Figure 4c-d further examine the error distribution with respect to structural characteristics. In Figure 4c, the test set data are grouped by primitive cell volume, revealing that CGCNN [19] and ALIGNN [21] struggle particularly with large-volume structures, often producing substantial or unphysical errors. This trend reflects the limitations of GNNs with finite receptive fields, which hinder their ability to capture long-range interactions essential for accurate C_v prediction in extended systems, consistent with the observation by Gong et al. [28]. In contrast, CLOUD significantly reduces the error across all volume bins, and CLOUD-DEBYE further improves upon this, achieving low and stable errors even for the largest structures. Figure 4d presents a heatmap of average error across crystal systems, where CLOUD-DEBYE consistently outperforms all GNNs and purely data-driven CLOUD. The error reduction is especially notable for complex systems such as orthorhombic, tetragonal, trigonal, and hexagonal crystals, underscoring the generalizability of physics-informed modeling across diverse symmetry classes. Therefore, the results highlight CLOUD’s superior accuracy and robustness for learning features relevant to C_v prediction, and also demonstrate the potential of combining it with physics-based models for real-world material applications.

After being trained on C_v data at T = 300K, CLOUD-DEBYE can be extrapolated to predict C_v at other temperatures. We evaluate the model performance on different chemical systems: Al_2O_3 and Li_2O are two oxides, CaTiO_3 belongs to perovskite, and GaN is a III-V compound which is broadly used as semiconductors. We vary the temperature from 0 K to near the melting points of the compounds. The predictions are compared to experimental results reported in previous literature [84–87] as well as phonon calculation results with Perdew-Burke-Ernzerhof (PBE) as the exchange-correlation functional from the Materials Data Repository (MDR) database [88]. As experiments are usually conducted under constant pressure, the experimental results are reported as constant-pressure heat capacity (C_p) which differs from C_v by C_p - C_v = \alpha^2 T V / \kappa_T where \alpha is the expansion coefficient and \kappa_T is the isothermal compressibility. However, C_p is usually a good approximation for C_v of solids under low temperature [56].

The prediction results by CLOUD-DEBYE are visualized in Figure 5. For all the four systems that are experimented with, the CLOUD-DEBYE’s prediction results match the experimental and DFT-calculated values well under different temperatures, even though the model is only trained with C_v under 300K. In addition, CaTiO_3 is not included in either of the training or test set, which reflects the superior generalization performance of our model. Furthermore, the integration of Debye model ensures the thermodynamic consistency of model predictions. According to Debye model, the heat capacity approaches

C_v \to 9Nk_B (T/\Theta_D)^3 \int_0^{x_D} x^2 dx = 3Nk_B (12)

as T \to \infty, which is expected (Dulong-Petit Law), since at high temperatures all phonon modes become activated. In the limit of T \to 0,

C_v \to 9Nk_B (T/\Theta_D)^3 \int_0^{\infty} x^4 e^x / (e^x-1)^2 dx \sim Nk_B (T/\Theta_D)^3 (13)

Both of the two constraints are satisfied by the model predictions.

The results above strongly suggest the advantage of combining deep learning model with physics-based model as the temperature dependence of C_v is directly embedded in the model, improving the data efficiency along the temperature dimension. Previous literature [11] has discussed the challenges of developing empirical [89, 90] or machine-learning models for temperature-dependent heat capacity, primarily due to the difficulty of obtaining thermochemical property data as a function of temperature. Gurunathan et al. [50] alleviated this problem with a DOS-mediated approach: the ALIGNN force field [21] is trained to predict phonon Density of States (pDOS) from which the phonon-based properties including C_v could be derived. In contrast, our model integrates Debye model into the learning framework so that it can be trained end-to-end with C_v data, learning the features to predict \Theta_D from the crystal structures while not requiring any DFT-calculated labels of pDOS or \Theta_D. We expect the model to serve as a differentiable framework that learns essential physical properties of crystal structures from data of high-quality and easy-accessibility and predicts material properties while observing physics consistency.

Looking ahead, the CLOUD–Debye framework exemplifies a broader approach: embedding differentiable physics models into machine learning for physics-consistent predictions. Beyond the heat capacity and phonon internal energy, the same principle could be applied to other crystal properties governed by well-established physics laws. Examples include learning the full phonon density of states to access multiple vibrational observables, deriving elastic moduli through Voigt–Reuss–Hill averaging of elastic tensors, and modeling thermal expansion via the Grüneisen parameter and bulk modulus. We view these as promising directions for extending CLOUD into a general framework for physics-informed property prediction, and we provide further discussion in the section S3.5 in the Supplementary Information.

Methods
Symmetry-Consistent Ordered Parameter Encoding

We design SCOPE, a symmetry-consistent string representation for crystal structures that integrates symmetry, equivalent sites, and constituting atoms to efficiently encode the structural and compositional information, eliminating the need for coordinate information or equivariant models. A SCOPE string for a crystal structure starts with the generator string which represents the basic symmetry operators of the corresponding space group. To build the generator string for each space group, all non-translational symmetry operations are represented by 14 generator matrices labeled a to n, while 10 uppercase letters encode translation components. Together, these generators define all symmetry operators and enable encoding all 230 space groups in a compact 4104-byte ASCII file [65].

Take the No. 35 space group as an example. Its generator string is given by ‘03aDDDbOOOjOOO0’. The ‘0’ at the beginning of the generator string indicates the inversion operator is not a generator. After that, ‘3’ indicates that there are three generator matrices. Every four letters correspond to a generator: the first lower-case letter determines which of the 14 matrices is to be used. and the subsequent 3 letters stand for the translation components in three dimensions. For instance, ‘aDDO’ corresponds to the symmetry operator:

g = (O | t) = ( (1 0 0 | 1) (0 1 0 | 2) (0 0 1 | 0) ) (14)

The last symbol in the generator string is ‘0’ which indicates that there is not any alternative choice of the origin for this space group.

We argue that generator strings are more informative than space group symbols, as space group symbols only provide a high-level classification of symmetry, whereas generator strings offer detailed insight into the specific symmetry operations that define the crystal’s structure. For example, Pmc2_1 (No. 36) and Amm2 (No. 38) only have one character in common in the space group symbols. However, their generator strings are ‘03aDDObOODjOOD0’ and ‘03aODDbOOOjOOO0’, respectively, suggesting that their symmetry generators are quite similar and only differ in the translation components. In fact, both space groups belong to the orthorhombic crystal system and exhibit mirror planes and two-fold rotation axes, but the subtle difference in the translation components (specifically in glide operations) reflects the underlying structural variation. This distinction is not apparent from the space group symbols alone, highlighting the value of generator strings in capturing more nuanced symmetry information. In contrast, Fmm2 (No. 42) has a face-centered lattice while Imm2 (No. 44) has a body-centered lattice, and the difference in lattice type significantly affects the symmetry operations and translations within the unit cell even though they both belong to Orthorhombic crystal system. Their space group names have 75% of the characters in common. Instead, their generator strings, ‘04aODDaDODbOOOjOOO0’ and ‘03aDDDbOOOjOOO0’, reflect more difference in the symmetry operators than their names.

Wyckoff positions in a SCOPE string are denoted by a symbol containing a number and a letter, e.g., 4d. The number denotes how many elements are on the corresponding orbit, also known as the “multiplicity” of the site. The letter is used to distinguish the different Wyckoff positions for a specific space group. The assignment rule for letters is the letters showing up earlier in the alphabet are given to sites with lower multiplicity.

Finally, the generator string for the space group, the Wyckoff symbols for atomic positions, and material compositions (e.g. NaCl, Al_2O_3, CaTiO_3, etc.) are integrated to form a coordinate-free representation of crystals while maintaining crucial structural information like symmetries, enabling efficient encoding of crystal structures so that large language models can be utilized for crystal representation learning. To construct the SCOPE representation, we parse CIF files using pymatgen [91] to extract structural and symmetry information. For each structure, we apply ‘SpacegroupAnalyzer’ to determine its space group number and Wyckoff positions, and consequently retrieve the generator strings corresponding to the space groups. Finally, the material composition is extracted from the structure’s fractional composition and appended to form the complete string. This results in a sequence of the form: [generator string] | [Wyckoff symbols] | [composition], efficiently capturing structural symmetry and stoichiometry without relying on atomic coordinates.

Datasets

Pre-train Dataset: for CLOUD pre-training, we collect ∼6.3M unique crystal structures from Open Databases Integration for Materials Design (OPTIMADE) [55] API. OPTIMADE integrates data from various major materials databases and provides an application programming interface (API) to make the data accessible to users. The contributing databases include: AFLOW [92], Alexandria [93], Computational Two-Dimensional Materials Database (C2DB) [94], Crystallography Open Database (COD) [95], Joint Automated Repository for Various Integrated Simulations (JARVIS) [96], Materials Cloud [97], Materials Platform for Data Science (MPDS) [98], Materials Project [72], Material-Property-Descriptor Database (MPDD) [99], Material Properties Open Database (MPOD) [100], NOMAD [101], Open Database of Xtals (odbx) [102], Open Materials Database (omdb) [103], and Open Quantum Materials Database (OQMD) [104]. Therefore, such comprehensive datasets enable effective pre-training on crystal representations. The original data we downloaded from OPTIMADE consists of ∼13M CIF files. To remove duplicate entries, we apply the deduplication rule from Antunes et al. [44]: structures are first grouped by identical chemical composition and space group number, and within each group we retain only the structure with the smallest volume per formula unit. This criterion follows the observation that structures with identical composition and space group but larger volumes are often redundant polymorphs or duplicates of less compact configurations. Applying this rule to our dataset removes ∼51.5% of the total structures, resulting in ∼6.3M structures for subsequent processing. Finally, we use the de-duplicated data as the pre-training set.

Fine-tune Datasets: we evaluate CLOUD across a range of materials datasets that test different aspects of model performance. These include:

• MatBench [70], a diverse suite of materials property prediction tasks based on DFT and experimental data. We focus on eight regression tasks that include structural information and are widely used to assess model accuracy across various physical and chemical properties.

• MatBench Discovery [71], which focuses on predicting the stability of novel materials using models trained on customized formation energy data. The evaluation set comprises approximately 257,000 out-of-distribution crystal structures generated through elemental substitution [73]. Summary results are provided in the main text, with detailed evaluation in the Supplementary Information.

• UnconvBench [75], a dataset of unconventional crystal systems, including low-dimensional materials, metal-organic frameworks (MOFs), and defected structures. It provides a challenging testbed for assessing model robustness beyond well-ordered bulk crystals. Key results are summarized in the main text, with full analysis in the Supplementary Information.

• Phonon thermodynamics dataset [28], which includes phonon internal energy (U) and constant-volume heat capacity (C_v) computed at 300K for Materials Project structures [72]. These properties are derived from phonon density of states (DOS) and reflect long-range structural periodicity and bonding strength. We use this dataset to examine CLOUD’s ability to capture global structural features relevant for thermodynamic property prediction.

We include the information about the fine-tuning datasets that are used in this work in Table 2. More information about the fine-tuning datasets can be found in section S2.1 in Supplementary Information.

Training Details

The CLOUD model is developed upon the transformer encoder architecture [34, 35]. Specifically, it employs Scaled Dot-Product Attention, which computes attention scores by measuring query-key alignment, scaling by the dimension of key, and applying softmax to derive weighted values:

Attention(Q,K,V) = Softmax(QK^T / \sqrt{d_k}) V (15)

Multi-head attention extends single attention by applying multiple linear projections to Q, K, and V, enabling parallel attention computations. The resulting outputs are concatenated and re-projected to form the final representation, allowing the model to capture diverse subspace information effectively.

CLOUD is first pre-trained via MLM on the dataset collected from OPTIMADE. The model is pre-trained in a self-supervised learning approach: without requiring the explicit supervision of labeled data, the model leverages intrinsic structure or patterns within the input data for learning. During pre-training, 15% of the tokens in a sequence are selected for potential modification, among which 80% are masked, 10% are replaced with randomly chosen vocabulary tokens, and the remaining 10% are left unchanged. The strategy, which is consistent with the original implementation proposed in literature [35], ensures the model learns meaningful contextual embeddings while biasing representations toward the actual observed tokens. The pre-training dataset is split into training and validation set by 80/20. AdamW optimizer is used with a learning rate of 5 \times 10^{-5} along with a cosine learning decay. The model is trained with batch size 2048 for a total of 50 epochs.

For the downstream task fine-tuning, we add a randomly initialized MLP after the encoder to output the property predictions. For MatBench, UnconvBench, and U and C_v datasets, five-fold cross-validation is performed and the model is evaluated according to the average score on the test set; for MatBench Discovery, the Materials Project data are split into the training and validation set by 95/5 and the model is evaluated according to the performance on the WBM test set. We perform hyperparameter tuning for each dataset separately for optimal model performance. We provide additional information on model implementation in section S2.2 in Supplementary Information.

To ensure a fair comparison on datasets from MatBench, we limit the prediction head for the three pre-trained language models (MatInFormer [31], SLICES-BERT where we train on SLICES [30] with the same architecture as CLOUD, and our CLOUD in this work) to be a single-layer MLP to ensure consistency with MatInFormer architecture implementation [31]. All the baseline results except for SLICES-BERT are obtained from the MatBench leaderboards [70]. When fine-tuned on MatBench Discovery for crystal structure stability classification, the model is first trained to predict the formation energy of the given structure. Based on the predictions for formation energies, the energy above the convex hull (E_hull), the distance to the convex hull spanned by competing phases in the same chemical system, is readily obtained to determine the stability of the materials. A material will be classified as stable if its E_hull is lower than a threshold. The stability threshold for E_hull, typically set to 0, could be dynamically adjusted in order to account for the systematic prediction shifts unique to each model. Therefore, we vary the threshold and take an record of the statistics such as true positive rate (TPR) and false positive rate (FPR) which are used for further analysis.

When fine-tuning the model on the phonon internal energy (U) and constant-volume heat capacity (C_v) datasets with the CLOUD-DEBYE framework, we integrate the pre-trained CLOUD with the Debye model. Specifically, CLOUD is trained to predict the Debye temperature (\Theta_D), a characteristic property of the crystal structure that governs the temperature dependence of U (T) and C_v (T) through the Debye model. This approach introduces an inductive bias grounded in solid-state physics, enabling the model to generalize more effectively across temperatures. Moreover, while experimental or high-fidelity computational labels for C_v are relatively accessible, direct labels for \Theta_D are more difficult to obtain. By leveraging the Debye model as an intermediate physical prior, CLOUD indirectly learns to predict \Theta_D from the structure and derives U(T) and C_v (T) accordingly. To enable end-to-end training with U or C_v as supervision targets, we adopt a differentiable implementation of the Debye model (Equations 7 and 8), where the integrals are evaluated using Gauss-Legendre quadrature in PyTorch to allow gradient backpropagation through the physical model.

Scaling Analysis

Transformer-based language models benefit from the implicit parallelism of the attention mechanism and exhibit predictable performance improvements as a function of scale. Empirical studies have shown that model performance, often measured via cross-entropy loss, follow power-law scaling laws with respect to both the number of model parameters and the amount of training data [39]. A canonical form of the scaling law introduced by Kaplan et al. [39] is:

L(N,D) = [N_c/N]^{\alpha_N} + [D_c/D]^{\alpha_D} (16)

where N is the number of tokens, D is the number of non-embedding model parameters, and N_c, D_c, \alpha_N, \alpha_D are dataset-dependent empirical constants. This formulation captures the trade-off between model capacity and training data volume under a fixed compute budget.

To better capture the behavior of loss in realistic settings and account for the irreducible entropy of the data distribution, Hoffmann et al. [40] propose a refined scaling law:

L(N,D) = A/N^\alpha + B/D^\beta + E (17)

where A, B, E, \alpha, and \beta are fitting parameters. The constant term E reflects the intrinsic entropy of the data distribution and defines the lower bound on achievable loss.

We evaluate the scaling performance of the proposed CLOUD model by systematically varying the number of pretraining tokens and the number of non-embedding parameters in the BERT encoder. We fit the empirical scaling behavior using the Hoffmann scaling law (Equation 17) to capture the interaction between model size and dataset scale.

The fitting process involves solving the following optimization problem:

min_{A,B,E,\alpha,\beta} \sum_i Huber_{\delta} (log \hat{L}(N_i, D_i) - log L_i) (18)

where \hat{L}(N_i, D_i) is the predicted loss from the scaling law and L_i is the observed loss for the (N_i, D_i) configuration. We use the Huber loss with threshold \delta = 10^{-3}, which provides robustness to outliers, as recommended in Hoffmann et al. [40]. The optimization is performed using the L-BFGS algorithm [77], with multiple random initializations to avoid suboptimal local minima.

This analysis allows us to estimate the compute-optimal scaling trajectory for CLOUD. By identifying the exponents \alpha and \beta, we can characterize how performance scales with additional compute and whether improvements are better obtained by increasing model size or dataset size. Similar to the findings of Hoffmann et al. [40], we observe that balanced scaling of model size and data yields optimal performance under fixed compute constraints, underscoring the importance of data collection alongside model architecture design.

Data Availability

The datasets used in this study are publicly available in previously published sources. The original crystal structure and property data can be found in the corresponding literature [28, 55, 70, 71, 75]. The processed datasets generated in this study and used for model training are available in the project GitHub repository at https://github.com/BattModels/CLOUD. No restricted-access or privacy-protected data were generated or used in this study.

Code Availability

All code developed and used in this study is publicly available in the project GitHub repository at https://github.com/BattModels/CLOUD, and the archived version used in this study is available at Zenodo (DOI: 10.5281/zenodo.17666827) [105]. The repository includes scripts for data processing, model pre-training, fine-tuning, and evaluation.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "In this work, Xu and coauthors reported a representation learning framework for crystal structures that encodes materials information from a symmetry aware text description of crystals without precise atomistic coordinates as input. I think this work is generally interesting, particularly because of its powerful ability to capture global information of crystals, and the competitive performance of the benchmark compared with other structure-partially-blind models. I would recommend publication of this work in Nature Communications, if my following comments and suggestions could be properly addressed.\n1. I think the main merit of this kind of model that do not require precise coordinates is to facilitate materials discovery in terms of thermodynamics. In other words, people can use space group and composition to propose new structures, and employ this type of model to evaluate stability without explicitly calculating coordinates. In this sense, the authors should find a way to demonstrate the ability of this model to evaluate convex hull of structures, maybe doing some test in the materials discovery benchmark.\n2. It would be great for the authors to analyze more about the intrinsic limitation of this model. One example I can think of is the inability to distinguish AB type compound that permute the occupancy of sites for the two elements.\n3. It would also be great to analyze the robustness of this model to small perturbation of crystal structures, because such perturbation would change the symmetry classification of the crystal, which might hugely alter the text description. But such perturbation can exist in everyday research, from sources like set up of geometry optimizer, inaccuracy of MLFF, and even the tolenrence of space group analyzer.\n4. In addition to the debye model, what are some other properties of crystal that can modeled in this thermodynamic consistent and differentiable way? The authors should include some discussion about it.", "reviewer_2": "This work introduces a transformer-based framework for crystal representation learning. The representation is SCOPE, a coordinate-free crystal string that concatenates space‑group generator symbols (rather than just space group numbers), occupied Wyckoff positions, and composition information. The model is pretrained on OPTIMADE structures and fine‑tuned on a wide suite of downstream materials property benchmarks (MatBench, MatBench‑Discovery, UnconvBench). The reported results show a broadly competitive performance compared to graph neural networks (GNNs) and better than composition-based models. Attention analysis suggest the model attends strongly to symmetry tokens, supporting interpretability claims. The authors also investigate the scaling‑law fits related to the data/model size to performance, offering guidance for future scaling. The authors further present a CLOUD‑DEBYE that predicts a Debye temperature and propagates it through a differentiable Debye model to obtain heat capacity and internal energy curves, improving accuracy while enforcing thermodynamic consistency. Overall, the work is comprehensive and timely, covering symmetry‑aware representations, foundation models, and physics‑informed ML for materials. However, the work can be better organized to show better novelty and clarity. I suggest that authors have a revision before the work can be accepted.\n1. One key innovation in this study is the replacement of space group numbers with string generator with symmetry operations. The concept seems appealing and makes sense intuitively, but an ablation study is missing. It is unclear how much gain comes from (a) generator‑level symmetry tokens vs space group numbers, (b) inclusion of Wyckoff positions. Comparisons to MatInFormer are confounded by dataset size, tokenizer scope, and model capacity differences.\n2. “Long‑range order/interactions” claim is underspecified and under‑validated. The manuscript attributes CLOUD’s advantage on UnconvBench and phonon thermodynamics to capturing long‑range effects but provides no operational definition or ablation test. Prior work shows that mainstream GNNs often fail to capture periodicity and long‑range information. The statement will be much stronger if the impact of the string generator of space group numbers on disordered materials. It will be nice if the authors can conduct a controlled ablation study on (a) SG‑number only, (b) SG Hermann–Mauguin tokens, (c) Full SCOPE w/o Wyckoff, (d) Symmetry‑only (mask composition), (e) Composition‑only for perfect and disordered crystals datasets.\nAside from the above two main points, there are a few minor points to be addressed as well\n3. Distinguish “long‑range order” (periodicity) from “long‑range interactions” (physical couplings such as elastic or electrostatic); the manuscript alternates between them.\n4. Clarify de‑duplication rule (“smallest volume per formula + SG”) and report how many items dropped.\n5. Some relevant work are not cited. A few generative symmetry-aware representations e.g. WyCryst, WyFormer, and SymmCD—demonstrate the importance of explicit symmetry and Wyckoff encodings for both generation and property prediction.", "reviewer_3": "The authors present a framework for predicting the properties of crystalline materials and compare it with multiple deep learning models, which is a hot topic. However, I have some concerns about the novelty of the framework and the performance comparison.The lack of application of some models in specific material systems also means that there is a lack of proof of the effectiveness of the models in practical applications.\n1. The authors perform deduplication on the OPTIMADE dataset, but more details are not shown.\n2. The authors pre-trained on the OPTIMADE dataset, evaluated the performance on the MatBench dataset, and compared it with other deep learning models (CGCNN, ALIGNN). Have CGCNN and ALIGNN also been pre-trained? If not, the performance comparison is unfair.\n3. Large-scale pre-training is very expensive. Have the authors considered the issue of computational time? training time, testing time, model parameters\n4. PMCGNN's performance on the JARVIS-DFT dataset is already higher than ALIGNN, Matformer, etc. The authors should also compare. https://pubs.acs.org/doi/10.1021/acs.jcim.4c01200\n5. The author also mentioned that the performance of the current model depends on a large amount of labeled data, but the author's pre-training-fine-tuning model also depends on downstream labeled data. The self-supervised model may better handle this problem: CDSSL, https://doi.org/10.48550/arXiv.2408.17255."}, "author_response": "Response to Reviews on Article “CLOUD: A Scalable and Physics-Informed Foundation Model for Crystal Representation Learning”\nChangwen Xu1, ShangZhu1, VenkatasubramanianViswanathan1,2*\n1Department of Mechanical Engineering, University of Michigan.\n2Department of Aerospace Engineering, University of Michigan.\n*Corresponding author(s). E-mail(s): venkvis@umich.edu;\nContributing authors: changwex@umich.edu; shangzhu@umich.edu;\nManuscript ID: NCOMMS-25-45609\nWe are very grateful for the reviewers’ comments, suggestions, and valuable input on our manuscript. There is no doubt that these comments are valuable and very helpful for improving our manuscript substantially. In what follows, we would like to answer the questions you mentioned and give a detailed account of the changes made to the original manuscript. We strongly believe that the revised manuscript should meet the standards for publication in Nature Communications.\nComment:In this work, Xu and coauthors reported a representation learning frameworkforcrystalstructuresthatencodesmaterialsinformationfromasymmetryaware text description of crystals without precise atomistic coordinates as input. I think this work is generally interesting, particularly because of its powerful ability to capture globalinformationofcrystals,andthecompetitiveperformanceofthebenchmarkcom- pared with other structure-partially-blind models. I would recommend publication of thisworkinNatureCommunications,ifmyfollowingcommentsandsuggestionscould be properly addressed.\nAuthor Reply: We sincerely thank the reviewer for the thoughtful and encouraging assessment of our work. We deeply appreciate the recognition of the novelty of our symmetry-aware representation, its ability to capture global crystal information, and its competitive benchmark performance. In the following, we will address each of the detailed comments and suggestions provided.\nComment 1: I think the main merit of this kind of model that do not require precise coordinates is to facilitate materials discovery in terms of thermodynamics. In other words, people can use space group and composition to propose new structures, and employthistypeofmodeltoevaluatestabilitywithoutexplicitlycalculatingcoordinates. In this sense, the authors should find a way to demonstrate the ability of this model to evaluate convex hull of structures, maybe doing some test in the materials discovery benchmark.\nAuthor Reply: We thank the reviewer for the comment. In terms of the materials discovery benchmark, we have included the benchmarking results on MatBench Discovery, which could be found in Figure 2b in the manuscript and Figure S1, Table S4, and Table S5 in the supplementary information. In summary, MatBench Discovery is designed to evaluate ML models for materials discovery, particularly for predicting thethermodynamicstabilityofmaterials.MLmodelsaretrainedtopredictformation energies of the material structure then compare them to the convex hull to determine the stability. Models are evaluated on WBM dataset which contains unrelaxed structures only and is out-of-distribution from the training set. Our model manifests competitive performance, with ROC-AUC score comparable to ALIGNN [1]. When a 2 threshold of 0 is used, our model shows comparable performance to Wrenformer [2], while we observe that CLOUD benefits from more negative thresholds which lead to higher precision and DAF metrics. Therefore, the empirical evidence suggests the potential of CLOUD for materials discovery only based on the structure template (space group and Wyckoff positions) and the composition of the material, reducing the needs for high-fidelity simulation for stability identification.\nIn addition to that, we have further examined the model’s ability for novel stable material discovery by directly simulating a prospective screening workflow on the WBM dataset. Specifically, we predict the energy above convex hull E from the hull model’sformationenergypredictions,enablingtherankingofcandidatestructuresby their predicted thermodynamic stability. Using this ranking, we carry out material stability evaluations under practical budget constraints. We assess (i) Precision@k, whichquantifiesthefractionoftruestablematerialsamongthetop-rankedcandidates, (ii) the Enrichment Factor, which measures how strongly the model concentrates true stables in the top fraction of its ranking compared to random selection, and (iii) Cost-to-M-stable, which reports the number of candidates that must be evaluated before discovering M stable materials. These complementary metrics simulate how anexperimentalistorhigh-throughputworkflowwouldemploythemodeltoprioritize follow-up calculations or experimental evaluations. We additionally pick one example of novel material that manifests thermodynamic stability from the WBM test data, comparing the model prediction with the existing convex hull to demonstrate the model’s capability in material discovery.\nTodemonstrateourmodel’scapabilityinmaterialdiscovery,wehavebenchmarked on MatBench Discovery. The related analysis is included in the section Transferrable and scalable crystal representation learning with CLOUD in Results:\nAside from predicting the material properties from the relaxed crystal structures, we are also interested in assessing CLOUD’s capability for structure stability prediction without knowingtherelaxedcrystalstructuresfromDFTcalculations.Therefore,wefine-tunethe pre-trainedmodelonMatBenchDiscovery[2]whichisdesignedtoevaluatemachinelearn- ingmodelsformaterialsdiscovery,particularlyforpredictingthethermodynamicstability of materials. The dataset aims to reflect practical challenges in the discovery process by requiringpredictionsbasedonunrelaxedcrystalstructures,whichavoidrelianceonexpen- sive DFT calculations. We fine-tune CLOUD on Materials Project [3] formation energy trainingdataandtestonWBMdataset[4]whichcontains∼257KOODcrystalstructures generated by systematically substituting elements in pre-existing structures from Mate- rials Project. Model performance is evaluated in terms of the ROC-AUC score, which is 3 calculated by varying the stability threshold, recording the false positive rate (FPR, the fraction of unstable structures being misclassified as stable) and true positive rate (TPR, thefractionofstablestructuresthatarecorrectlyidentified),plottingtheROC(Receiver operatingcharacteristic)curve,andcalculatingtheareaundercurve(AUC)score.Incon- trasttodirectlypredictingtheformationenergyforrelaxedstructureswhereCLOUDgives a larger prediction error than many structure-based GNNs, CLOUD manifests compara- ble classification performance despite not using atomic coordinates, given that the model achievestheAUCscoreof0.81whichisclosetothatofALIGNN[1]andhigherthanthat of CGCNN [5] and MEGNet [6], as shown in Figure 1. These results highlight CLOUD’s potentialforacceleratingmaterialsdiscoverybyeffectivelyscreeningunrelaxedstructures and directly predicting the stability of the relaxed ones, thereby minimizing reliance on costly DFT calculations.\nFig. 1: Receiver operating characteristic (ROC) curves for models evaluated on Mat- Bench Discovery for structure stability classification. False positive rate (FPR) and true positive rate (TPR) are the fraction of nonstable/stable structures classified as stable. Area under curve (AUC) scores are calculated and presented in descending order.\nMorerelateddiscussiononMatBenchDiscoveryresultswasincludedinthesection S3.2 MatBench Discovery in the Supplementary Information:\nWe train the model on Materials Project [3] formation energy data from the v2022.10.28 MP release, then make predictions for the unrelated structures in the WBM dataset [4] which are generated via elemental substitution of MP source structures so that the gen- erated structures are not included in the training set. When evaluating the models on 4 MatBench Discovery, we record the classification results under varied thresholds and plot the receiver operating characteristic (ROC) curve for CLOUD and the major structure- based models that are trained on MatBench. We plot the ROC curves with more models included in Figure S1 compared to Figure 2c. Though not comparable to machine- learnedinteratomicpotentials(MLIPs)atthisstage,CLOUDstillshowssolidperformance compared with CGCNN, MEGNet, Wrenformer, etc.\nWe further examine the other metrics listed in Table S4. The regression metrics of CLOUD are again second to ALIGNN while outperforming CGCNN and MEGNet. The classification metrics for CLOUD are close to those of Wrenformer under the stability threshold of 0, as shown in Table S4.\nWelisttheclassificationresultsbyCLOUDunderdifferentstabilitythresholdsinTable S5. Note that the true labels for the test data are derived with the threshold of 0, con- sistent with the benchmark setting [2], while the dynamic threshold applies to the model prediction.MorenegativethresholdswillresultinhigherprecisionforCLOUDandsubse- quentlyhigherDAF,whichisalsoobservedformodelsthataremoreoptimisticinstability predictions like CHGNet [2].. However, the trade-off across metrics leads to decreased F1 score and TPR when a negative threshold is used.\nTo further demonstrate our model’s potential in material discovery workflows, we have added the following comment in the section Transferrable and scalable representationlearningwithCLOUDinResultsatpage10intherevisedmanuscript:\nTofurthershowcasethecapabilityofCLOUDinmaterialdiscovery,weconductaprospec- tive discovery-style evaluation. We use CLOUD’s predicted formation energies to derive E , rank candidates by predicted stability, and quantify discovery efficiency using Pre- hull cision@k, Enrichment Factor, and Cost-to-M-stable metrics. These metrics emulate a realistic screening workflow under finite resource constraints. As shown in Figure S2a–c, CLOUD substantially improves efficiency by reducing the number of DFT calculations requiredcomparedtorandomscreening.InacasestudyontheAl–Fesystem,CLOUDsuc- cessfully identifies Al Fe as stable, despite its absence from the Materials Project convex 2 hull,inagreementwithresultsfromtheWBMdataset(FigureS2d).Theseresultshighlight the strong potential of CLOUD to accelerate the discovery of novel stable materials.\nMeanwhile, We have added the following results and discussion in the S3.2 MatBench Discovery in the Supplementary Information:\nInordertoshowcasethepotentialofCLOUDinmaterialdiscoverytasks,wefurthercarry out a prospective discovery-style evaluation to directly demonstrate CLOUD’s ability to screen for new stable materials. On the WBM dataset which contains material structures that are out-of-distribution (OOD) to the training data, we calculate the energy above hull (E ), so far the same as the evaluation for MatBench Discovery. We rank all test hull candidatesbytheirpredicteddistancetohull,Epred,inascendingorder.Fixthestability hull 5 threshold τ (e.g., 0, 0.02, or 0.05 eV/atom) and define binary labels (cid:104) (cid:105) s = 1 Etrue ≤τ ∈{0,1}, (1) i hull,i for i = 1,...,N. Let π denote the permutation that sorts candidates by Epred (best hull first),soπ(1)isrank-1,etc.Wethenquantifydiscoveryefficiencyunderrealisticscreening budgets using three complementary metrics:\n• Precision@k: the fraction of true stables in the top-k suggestions,\nk\n1 (cid:88)\nP@k = s . (2)\nk π(j)\nj=1\n• Enrichment Factor at top α% (with k = ⌊αN⌋): the fold improvement over random screening,\nk N\n1 (cid:88) 1 (cid:88)\nEF = s / s (3)\nα k π(j) N i\nj=1 i=1\n• Cost-to-M-stable:thenumberofevaluationsneededtouncoverM truestables when following the ranked list,\n \n (cid:88) m \nCost-to-M = min m : s ≥ M . (4)\nπ(j)\n \nj=1\nFor comparison to the theoretical optimum, we also report the cost factor\nCost-to-M/M ≥1, where 1 corresponds to the ideal case.\nThe workflow resembles the true material discovery process where the most promising candidates,ratherthanthewholepool,aresortedoutforfinalevaluationwiththeoretical calculation or experiments. As shown in Figure 2a, the precision for stability screening remains high even up to the top-1000 candidates and reaches nearly 100% for k ≤ 200, showing that the model consistently prioritizes true stable materials at the top of the rankedlist.Theenrichmentfactorfurtherdemonstratesthatthemodelconcentratesstable compoundseffectively,withasmuchas∼10-foldenrichmentoverrandomselectionwithin the top fraction of candidates selected by the model (Figure 2b). The cost factor curves reveal that CLOUD identifies a target number of stable materials with only about 1.1– 1.3 times the evaluations required by the theoretical optimal (Figure 2c), highlighting its efficiency under practical screening budgets.\n1.0\n0.8\n0.6\n0.4\n0.2\n0.0\n-0.05 -0.02 0.00 0.02 0.05\nStability threshold\nnoisicerP\na Precision@k\n5%\n2%\nP@50\nP@100\nP@200\nP@500 1%\nP@1000\n-0.05 -0.02 0.00 0.02 0.05\nStability threshold\ndetaulave\n%\n-poT\nb Enrichment Factor\n8.0 5.2 3.5 2.3 1.8\n9.7 5.9 3.9 2.5 1.9\n10.9 6.3 4.1 2.6 1.9\n1.5\n1.4\n1.3\n1.2\n1.1\n1.0\n0.9\n0 200 400 600 800 1000\nTarget # of stable materials (N)\n)N\n/ slave(\nrotcaf\ntsoC\nc Cost-to-M-stable\n0.1\n=-0.05\n=-0.02 0.0\n=0.00\n=0.02 0.1\n=0.05\n0.2 Ideal (×1)\n0.3\n0.4\n0.0 0.2 0.4 0.6 0.8 1.0\nFraction Fe\n)mota/Ve(\nygrene\nnoitamroF\n10\n8\n6\n4\n2\nd Al-Fe convex hull + discovery\nAl Stable Al2Fe (true) Fe\nUnstable Al2Fe (pred)\nAl2Fe Ehull true = -0.026 eV\nEhull pred = -0.127 eV\n)modnar\n×(\nFE\nFig. 2: Evaluation of CLOUD in materials discovery tasks on the WBM dataset.\nWe vary the stability threshold τ ∈ {−0.05,−0.02,0.0,0.02,0.05} and calculate the following metrics: a Precision@k: fraction of true stables recovered within the top- k predictions across stability thresholds; b Enrichment Factor: fold enrichment over randomscreeninginthetopα%oftherankedlist;cCost-to-M-stable:numberofeval- uations required to uncover M stable materials, normalized by the ideal case (dashed line). d Case study: Al–Fe convex hull. CLOUD correctly identifies the intermetallic Al Fe as a new stable candidate.\n2\nAuthor Reply: We thank the reviewer for this thoughtful comment. We agree that discussing the intrinsic limitations of our model is important for a comprehensive understanding of its scope and potential extensions.\nOne limitation arises from the ordering of atoms in the input structure. Although ourSCOPErepresentationistheoreticallyinvarianttosymmetryoperations,thespe- cific assignment of elements to Wyckoff positions is determined by the atom order provided in the structure file. Since language models are inherently sensitive to token orderduetopositionalembeddings,thecurrentSCOPEencodingisnotstrictlyinvari- ant to permutations of the atomic list. To address this limitation, a natural remedy is to impose a canonical ordering of Wyckoff positions (e.g., sorting by site symme- try, starting from 1a, and permuting the elements accordingly), thereby removing dependence on input atom order.\nA second limitation concerns partial occupancies. In SCOPE, each Wyckoff posi- tionisexplicitlyassociatedwithitsoccupyingelement,whichensuresthatcompounds with different element–Wyckoff assignments (e.g., AB-type structures occupying dis- tinct Wyckoff sites) are distinguishable. However, when multiple elements partially occupythesameWyckoffposition,SCOPEdoesnotencodewhichspecificsymmetry- equivalentsitesareoccupiedbywhichelement.Instead,onlythecombinedoccupancy is represented. While such cases are relatively rare in the benchmark datasets used in thiswork,thislimitationbecomesmorerelevantinmodelingdisorderandsite-mixing phenomena.WeviewtheextensionofSCOPEtoexplicitlyhandlepartialoccupancies as an interesting direction for future work.\nFinally, we also provide a broader discussion of limitations and future work in SectionS1oftheSupplementaryInformation.Theseincludethefinitesizeofavailable pretraining data and the absence of certain structural details in SCOPE (such as free positional parameters within Wyckoff sites).\nWe have included the discussion on the limitation of shortage of pre-training data and missing some structural information in the section S1 Limitations and Future Work in the Supplementary Information:\nWhile CLOUD shows competitive performance in crystal representation learning and property prediction, it does encounter limitations and potential improvements for future explorations.Onelimitationisthatourmodelrequiresalargeamountofpre-trainingdata 8 to achieve SOTA performance in many of the downstream tasks according to the scaling law we fitted, whereas such a dataset has not yet been built. Currently, the largest crys- tal database we can build is the one collected from OPTIMADE which includes ∼6.3M unique crystal structures. Compared to molecular science in which ∼50 billion unlabeled SMILESofmoleculescanbeeasilyobtained[7–9],obtainingsuchalargenumberofcrys- talstructuresischallengingdespitethelargedesignspaceofcrystals.Notallcombinations ofatomscanformvalidcrystalstructures[10]asmanypotentialconfigurationsareinher- entlyunstableorchemicallyincompatible.Recently,alotofmachinelearningmodelshave beenproposedforcrystalgeneration[11–14],however,noneofthemhavepushedthedata size to the order of billions. Recently, OMat24 dataset [15] was released which contains over 100 million DFT calculations for solid-state materials, which is a huge contribution to the research community to advance materials science development with AI. In future work,weplantocomprehensivelyintegrateexistingdatasetslikeOMat24anddiversegen- eration approaches for crystal structure generation, ranging from rule-based methods to data-drivenapproaches.Wewillexaminehowtheexpanded,syntheticdatasetimpactsthe model performance in both pre-training scaling and downstream predictive accuracy.\nAnotherlimitationliesintherepresentationitself:therepresentationwedesignisbuilt on symmetry, whereas symmetry does not capture everything. For instance, while the spacegroup,Wyckoffpositions,andstoichiometryremainunchangedduringrelaxation,the atomic positions and cell parameters undergo significant variations. As a result, CLOUD is unable to leverage the MPtrj dataset, which comprises over 1.5 million crystal struc- turesfromrelaxationtrajectories,unlikemanyothermodelsontheleaderboard[2,16].In addition, the representation does not work for amorphous materials which lack the long- range order and symmetry that crystalline materials exhibit, thus limiting the versatility ofCLOUDinhandlingabroaderrangeofmaterials.Furthermore,someWyckoffpositions include free variables, meaning the representation encodes an ensemble of materials. As a result,itmaynotdistinguishbetweenmaterialswithdifferentphysico-chemicalproperties thatsharethesameprototypebuthavedifferentatomicpositions.Infuturework,weplan tosystematicallyinvestigatethedesignofstringrepresentationsforcrystalswhichenables effective encoding of more information besides symmetry and composition, while build- ing the multi-modal framework for crystalline materials as an alternative to incorporate comprehensive structural information.\nPlus, we have also mentioned the limitations and future work concisely in the Discussion section in the manuscript:\nSeveralpromisingdirectionsremainforfutureexploration.DespiteCLOUD’sstrongscal- ing performance with respect to data size, the availability of high-quality crystal data for pretraining remains limited. Therefore, constructing a large-scale, high-quality crystal structure database would be instrumental in advancing scientific foundation models for crystals.Moreover,thecurrentdesignofSCOPEcanbefurtherrefined.Currently,SCOPE 9 represents an ensemble of crystal structures that share the same space group, Wyckoff positions, and composition but differ in atomic coordinates due to the presence of free parameters in many Wyckoff positions.\nInadditiontothelimitationswehavealreadydiscussedintheoriginalmanuscript andtheSupplementaryInformation,wehaveaddedcontentrelatedtothelimitations we havediscussed intheresponseinto thesecond-to-lastparagraphof Discussion in the manuscript:\nAdditionally, we note two current limitations of the SCOPE representation. First, since Wyckoff labels depend on the input atom order, the representation is not strictly permutation-invariant, though this can be addressed by adopting a canonical ordering of Wyckoffsites. Second,SCOPEcannot fullycapture partialoccupancies orsite mixing,as it encodes only the combined occupancy rather than element-specific assignments. Both aspects highlight promising directions for future work. We have provided a more detailed andin-depthdiscussiononlimitationsandfutureworkinS1LimitationsandFutureWork in the Supplementary Information.\nMeanwhile, we have also added a detailed discussion into the second paragraph in the S1 Limitations and Future Work section in the Supplementary Information:\nInaddition,twofurtherlimitationsdeservemention.First,althoughtheSCOPErepresen- tationisinvarianttosymmetryoperationsbydesign,itisnotinvarianttopermutationsof the atom list in the input structure. Since Wyckoff labels are assigned in the order atoms areprovided,permutingtheatomsinaCIFcanleadtodifferentsequences,whichisunde- sirablegiventhatlanguagemodelsaresensitivetotokenorder,andinfutureworkweplan tosolveitbyenforcingacanonicalorderingofWyckoffsites(e.g.,sortedbysitesymmetry, startingfrom1a)andpermutingelementsaccordingly.Second,whileSCOPEdistinguishes compoundswithdifferentelement–Wyckoffassignments,itcannotfullyresolvecaseswhere multipleelementspartiallyoccupythesameWyckoffposition.Insuchsituations,onlythe combinedoccupancyisencoded,ratherthanthespecificsiteassignmentsofeachelement. This limitation becomes important for modeling disorder and site mixing, and extending SCOPE in this direction represents a promising direction for future work.\nComment 3: It would also be great to analyze the robustness of this model to small perturbation of crystal structures, because such perturbation would change the symmetry classification of the crystal, which might hugely alter the text description. But such perturbation can exist in everyday research, from sources like set up of geometry optimizer, inaccuracy of MLFF, and even the tolenrence of space group analyzer.\nAuthor Reply: We thank the reviewer for this valuable suggestion. We agree that robustness to small perturbations—such as those arising from geometry opti- mization settings, inaccuracies in ML force fields, or the thresholds of symmetry analyzers—is an important practical consideration for any structure-based representation.\nInterestingly, a recent study by Siron et al. [17] — which cites our workshop paper[18]—hasalreadyperformedasystematicrobustnessanalysisthatincludesthe SCOPErepresentationasoneofthebenchmarkedmethods.Intheirevaluation,multi- plestructureidentificationandfingerprintingapproaches(includingSCOPE,SLICES, Pymatgen’s StructureMatcher, and their proposed BAWL hashing) were tested for matchingperturbedcrystalstructures.PerturbationsincludedGaussiannoiseonfrac- tional atomic coordinates and on lattice vectors, as well as symmetry operations and translations. Their results (Table 1 and Figure 3) show that CLOUD maintains mod- erate robustness to small perturbations in atomic positions and lattice strain, though methods explicitly designed for perturbation invariance (e.g., BAWL) exhibit higher tolerance at large noise levels.\nTable 1: Comparison of different methods matching pairs of disordered structures, adapted from Siron et al. [17]. The success rate is reported on all the combinations of pairwise similarity for different chemical formulas.\nBAWL Short-BAWL EqV2-sim PDD Pymatgen SCOPE SLICES\n0.14±0.18 0.30±0.33 0.61±0.38 0.56±0.46 0.00±0.01 0.46±0.33 0.00±0.00\nThis independent analysis provides evidence that our SCOPE representation is reasonably stable under minor structural distortions, even though such distortions may change the exact space group classification and thus the tokenized description. We will cite these results in the revised manuscript and add a short discussion in the Methods section to acknowledge current robustness performance and potential future work.\nIn addition to the existing analysis on the robustness on the SCOPE representa- tion, we have also investigated in the robustness of model predictions with respect to differentnoiselevelsonatomicpositionsandlatticevectors.Similartotheanalysisby Sironetal.[17],weappliedGaussiannoiseofdifferentlevelsofperturbationstomimic atomic coordinates noise or the lattice strain, and report how the model performance varies with the standard deviation σ of the noise.\nWe have added the following results and discussions on representation and model robustness under small perturbation of crystal structures in the section S3.6 Robust- ness Analysis in the revised Supplementary Information:\nIn addition to the property prediction results, we investigate the robustness of building SCOPErepresentationsandCLOUDmodelpredictionsundercontrolledperturbationsof the input crystal structures. Such perturbations are practically unavoidable in realistic materialsresearch,arisingfromgeometryoptimizationtolerances,noiseinmachine-learned force fields, or the thresholds used by symmetry analyzers. Because the SCOPE repre- sentation relies on symmetry parsing, even small distortions can alter the assigned space group or Wyckoff positions. We therefore probe how predictive accuracy and parsing stability degrade under systematically applied noise to atomic coordinates and lattice vectors, providing a quantitative measure of the sensitivity of our models to structural perturbations.\nA recent independent study by Siron et al. [17] provides an external benchmark of this issue. They systematically compared multiple structure identification and fingerprinting approaches—including the SCOPE represention for CLOUD, Pymatgen’s Structure- Matcher, SLICES, and their proposed BAWL hashing—for matching perturbed crystal structures.PerturbationsincludedGaussiannoiseonfractionalatomiccoordinatesandon lattice vectors, as well as symmetry operations and translations. Their results show that SCOPE (and thus CLOUD) maintains moderate robustness to small perturbations, out- performingotherrepresentationssuchasSLICES,whilespecializedperturbation-invariant methods such as BAWL achieve higher tolerance at larger noise levels. This analysis establishes that SCOPE is reasonably stable under minor structural distortions.\n2\nthis perturbs both lattice lengths and angles (normal and shear). For every perturbed structure we regenerate SCOPE and evaluate fixed checkpoints (no retraining) over a log gridσfrom10−4 uptoatmost100.Alongsidepredictionerror(Figure4a-c)wereportthe symmetry-parse success rate (Figure 4d). On the gvrh dataset from MatBench [19], the MAE for unperturbed structures is 0.0873±0.0028. Fractional coordinate noise induces a monotonic increase—0.1004 at σ = 3×10−4 (+15%) and 0.1201 at 10−3 (+37.6%)— whereaslatticestrainshowsmuchlesssignificantimpactonmodelperformanceupto10−2 (MAE≈0.090,≤+4%).Parsingsucceedsforallsmall-σ cases(rate=1.0through10−2);\ndegradation appears only under heavy stress (e.g., at σ = 0.1 the success rate drops to 0.774 for coordinate noise and 0.987 for strain, and to 0.075 / 0.399 at σ = 0.3). On the Cv dataset, CLOUD is markedly more sensitive to coordinate noise: 8.92→11.00 at 10−4 (+23%),19.71at3×10−4 (+121%),and27.04at10−3 (+203%);straingrowsmoregently (12.29 at 10−3, 13.66 at 10−2). Symmetry parse success remains ≈ 0.997 up to 10−2, indicatingthatourmethodofconstructingSCOPEforcrystalstructuresisrobusttosmall perturbation levels. Meanwhile, we need to point out that the symmetry parsing success doesnotnecessarilyguaranteethatthetruesymmetryfortheunperturbedstructurecould be recovered. Adding the physics head in CLOUD–DEBYE both lowers baseline error (3.19±0.04) and dampens sensitivity: under coordinate noise at σ = 10−3 the MAE is 5.27 (versus 27.04 for CLOUD), and under strain at σ = 10−2 it is 4.66 (versus 13.66).\nOverall,theseresultsdemonstratethatCLOUDisreasonablyrobusttomodeststructural perturbations, with lattice strain having limited effect and coordinate noise being more impactful, while the physics-informed CLOUD–DEBYE variant substantially mitigates sensitivity and enhances stability across perturbation regimes.\nCommet 4: In addition to the debye model, what are some other properties of crys- tal that can modeled in this thermodynamic consistent and differentiable way? The authors should include some discussion about it.\nAuthor Reply: We thank the reviewer for the comment. In addition to the Debye heat capacity, we have extended our CLOUD-Debye framework to predict the vibra- tional entropy S . This quantity arises consistently from the Debye model and is a vib distinct thermodynamic observable beyond C .\nv\nBeyond the Debye model, the same framework can naturally be extended to other properties governed by well-established physics laws. In principle, our differentiable physics framework is generalizable as long as there exist enough training data and the physics model could be implemented in a differentiable way. In such a framework, 13\n0.04\n0.03\n0.02\n0.01\n0.00\n10 4 10 3 10 2 10 1 10 0\nEAM\na\n20\nStrain 10\n5\n0\n10 4 10 3 10 2 10 1\nEAM\nb\nCoord noise\nStrain\n2.0\n1.5\n1.0\n0.5\n0.0\n10 4 10 3 10 2 10 1\nEAM\nc\n1.0\n0.8\n0.6\n0.4\nCoord noise 0.2\nStrain\n0.0\n10 4 10 3 10 2 10 1 10 0\netar\nsseccus\nesraP\nd\ngvrh coord\ngvrh strain\nCv coord\nCv strain\nFig. 4: Model performance robustness analysis. Gaussian noise with varying σ’s are added to the atomic positions or lattice strain applied to lattice vectors. We show the ∆MAE due the noise for a gvrh with CLOUD model; b Cv with CLOUD model; c Cv with CLOUD-DEBYE model. In addition, the symmetry parsing successful rate under different noises are presented in d.\nthe differentiable physics model provide constraints as well as guidance to enforce physics consistency, while the complex mapping will be learned by the machine learn- ing model in a data-driven manner. We have added the discussion on future work of extending the framework to a wider range of physics models, investigating the possibility of integrating CLOUD with the (i) pDOS for more general and accurate predictions of phonon-related properties; (ii) predicting elastic modulus from elastic tensors;(iii)predictingthermalexpansioncoefficientthroughtheGru¨neisenparameter γ and temperature-dependent bulk modulus. These cases demonstrate the generality of our approach in coupling learned representations with differentiable physical laws.\nWe have added one concise paragraph at the end of the section CLOUD-DEBYE: Integration of CLOUD in a differentiable and physics-consistent materials thermodynamics framework in Results:\nThis paragraph forms a brief outlook for extending the CLOUD-DEBYE frame- work we proposed.\nMeanwhile, we have also added the following discussion in S3.5 Extension of CLOUD-DEBYE: Towards Broader Applications in the revised Supplementary Infor- mation to demonstrate the generalizability of our differentiable physics approach in materials modeling.\nWe have extended our CLOUD-DEBYE framework for predicting vibrational entropy (S ) to showcase its broad application in many other property prediction tasks. Table 2 vib summarizes the results for predicting S . Similar to our implementations for predict- vib ing constant-volume heat capacity Cv and phonon internal energy U, the train and test data are under 300K, and we use MAE/MAD as the metric for comparison. CLOUD- DEBYEoutperformsthetwodescriptor-hybridizedGNNs(de-CGCNNandde-ALIGNN) and CLOUD, which is consistent with the observations for Cv and U.\nTable 2: Results on predicting S . All\nvib\nresults are shown in the mean of test\nMAE/MAD along with the standard devi-\nation in brackets. The best results are in\nboldandthesecond-bestresultsareunder-\nlined.\nModel S (↓)\nvib\nde-CGCNN 0.154(0.003)\nde-ALIGNN 0.195(0.005)\nCLOUD 0.212(0.009)\nCLOUD-DEBYE 0.121 (0.004)\na b\nc d\nFig. 5: S predictions by CLOUD-DEBYE for a Al O , b Li O, c CaTiO , and\nvib 2 3 2 3\nd GaN. The CLOUD-DEBYE model, fine-tuned on S data at 300K, accurately\nvib\nextrapolates heat capacity across a wide temperature range, approaching the melting point of each material. Model predictions are compared against experimental and DFT-calculated results.\nOur experimental results present a compelling demonstration of the significance of inte- grating CLOUD with physical models for physics-consistent predictions. Moreover, our framework of embedding physics-based, differentiable laws into machine learning models is not limited to the Debye model. In principle, this approach can be generalized to any physical system, provided that the underlying physics can be implemented in a differen- tiablemannerandsufficientdataareavailabletolearnthecomplexmappingsthatremain intractable in the physics model. In the following, we will briefly discuss the possible extensions of the CLOUD-physics framework which can serve as future work:\n16\n• Phonon-related properties from pDOS While the Debye model provides a quadratic approximation to the phonon spectrum and thus misses optical branches and fine spectral features, these can be captured by predicting the full phonon density of states (pDOS). In this setting, the ML model predicts the pDOSg(ω)directlyfromthecrystalrepresentation,whichisthenfedintoanalytic statistical–mechanics relations for vibrational thermodynamics. For example, the constant–volume heat capacity C is given by\nv\n(cid:34) (cid:35)\n(cid:90) ∞ 1 ℏω\nC (T)= ℏωg(ω) e ℏω/kBT dω, (5)\nv 0 (cid:0) eℏω/kBT −1 (cid:1)2 k B T2\nand the vibrational entropy by\n(cid:90) ∞ (cid:20) ℏω/k T (cid:16) (cid:17)(cid:21)\nS (T)=k g(ω) B −ln 1−e−ℏω/kBT dω. (6)\nvib B\n0\neℏω/kBT −1\nTherefore, the mapping from pDOS to thermodynamic observables is fully dif- ferentiable, allowing gradients to backpropagate through the integrals, and (ii) a single ML–predicted g(ω) provides access to multiple vibrational quantities (e.g. C , S , free energy), enabling multitask training and richer physics constraints.\nv vib\n• Elastic modulus: Elastic properties can be described by closed-form averages of the elastic tensor C . The Voigt bounds depend linearly on C:\nij\nV 9 11 22 33 12 23 13\nV 15 11 22 33 12 23 13 44 55 66\nThe Reuss bounds depend on the compliance tensor S =C−1:\nR S +S +S +2(S +S +S )\n11 22 33 12 23 13\nR 4(S +S +S )−4(S +S +S )+3(S +S +S )\n11 22 33 12 23 13 44 55 66\nTherefore, the Voigt–Reuss–Hill averages are\nK = 1(K +K ), G = 1(G +G ) (11)\nVRH 2 V R VRH 2 V R\npropagate K or G back to update the model parameters.\nVRH VRH\n• Thermal expansion coefficient: The thermal expansion coefficient relates to\nthermodynamic quantities through the Gru¨neisen parameter γ and temperature- dependent bulk modulus [25]:\nγC (T)\nα(T)= v (12)\nK (T)V\nT molar\nwhere C (T) is the heat capacity from the Debye model and the bulk modulus\nv\nvaries with temperature according to:\nK (T)=K (1−α T) (13)\nT 0 K\nThe linear thermal expansion follows:\n∆L (cid:90) T\n= α(T′)dT′ (14)\nL\nThemodelpredictsthezero-temperaturebulkmodulusK ,itstemperaturecoeffi-\n0\ncientα ,andtheGru¨neisenparameterγ.Sincealloperationsaresmoothandthe\nK\nthermalexpansioncoefficientdependsdifferentiablyonthesepredictedquantities,\ngradients flow seamlessly from experimental thermal expansion measurements\nback to optimize the model parameters for representation learning.\nTogether,theseexamplesillustratethebreadthofourframework:CLOUD–pDOScaptures\nfull vibrational spectra for Cv and S\nvib\n, VRH averaging provides smooth and bounded\nmappingsforelasticmodulus,andthermodynamicrelationsyieldtemperature-dependent\nthermalexpansion.ThesecaseshighlightthebroaderpromiseofcombiningCLOUDwith\nphysical laws wherever sufficient training data are available.\nReviewer #2 Comments\nComment: This work introduces a transformer-based framework for crystal repre- sentation learning. The representation is SCOPE, a coordinate-free crystal string that concatenates space-group generator symbols (rather than just space group numbers), 18 occupied Wyckoff positions, and composition information. The model is pretrained on OPTIMADE structures and fine-tuned on a wide suite of downstream mate- rials property benchmarks (MatBench, MatBench-Discovery, UnconvBench). The reported results show a broadly competitive performance compared to graph neu- ral networks (GNNs) and better than composition-based models. Attention analysis suggest the model attends strongly to symmetry tokens, supporting interpretability claims. The authors also investigate the scaling-law fits related to the data/model size to performance, offering guidance for future scaling. The authors further present a CLOUD-DEBYE that predicts a Debye temperature and propagates it through a dif- ferentiable Debye model to obtain heat capacity and internal energy curves, improving accuracy while enforcing thermodynamic consistency. Overall, the work is comprehen- sive and timely, covering symmetry-aware representations, foundation models, and physics-informedMLformaterials.However,theworkcanbebetterorganizedtoshow better novelty and clarity. I suggest that authors have a revision before the work can be accepted.\nAuthor Reply: We sincerely thank the reviewer for the comprehensive and con- structive evaluation of our work. We are grateful for the recognition of the breadth and timeliness of our study, including the development of the SCOPE representation, the CLOUD and CLOUD-DEBYE models, and the exploration of interpretability and scaling laws. We also appreciate the reviewer’s thoughtful suggestion regarding organization and clarity, and we have carefully revised the manuscript to better highlight the novelty of our work and improve the overall presentation.\nComment 1: One key innovation in this study is the replacement of space group numbers with string generator with symmetry operations. The concept seems appeal- ing and makes sense intuitively, but an ablation study is missing. It is unclear how much gain comes from (a) generator-level symmetry tokens vs space group numbers, (b) inclusion of Wyckoff positions. Comparisons to MatInFormer are confounded by dataset size, tokenizer scope, and model capacity differences.\nAuthor Reply: We thank the reviewer for the comment. We agree that the design choiceofourSCOPErepresentationisintuitivelyplausiblewhilewecanprovidemore emprirical evidence to support the claim. We also agree that the direct comparison to MatInFormer, though partly suggesting the benefit of using generator strings com- pared to space group numbers or Hermann-Mauguin notations, is not comprehensive enough given the different pre-training data and model capacity. Therefore, we have designed different variants of SCOPE:\n19\n• Varying space group notations: index (1-230), Hermann-Mauguin notations, or Generator strings.\n• Inclusion of Wyckoff positions (WPs): using WPs in SCOPE or not.\nIn addition, we experimented with another variant which has composition-only for comparison, excluding space group and Wyckoff position information. Then we pre- train the model with all those variants and fine-tune them on two representative datasets from MatBench [19]: dielectric and e form, which have 4764 and 132752 data points respectively. We use the empirical result to demonstrate the effectiveness of using generator strings as the space group notations and the inclusion of Wyckoff positions for improved model predictive accuracy.\nWe have pre-trained CLOUD for each of the SCOPE variants and fine-tuned the models on dielectric and e form datasets. As shown in Table 3, across both datasets, including WPs consistently improves performance, regardless of the space group notation used. Meanwhile, generator strings with WPs achieve the best perfor- mance, yielding the lowest MAE on both dielectric (0.3038) and e form (0.0542). Using generator strings without WPs ranks second on dielectric, suggesting that generator-level symmetry tokens alone provide a tangible gain over index or HM notations.Thecomposition-onlyvariantperformsworst,confirmingthevalueofstruc- tural symmetry information. These results quantitatively validate the benefit of both generator-level symmetry tokens and the inclusion of WPs, with their combination offering the largest improvement.\nAll the above results have been merged into the section S3.7 Ablation Studies in the revised Supplementary Information, which is attached after the response to Comment 2.\nComment 2: “Long-range order/interactions” claim is underspecified and under-validated.ThemanuscriptattributesCLOUD’sadvantageonUnconvBenchand phonon thermodynamics to capturing long-range effects but provides no operational definition or ablation test. Prior work shows that mainstream GNNs often fail to captureperiodicityandlong-rangeinformation.Thestatementwillbemuchstrongerifthe impact of the string generator of space group numbers on disordered materials. It will be nice if the authors can conduct a controlled ablation study on (a) SG-number only, (b) SG Hermann–Mauguin tokens, (c) Full SCOPE w/o Wyckoff, (d) Symmetry-only (maskcomposition),(e)Composition-onlyforperfectanddisorderedcrystalsdatasets.\nAuthor Reply: We thank the reviewer for the comment. To address the concern, we have conducted a controlled ablation study on different space group notations (index, 20 Hermann–Mauguin, generator string) and on the inclusion or exclusion of Wyckoff positions,similartotheablationstudiesdescribedinourresponsetoComment1.We didnotincludethesymmetry-onlyvariantintheablationstudy,thoughsuggestedby thereviewer,assymmetrywithoutcompositiondoesnotdefineaphysicallymeaningful materialandthuswouldnotprovideaninterpretablebaseline.Webenchmarkedthese variants of SCOPE on two datasets: the defected dataset from UnconvBench [26], which evaluates the representation’s ability to model defected materials, and the C v dataset from Gong et al. [27], which probes the capacity to capture long-range inter- actions. The experimental results are summarized in Table 3 and demonstrate the relative contributions of different symmetry encodings and the inclusion of Wyckoff positions.\nFor the defected dataset, we observe that the SCOPE variants based on the generator string outperform those based on space group index or Hermann–Mauguin notation,indicatingthatthegeneratorstringindeedprovidesaricherencodingofsym- metry. Interestingly, the inclusion of Wyckoff positions slightly degrades performance on defected structures. This trend is likely due to how the dataset was constructed: vacancieswereintroducedbyremovingatomswithuniqueWyckoffpositions[28],but the specific defect site is not annotated. As a result, including Wyckoff positions may introduce noise rather than additional information in this task. For the C dataset, v the composition-only baseline performs poorly, confirming that composition alone is insufficient to model the structure-property relationship. However, we see a substan- tialimprovementoncesymmetryinformationisincorporated,regardlessofthespecific notation, and the generator string yields the best overall results. This demonstrates that symmetry encodings play a central role in enabling the model to capture ther- modynamicpropertieslinkedtolong-rangeinteractions.Overall,theseresultssupport our design choice of using generator strings while highlighting potential limitations when Wyckoff information is incomplete or ambiguous.\nCombining our response to Comment 1 & 2 from the reviewer, we have added the following results and discussion regarding the ablations studies in the section S3.7 Ablation Studies in the Supplementary Information:\nTo systematically evaluate the design choices underlying our SCOPE representation, we conduct a series of ablation studies across multiple datasets. Specifically, we vary (i) the notation used to encode space group symmetry—numerical index (1–230), Her- mann–Mauguin (HM) symbols, or generator strings (GS) of symmetry operations—and 21 (ii)whetherWyckoffpositions(WPs)areincluded.Wealsoexaminedacomposition-only baseline,wheresymmetryinformationisignoredentirely,toassesstheaddedvalueofstruc- turalencodings.Foreachvariant,wepre-trainCLOUDandfine-tunetheresultingmodels on four representative benchmark datasets: dielectric and e form from MatBench [19], the defected dataset from UnconvBench [26], and the Cv dataset from Gong et al. [27]. These tasks span both conventional property prediction benchmarks and settings where long-range interactions or structural disorder are especially relevant.\nTheresults,summarizedinTable3,yieldseveralconsistenttrends.First,acrosstheMat- Bench tasks, including WPs consistently improves accuracy regardless of the space group notation, with generator strings plus WPs achieving the best performance overall. This indicates that both generator-level tokens and fine-grained Wyckoff position information arecrucialinmodelingthematerialstructures.Meanwhile,forthedefecteddataset,gen- erator strings again provide the strongest symmetry encoding, but the inclusion of WPs slightlyreducesperformance.Thiscanbeattributedtothedatasetconstruction:vacancies areintroducedatuniqueWyckoffsiteswithoutexplicitannotationofdefectlocation[28], making WP tokens a source of noise rather than signal. Finally, in the Cv task, the composition-onlybaselineperformspoorly,underscoringtheinsufficiencyofstoichiometry alone for material representation learning. Incorporating symmetry information markedly improvesaccuracy,withgeneratorstringsdeliveringthelowesterrors.Together,thesefind- ings quantitatively validate the benefit of generator-level symmetry encoding while also clarifying the contexts in which Wyckoff positions provide complementary information.\nTable 3: Ablation study on the SCOPE representation across four datasets with dif- ferent variants of SCOPE design. Results are reported as the mean test MAE for dielectric, e form, and defected and MAE/MAD for C (and standard devia- v tion) from five-fold cross-validation. Bold indicates the best performance within each dataset, and underline denotes the second-best.\nSGnotation WP dielectric(↓) eform(↓) defected(↓) Cv (↓)\nN/A ✗ 0.3438(0.0916) 0.0845(0.0007) 0.8544(0.1894) 0.48(0.01)\n✗ 0.3551(0.0974) 0.0599(0.0006) 0.7769(0.1946) 0.17(0.00)\nIndex\n✓ 0.3345(0.0806) 0.0568(0.0004) 0.9360(0.1897) 0.17(0.01)\n✗ 0.3333(0.0751) 0.0604(0.0003) 0.8188(0.1890) 0.18(0.00)\nHM\n✓ 0.3290(0.0894) 0.0569(0.0008) 0.8257(0.1539) 0.16 (0.01)\n✗ 0.3289(0.0771) 0.0601(0.0006) 0.7482 (0.1648) 0.17(0.01)\nGS\n✓ 0.3038 (0.0782) 0.0542 (0.0007) 0.8197(0.1093) 0.16 (0.00)\nComment 3: Distinguish “long-range order” (periodicity) from “long-range interac- tions” (physical couplings such as elastic or electrostatic); the manuscript alternates between them.\nWe have corrected the usage of terminologies in the revised manuscript accordingly:\n• In Introduction, the original text was ‘Figure 1d) and assess its ability to predict phonon internal energy and heat capacity, two properties that rely on long- range dependencies of crystal structures.’ We have changed ‘dependencies’ to ‘interactions’.\n• In the Transferable and scalable crystal representation learning with CLOUD section in Results, the original text was ‘... enabling it to effectively models long-rangedependenciesandcontextualrelationshipswithinsequences.’Wehave changed ‘dependencies’ to ‘interactions’.\n• In the CLOUD-DEBYE: Integration of CLOUD in a differentiable and physics- consistent materials thermodynamics framework section in Results, the original text was ‘... , which restrict their receptive fields and impair their ability to learn long-range dependencies.’ We have changed ‘dependencies’ to ‘orders’.\n• In the CLOUD-DEBYE: Integration of CLOUD in a differentiable and physics- consistent materials thermodynamics framework section in Results, the original text was ‘Unlike GNNs, they capture interactions between all tokens without relying on sequential message passing, thereby enabling the modeling of long- range dependencies.’ We have changed ‘dependencies’ to ‘orders’.\n• In the CLOUD-DEBYE: Integration of CLOUD in a differentiable and physics- consistent materials thermodynamics framework section in Results, the original text was ‘The results manifest CLOUD’s capability of learning long-range dependencies in crystal structures as well as the significance of introducing ther- modynamic consistency into the material machine learning framework.’ We have changed ‘dependencies’ to ‘interactions’.\n23\n• InDiscussion,theoriginaltextwas‘Moreover,CLOUDiscapableoflearninglong- rangestructuraldependencies.’Wehavechanged‘dependencies’to‘information’.\nComment 4: Clarify de-duplication rule (“smallest volume per formula + SG”) and report how many items dropped.\nAuthorReply:Wethankthereviewerforpointingouttheneedforclarification.Our deduplicationmethod,adaptedfromAntunesetal.[13],groupsstructuresbyidentical chemical composition and space group number, and retains only the entry with the smallest volume per formula unit in each group. This approach, which preferentially selects the most compact configuration, removes approximately 51.5% of entries from the raw dataset, leaving ∼6.3 M structures for pre-training. The revised manuscript now includes this description and the corresponding statistics in the Methods section for clarity.\nWe have added the following clarification in the first paragraph of the Datasets section in Methods:\nTo remove duplicate entries, we apply the deduplication rule from Antunes et al. [13]: structuresarefirstgroupedbyidenticalchemicalcompositionandspacegroupnumber,and withineachgroupweretainonlythestructurewiththesmallestvolumeperformulaunit. Thiscriterionfollowstheobservationthatstructureswithidenticalcompositionandspace group but larger volumes are often redundant polymorphs or duplicates of less compact configurations. Applying this rule to our dataset removes ∼51.5% of the total structures, resulting in ∼6.3M structures for subsequent processing.\nComment 5: Some relevant work are not cited. A few generative symmetry-aware representations e.g. WyCryst, WyFormer, and SymmCD—demonstrate the impor- tance of explicit symmetry and Wyckoff encodings for both generation and property prediction.\nAuthor Reply: We thank the reviewer for the comment. We acknowledge that the reviewer mentioned very important and inspiring literature which are highly rel- evant to our work. Therefore, we have added discussion about the importance of explicit symmetry and Wyckoff encodings and cited the relevant work in our revised manuscript.\nWe have added the following discussion into the revised manuscript, which could be found in the first and last paragraph in the Symmetry-Consistent Ordered Param- eter Encoding (SCOPE) crystal representation section in Results:\n“Meanwhile, recent work has demonstrated the importance of explicitly incorporat- ing crystallographic symmetry and Wyckoff information into generative and predictive models. For example, WyCryst [29] introduces a Wyckoff-based representation that pre- serves symmetry during generation and relaxation; WyFormer [30] encodes crystals as 24 fused element–Wyckoff tokens to enable symmetry-conditioned autoregressive generation; SymmCD[31]decomposescrystalsintoasymmetricunitsandsymmetryoperations,ensur- ingsymmetry-preservingdiffusionevenforrarespacegroups.Theseworksunderscorethat symmetry-aware, Wyckoff-based encodings serve as effective and expressive representa- tions for crystal structures, substantially improving the physical validity and diversity of generated crystals.”\n“This provides a more granular and physically complete encoding than categorical site- symmetry tokens [30] or binary symmetry matrices [31]. Unlike WyCryst’s one-hot representationwithdegreesoffreedom[29],SCOPEcapturesthecompletesetofsymmetry generators in a compact symbolic form, enabling fine-grained distinction between closely relatedspacegroupsandnaturallyaligningwithmaskedlanguagemodelingforlarge-scale transformer pretraining.”\nReviewer #3 Comments\nComment: The authors present a framework for predicting the properties of crys- talline materials and compare it with multiple deep learning models, which is a hot topic. However, I have some concerns about the novelty of the framework and the performance comparison. The lack of application of some models in specific material systems also means that there is a lack of proof of the effectiveness of the models in practical applications.\nAuthor Reply: We sincerely thank the reviewer for the thoughtful assessment of our work and for recognizing the relevance of crystal property prediction as a timely and important research direction. We appreciate the reviewer’s concerns regarding the novelty of our framework, the performance comparisons, and the demonstration of practical applications. We have carefully addressed these points in our revision by clarifying the unique contributions of our approach, providing additional comparisons to various baselines, and expanding the discussion of applications to better support the effectiveness of the proposed model.\nAuthor Reply: We thank the reviewer for pointing out the need for clarification. The OPTIMADE dataset indeed contains duplicate structures that can arise from multiple database sources reporting the same material. To address this, we adopted the deduplication strategy from Antunes et al. [13], which identifies duplicates by groupingstructureswiththesamechemicalcompositionandspacegroupnumber,then retaining only the one with the smallest volume per formula unit. This step removed approximately51.5%ofentriesfromtherawOPTIMADEcollection,yielding∼6.3M unique structures for training.\nWe have added the following clarification in the first paragraph of the Datasets section in Methods:\nTo remove duplicate entries, we apply the deduplication rule from Antunes et al. [13]: structuresarefirstgroupedbyidenticalchemicalcompositionandspacegroupnumber,and withineachgroupweretainonlythestructurewiththesmallestvolumeperformulaunit. Thiscriterionfollowstheobservationthatstructureswithidenticalcompositionandspace group but larger volumes are often redundant polymorphs or duplicates of less compact configurations. Applying this rule to our dataset removes ∼51.5% of the total structures, resulting in ∼6.3M structures for subsequent processing.\nComment 2:Theauthorspre-trained onthe OPTIMADE dataset,evaluatedtheper- formance on the MatBench dataset, and compared it with other deep learning models (CGCNN, ALIGNN). Have CGCNN and ALIGNN also been pre-trained? If not, the performance comparison is unfair.\nAuthor Reply: We thank the reviewer for the comment. To clarify, the results by CGCNN and ALIGNN reported in our manuscript are not from pre-trained mod- els but trained from scratch in the supervised learning approach instead. While pre-training has become a well-established paradigm for transformer-based language models due to their natural compatibility with self-supervised tasks (e.g., masked language modeling which is used for pre-training CLOUD) and scalability across mil- lions of unlabeled structures, large-scale self-supervised pre-training for crystal GNNs remains non-trivial.\nWe emphasize two reasons why this comparison remains fair. Most available pre-trainedGNNs,suchasALIGNN,relyonsupervisedpre-trainingwithlabeledprop- erties(e.g.,energyorforcelabels,ortransferringamodeltrainedforonepropertytoa similar property), which is fundamentally different from our self-supervised approach 26 and requires costly labeled data. Second, recent attempts to bring self-supervised pre-training to GNNs, such as Crystal Twins [32], have shown only modest improve- mentsoversupervisedbaselines.CrystalTwinsadoptsatwin-networkframeworkwith CGCNN encoders, using data augmentations (random perturbations, atom masking, and edge masking) and objectives such as Barlow Twins or SimSiam losses. While the method shows some performance gains over a purely supervised CGCNN, these improvements are not substantial, nor consistent across all benchmark datasets, and remain far from surpassing state-of-the-art models such as ALIGNN. This highlights both the promise and the current difficulty of making GNN pre-training competitive at scale.\nTo further probe this issue, we also conducted self-supervised pre-training of CGCNN using the same ∼6.3M crystal structures employed for CLOUD, in con- trast to the 428K one used in the Crystal Twins paper [32]. In practice, the training proved prohibitively slow with the large dataset: after eight hours of computation on NVIDIAH100GPU,only 12%ofthepre-trainingstepsofoneepochwerecompleted, whereas CLOUD requires on the order of ten minutes per epoch under compara- ble conditions. This dramatic efficiency gap underscores the scalability advantage of transformer-basedmodelsforlarge-scalepre-training.Weacknowledgethatourimple- mentation may not be fully optimized—whether due to our own unfamiliarity with theCrystalTwinsframeworkorlimitationsinthecurrentcodebase.Nevertheless,the evidenceindicatesthatself-supervisedpre-trainingforGNNsremainsfarfrommature or practically effective at scale. Indeed, even in the Crystal Twins public repository (https://github.com/RishikeshMagar/Crystal-Twins/issues/2),theauthorsnotethat attempting ALIGNN pre-training leads to prohibitive compute and memory require- ments due to the construction of additional line graphs, and that their preliminary experiments on smaller datasets yielded negligible benefits.\nWe further support this conclusion with quantitative results on the MatBench GVRH dataset, summarized in Table 4. Crystal Twins yields only minor improve- mentsforCGCNN(from0.089to∼0.086–0.087),andourownlarge-scalepre-training resultsinevennegativegains,probablybecauseoftheunder-optimization.Incontrast, CLOUD shows a clear and consistent benefit from pre-training (0.117 → 0.087), with improvementsthatarebothstatisticallysignificantandcomparableinmagnitudetoor exceeding those of Crystal Twins. Even though the accuracy by pre-trained CGCNN and CLOUD is comparable, CGCNN requires atomic coordinates as input while our representation is coordinate-free, highlighting the advantage of CLOUD in effective representation learning. Note that gvrh is one of the few tasks from MatBench that 27 CLOUDisnotsignificantlyoutperformingthenon-pretrainedorpre-trainedCGCNN. Asamatteroffact,basedontheresultsfromCrystalTwinspaper[32],CLOUDout- performs the pre-trained CGCNN on four out of the eight regression tasks while even the non-pretrained CLOUD outperforms the pre-trained CGCNN on two tasks.\nTaken together, these findings suggest that while GNNs are strong supervised learners, their current pre-training strategies are less scalable and less effective com- pared to transformer-based language models like CLOUD. We therefore believe that thecomparisonisfairandmeaningful,asitreflectsthestateoftheartineachmodeling paradigm.\nWe have provided the following results and discussions in the section S3.8 Train- ing Efficiency and Scalability in the revised Supplementary Information:\nForbaselinecomparisons,wefollowtheMatBenchprotocolwhereCGCNNandALIGNN are trained from scratch in a supervised manner, without pre-training. While ALIGNN providespre-trainedcheckpoints,theseareobtainedviasupervisedpre-trainingonlabeled data such as energies and forces, in contrast to the self-supervised strategy adopted for CLOUD. Attempts to develop self-supervised pre-training for GNNs, such as the Crystal Twinsframework[32],havedemonstratedonlymodestandinconsistentimprovementsover supervisedCGCNNandremainlesseffectivethanALIGNN.WefurtherapplytheCrystal Twinsframeworktopre-trainCGCNNonthesame∼6.3McrystalstructuresasCLOUD, but find that the pre-training process is prohibitively slow: after eight hours, only ∼12% of the training steps in one epoch are completed, compared to ∼10 minutes per epoch for CLOUD. We acknowledge that our implementation may not be fully optimized, but the evidence nevertheless indicates that self-supervised pre-training for GNNs is still far from mature. The comparative results on the MatBench GVRH dataset (Table 4) further supportthispoint:whileCGCNNwithCrystalTwinspre-trainingachievesonlymarginal gains over its supervised counterpart, CLOUD demonstrates a substantial and consistent improvementwhenpre-trained.Theseobservationshighlightboththelimitedeffectiveness andpoorscalabilityofGNN-basedpre-training,underscoringwhytransformer-basedmod- els like CLOUD are particularly well-suited for large-scale self-supervised learning with crystals.\nComment 3: Large-scale pre-training is very expensive. Have the authors considered the issue of computational time? training time, testing time, model parameters.\nAuthor Reply: We thank the reviewer for the comment. For the results reported in the manuscript, we performed pre-training on the deduplicated OPTIMADE dataset (∼6.3M sequences). On two NVIDIA H100 GPUs, the base CLOUD model (∼110M parameters) required only ∼10 minutes per epoch, which is areasonablecostcompared 28\nTable 4: Comparison of model performance on the MatBench\ngvrhdataset.Resultsarereportedasthemeanoffive-foldMAE\nwith standard deviations. Results for CGCNN (supervised) and\nCrystal Twins variants are taken from [32]; CLOUD and repli-\ncated CGCNN-CT results are from this work.\nModel Pre-training MAE±std\nCGCNN[32] None(supervised) 0.089±0.001\nCGCNN(CrystalTwins)[32] SSL(CTBarlow) 0.086±0.004\nCGCNN(thiswork) SSL(6.3M,CTcode) 0.092±0.004\nCLOUD(thiswork) None(supervised) 0.117±0.003\nCLOUD(thiswork) SSL(6.3M,MLM) 0.087±0.003\nto other large-scale machine learning models for materials science. More importantly, large-scale pre-training with CLOUD is computationally efficient due to the excellent scalability of transformer-based language models. With access to larger computing resources, the wall-clock training time can be further reduced with minimal loss of training efficiency. Indeed, benchmarking on Aurora, the exascale supercomputer at Argonne National Laboratory, confirmed that CLOUD maintains high parallel efficiency up to 64 nodes, even without specialized adaptations for large-scale compute. The results will be added to our revised manuscript. Once pretrained, the model can be reused across all fine-tuning tasks, thereby amortizing the one-time pre-training cost. As the number of downstream applications increases, the relative expense of pre-training becomes almost negligible.\nTofurtherassesscomputationalefficiencyinpracticalsettings,webenchmarkfine- tuning and inference costs relative to common GNN baselines on two representative benchmarks (Jarvis and MatBench FE). Table 5 shows that CLOUD achieves shorter training times per epoch and significantly faster inference, despite its larger model parameter count. These results demonstrate that CLOUD is not only scalable in pre- trainingbutalsocomputationallyefficientandlightweightindownstreamdeployment, indicating its practicality for large-scale discovery pipelines.\nFurthermore, compared to graph-based models that require explicit 3D atomic coordinates or costly DFT relaxations, CLOUD substantially lowers the end-to-end computational burden for large-scale discovery workflows. Most state-of-the-art modelsformaterialstructuremodelingarebasedongraphneuralnetworks(GNNs),which faceintrinsicscalingchallengessuchasoversmoothingandoversquashing[33].Inaddi- tion,distributedtrainingofGNNsrequiresgraphpartitioning,whereboundarynodes 29 and edges induce significant cross-device communication overhead. In contrast, trans- formers benefit from well-established parallelization strategies that scale efficiently with modern hardware.\nOverall,whilepre-trainingCLOUDrequiresanon-trivialbutone-timeinvestment, the efficiency of training at scale, and reusability across tasks, justifies the cost. We therefore view pre-training as a worthwhile tradeoff for enabling efficient and scalable downstream applications in materials science tasks.\nWe have added the following results and discussions in the section S3.8 Training Effiency and Sccalability in the revised Supplementary Information:\nThe computational cost of large-scale pre-training is an important consideration for prac- tical application of foundation models. For the results reported in this work, we pre-train CLOUDonthededuplicatedOPTIMADEdataset(∼6.3Msequences).UsingtwoNVIDIA H100 GPUs, the base model (∼110M parameters) requires ∼10 minutes per epoch. Such a cost is modest compared to other large-scale machine learning models in materials sci- ence and is further mitigated by the fact that pre-training is a one-time investment. Once obtained, a pretrained checkpoint can be reused across diverse fine-tuning tasks, consequently amortizing the initial cost over a broad range of downstream applications.\nToevaluatethescalingefficiencyofCLOUDontheAuroraHPC,theexascalesupercom- puter at Argonne National Laboratory, we measure the model throughput while scaling the number of nodes for two versions of the CLOUD model: a 110M parameter baseline and a 1B parameter larger model. The weak scaling setup fixes the local batch size per GPU,allowingustoisolatetheeffectofcommunicationoverheadastheclustersizegrows. As shown in Figure 6a, both models sustain high efficiency up to 128 nodes, maintain- ingefficiencyhigherthan60%.Notably,the1Bmodelconsistentlyoutperformsthe110M model in scalability, particularly at higher node counts, benefiting from a larger compute footprint that better amortizes communication overhead. We evaluate the strong scaling behavioroftheCLOUDmodelonAuroraforthetwoversionsofmodel.Inthissetup,the totaltrainingworkloadremainsconstantasthenumberofnodesincreases,andwerecord the reduction in runtime to compute speedup. As shown in Figure 6b, the 110M model shows near-linear speedup up to 64 nodes and reasonable scaling beyond. The 1B model, evaluated from 32 to 128 nodes, demonstrates consistent gains, though with diminishing returns at higher scales due to communication overhead. The log-log plots highlight the efficiency of distributed training, validating our stack’s readiness for exascale foundation model training.\nFurthermore, compared to graph-based models that require explicit 3D atomic coordi- nates or costly DFT relaxations, CLOUD substantially lowers the end-to-end computa- tionalburdenforlarge-scalediscoveryworkflows.Moststate-of-the-artmodelsformaterial 30\n100\n80\n60\n40\n20\n0\n1 2 4 8 16 32 64 128 256\nNumber of Nodes\n)%(\nycneiciffE\n25\n24\n23\n22\nCLOUD-1B\nIdeal Scaling 20\n8 16 32 64 128 256\nNumber of Nodes\n(a) Scaling Efficiency vs. Number of Nodes\nfor training CLOUD with 110M and 1B\nparametermodelsonAurorausingourtrain-\ningstack.Efficiencyiscomputedunderweak\nscaling by fixing the local batch size per\nGPU. Both models maintain high parallel\nefficiencyupto128nodes(higherthan60%),\nwith the 1B model exhibiting better scala-\nbility at large node counts due to improved\ncompute-to-communication ratio.\npudeepS\nCLOUD-110M\nCLOUD-1B\nIdeal (110M)\nIdeal (1B)\n(b) Scaling Speedup vs. Number of Nodes\nfor training CLOUD with 110M and 1B\nparametermodelsonAurorausingourtrain-\ning stack. Speedup is computed relative to\nthe smallest measured node count (8 for\n110M, 32 for 1B), with ideal scaling shown\nas dashed lines. Both the two models mani-\nfest near-perfect strong scaling behavior up\nto64nodesfor110Mmodelor128nodesfor\n1B model.\nFig. 6: Scalability analysis for CLOUD on Aurora [34].\nstructuremodelingarebasedongraphneuralnetworks(GNNs),whichfaceintrinsicscal- ing challenges such as oversmoothing and oversquashing [33]. In addition, distributed training of GNNs requires graph partitioning, where boundary nodes and edges induce significant cross-device communication overhead. In contrast, transformers benefit from well-established parallelization strategies that scale efficiently with modern hardware.\nIn addition to pre-training efficiency, the practicality of CLOUD also depends on the cost of fine-tuning and inference. To further assess computational efficiency in practical settings, we compared training and inference times of CLOUD against widely used GNN baselines (ALIGNN [1], Matformer [35], MEGNet [6], M3GNet [36], TensorNet [37], and SO3Net [38]) on representative property prediction tasks (Jarvis [39] and MatBench [19] formation energy datasets). The results for ALIGNN and Matformer on Jarvis formation energytaskareobtainedfromYanetal.[35]whilethenumbersfortheMatBenchforma- tion energy task are collected from implementations with MatGL [40]. As summarized in Table 5, CLOUD achieves faster training time per epoch and significantly shorter infer- ence time despite being larger in model parameter count. For example, on the Jarvis FE task, CLOUD requires only 34 seconds per epoch and 19 seconds for full inference, out- performingALIGNNandMatformerwhileofferinganorder-of-magnitudefasterinference\nthroughput. Similarly, on the MatBench FE benchmark, CLOUD completes inference in\nonly7secondsonthetestset,thefastestamongalltestedmodels.Theseresultshighlight 31 that the transformer-based architecture not only scales efficiently during pre-training but also provides lightweight fine-tuning and inference in downstream applications.\nTable 5: Training and inference time comparison on material\nproperty prediction tasks. Reported numbers include training\ntime per epoch, inference time on the full dataset (for Jarvis\nformationenergytask)orthetestset(forMatBenchformation\nenergy task), and model size.\nModel Time/epoch Inference Parameters\nALIGNN(JarvisFE) 327s 156s 15.4M\nMatformer(JarvisFE) 64s 59s 11.0M\nCLOUD(JarvisFE) 34s 19s 110M\nMEGNet(MatBenchFE) N/A 11.137s 168K\nM3GNet(MatBenchFE) N/A 20.089s 228K\nTensorNet(MatBenchFE) N/A 13.694s 2.8M\nSO3Net(MatBenchFE) N/A 32.601s 341K\nCLOUD(MatBenchFE) 59.892s 7.003s 110M\nComment 4: PMCGNN’s performance on the JARVIS-DFT dataset is already\nhigher than ALIGNN, Matformer, etc. The authors should also compare.\nhttps://pubs.acs.org/doi/10.1021/acs.jcim.4c01200\nAuthor Reply: We thank the reviewer for the comment. We acknowledge that\nPMCGNN is indeed an important model with superior performance on many bench- marktasks.Therefore,wehavefine-tunedCLOUDonJARVIS-DFTdataset,following\nthe same data splitting strategies as described in Feng et al. [41].\nWe have fine-tuned CLOUD on the five datasets from JARVIS-DFT: formation\nenergy, bandgap (OPT), total energy, bandgap (MBJ), and E . We follow the\nhull\nsame dataset version and data-splitting strategy as described in Feng et al. [41]\nand report the results in Table 6. Our model CLOUD achieves competitive pre-\ndictive performance, consistently outperforming CGCNN across the five tasks and\nexhibits the state-of-the-art (SOTA) performance on the E dataset. For 4 out\nhull\nof the 5 tasks, PMCGNN stands out with the lowest MAE, reflecting its tailored\ngraph-basedinductivebias.Nevertheless,CLOUDreachescomparableaccuracywhile\nbeingcoordinate-freeandsymmetry-consistentbydesign,highlightingitscomplemen-\ntary strengths. In particular, achieving the best performance on E demonstrates\nhull\nCLOUD’s capability for stability predictions, which are of central importance for\nmaterials discovery.\nbaselines on JARVIS Dataset. The results for CGCNN, ALIGNN, Matformer, and\nPMCGNN are directly obtained from Feng et al. [41].\nformationenergy bandgap(OPT) totalenergy bandgap(MBJ) E\nModel hull\nev/atom ev ev/atom ev eV\nCGCNN 0.063 0.20 0.078 0.41 0.17\nALIGNN 0.0331 0.142 0.037 0.31 0.076\nMatformer 0.0325 0.137 0.035 0.30 0.064\nPMCGNN 0.0278 0.122 0.029 0.25 0.040\nCLOUD 0.0606 0.159 0.069 0.31 0.030\nWe have also added the results and related discussion in S3.4 JARVIS-DFT\nsection in the revised Supplementary Information:\nAdditionally, we have also fine-tuned CLOUD on JARVIS-DFT-2021 3D dataset [39].\nWe chose five crystal property prediction tasks: formation energy, bandgap (OPT), total\nenergy, bandgap(MBJ), and E for benchmarking. We follow the same data-splitting\nhull\nstrategy as described in Feng et al. [41], and compare our model performance with the\nfollowing baselines: CGCNN [5], ALIGNN [1], Matformer [35], and PMCGNN [41]. As\nshown in Table 6, our model CLOUD achieves competitive predictive performance, con-\nsistently outperforming CGCNN across the five tasks and exhibits the state-of-the-art\n(SOTA) performance on the E dataset. For 4 out of the 5 tasks, PMCGNN stands\nhull\nout with the lowest MAE, reflecting its tailored graph-based inductive bias. Nevertheless,\nCLOUDreachescomparableaccuracywhilebeingcoordinate-freeandsymmetry-consistent\nbydesign,highlightingitscomplementarystrengths.Inparticular,achievingthebestper-\nformance on E demonstrates CLOUD’s capability for stability predictions, which are\nhull\nof central importance for materials discovery.\nPlus, we have also cited this excellent paper as one of the pioneering work in this\nfield in the second paragraph of the Introduction in the revised manuscript:\nPioneeringworkhassignificantlyimprovedthepredictionaccuracyforcrystals[1,5,41–44].\nComment 5: The author also mentioned that the performance of the current model\ndepends on a large amount of labeled data, but the author’s pre-training-fine-tuning\nmodel also depends on downstream labeled data. The self-supervised model may better\nhandle this problem: CDSSL, https://doi.org/10.48550/arXiv.2408.17255.\nAuthor Reply: We thank the reviewer for the comment. CDSSL [45] is an effec-\ntive pre-training approach for predictive models, in particular, graph neural networks\n33\n(GNNs) with a pretext task based on recovering valid material structures when given\nperturbed versions of these structures.\nMeanwhile, we would like to clarify that CDSSL is, like CLOUD, a pre-\ntrain–finetune framework in which labeled data is still required for downstream\nproperty prediction. Both methods aim to improve data efficiency when fine-tuning\non limited labeled data, rather than eliminating the need for labeled data entirely. In\ncomparison, CDSSL pre-trains structure-based graph neural networks via a geomet-\nric denoising objective, while CLOUD pre-trains a transformer-based model on our\nSCOPE representation using masked language modeling. CDSSL has shown improve-\nments for GNNs by learning local geometric information, whereas CLOUD learns\nglobal symmetry and compositional patterns.\nIn addition, These approaches are complementary, and future work could explore\ncombining SCOPE-based global symmetry encoding with CDSSL-style geometric\ndenoising to enhance data efficiency further.\nWe have added the citation for the CDSSL paper and included the related discus-\nsion in the third paragraph in the Introduction section:\nAside from language-model-style pre-training, alternative self-supervised strategies have\nbeen explored in the crystal domain, such as Crystal Denoising Self-Supervised Learning\n(CDSSL)[45],whichpre-trainsgraphneuralnetworksbyperturbingatomicpositionsand\nlearning to recover the original interatomic distances.\nReferences\n[1] Choudhary, K., DeCost, B.: Atomistic line graph neural network for improved\nmaterials property predictions. npj Computational Materials 7(1), 185 (2021)\n[2] Riebesell, J., Goodall, R.E., Jain, A., Benner, P., Persson, K.A., Lee, A.A.: Mat-\nbench discovery–an evaluation framework for machine learning crystal stability\nprediction. arXiv preprint arXiv:2308.14920 (2023)\n[3] Jain,A.,Ong,S.P.,Hautier,G.,Chen,W.,Richards,W.D.,Dacek,S.,Cholia,S.,\nGunter,D.,Skinner,D.,Ceder,G.,etal.:Commentary:Thematerialsproject:A\nmaterials genome approach to accelerating materials innovation. APL materials\n1(1) (2013)\n[4] Wang, H.-C., Botti, S., Marques, M.A.: Predicting stable crystalline compounds\n34\nusing chemical similarity. npj Computational Materials 7(1), 12 (2021)\n[5] Xie,T.,Grossman,J.C.:Crystalgraphconvolutionalneuralnetworksforanaccu-\nrate and interpretable prediction of material properties. Physical review letters\n120(14), 145301 (2018)\n[6] Chen, C., Ye, W., Zuo, Y., Zheng, C., Ong, S.P.: Graph networks as a universal\nmachine learning framework for molecules and crystals. Chemistry of Materials\n31(9), 3564–3572 (2019)\n[7] Grygorenko, O.O., Radchenko, D.S., Dziuba, I., Chuprina, A., Gubina, K.E.,\nMoroz,Y.S.:Generatingmultibillionchemicalspaceofreadilyaccessiblescreening\ncompounds. iScience 23(11), 101681 (2020) https://doi.org/10.1016/j.isci.2020.\n101681\n[8] Bellmann, L., Penner, P., Gastreich, M., Rarey, M.: Comparison of\ncombinatorial fragment spaces and its application to ultralarge make-\non-demand compound catalogs. Journal of Chemical Information and\nModeling 62(3), 553–566 (2022) https://doi.org/10.1021/acs.jcim.1c01378\nhttps://doi.org/10.1021/acs.jcim.1c01378. PMID: 35050621\n[9] Patel,H.,Ihlenfeldt,W.-D.,Judson,P.N.,Moroz,Y.S.,Pevzner,Y.,Peach,M.L.,\nDelann´ee, V., Tarasova, N.I., Nicklaus, M.C.: Savi, in silico generation of billions\nof easily synthesizable compounds through expert-system type rules. Scientific\nData (2020)\n[10] Gavezzotti,A.:Arecrystalstructurespredictable?Accountsofchemicalresearch\n27(10), 309–314 (1994)\n[11] Merchant, A., Batzner, S., Schoenholz, S.S., Aykol, M., Cheon, G., Cubuk, E.D.:\nScaling deep learning for materials discovery. Nature, 1–6 (2023)\n[12] Zeni, C., Pinsler, R., Zu¨gner, D., Fowler, A., Horton, M., Fu, X., Shysheya, S.,\nCrabb´e,J.,Sun,L.,Smith,J.,etal.:Mattergen:agenerativemodelforinorganic\nmaterials design. arXiv preprint arXiv:2312.03687 (2023)\n[13] Antunes,L.M.,Butler,K.T.,Grau-Crespo,R.:Crystalstructuregenerationwith\nautoregressive large language modeling. arXiv preprint arXiv:2307.04340 (2023)\n35\n[14] Gruver,N.,Sriram,A.,Madotto,A.,Wilson,A.G.,Zitnick,C.L.,Ulissi,Z.:Fine-\ntunedlanguagemodelsgeneratestableinorganicmaterialsastext.arXivpreprint\narXiv:2402.04379 (2024)\n[15] Barroso-Luque,L.,Shuaibi,M.,Fu,X.,Wood,B.M.,Dzamba,M.,Gao,M.,Rizvi,\nA.,Zitnick,C.L.,Ulissi,Z.W.:Openmaterials2024(omat24)inorganicmaterials\ndataset and models. arXiv preprint arXiv:2410.12771 (2024)\n[16] Deng, B., Zhong, P., Jun, K., Riebesell, J., Han, K., Bartel, C.J., Ceder, G.:\nChgnet as a pretrained universal neural network potential for charge-informed\natomistic modelling. Nature Machine Intelligence 5(9), 1031–1041 (2023)\n[17] Siron,M.,Djafar,I.,Fayet,E.,Rossello,A.,Ramlaoui,A.,Duval,A.:Lemat-bulk:\naggregating, and de-duplicating quantum chemistry materials databases. In: AI\nfor Accelerated Materials Design-ICLR 2025 (2025)\n[18] Xu, C., Zhu, S., Viswanathan, V.: Cloud: A scalable scientific foundation model\nforcrystalrepresentationlearning.In:Neurips2024WorkshopFoundationModels\nfor Science: Progress, Opportunities, and Challenges (2024)\n[19] Dunn, A., Wang, Q., Ganose, A., Dopp, D., Jain, A.: Benchmarking materials\nproperty prediction methods: the matbench test set and automatminer reference\nalgorithm. npj Computational Materials 6(1), 138 (2020)\n[20] Furukawa, G.T., Douglas, T.B., McCoskey, R.E., Ginnings, D.C.: Thermal prop-\nerties of aluminum oxide from 0 to 1200 k. Journal of research of the National\nBureau of Standards 57(2), 67–82 (1956)\n[21] Stull,D.R.:JANAFThermochemicalTables...vol.1.Clearinghouse,Washington,\nD.C. (1965)\n[22] Woodfield, B.F., Shapiro, J.L., Stevens, R., Boerio-Goates, J., Putnam, R.L.,\nHelean, K.B., Navrotsky, A.: Molar heat capacity and thermodynamic functions\nfor catio3. The Journal of Chemical Thermodynamics 31(12), 1573–1583 (1999)\n[23] Leitner, J., Strejc, A., Sedmidubsky, D., Ruzicka, K.: High temperature enthalpy\nand heat capacity of gan. Thermochimica acta 401(2), 169–173 (2003)\n[24] Loew, A., Sun, D., Wang, H.-C., Botti, S., Marques, M.A.: Universal\nmachine learning interatomic potentials are ready for phonons. arXiv preprint\n36\narXiv:2412.16551 (2024)\n[25] Gru¨neisen, E.: Theorie des festen zustandes einatomiger elemente. Annalen der\nPhysik 344(12), 257–306 (1912)\n[26] Wang, H., Sun, J., Liang, J., Zhai, L., Tang, Z., Li, Z., Zhai, W., Wang, X.,\nGao, W., Gong, S.: Crystograph: A comprehensive predictive model for crystal\nmaterial properties and the benchmark. Battery Energy, 70004 (2024)\n[27] Gong, S., Yan, K., Xie, T., Shao-Horn, Y., Gomez-Bombarelli, R., Ji, S., Gross-\nman, J.C.: Examining graph neural networks for crystal structures: limitations\nandopportunitiesforcapturingperiodicity.ScienceAdvances9(45),3245(2023)\n[28] Choudhary, K., Sumpter, B.G.:Can a deep-learning modelmakefast predictions\nof vacancy formation in diverse materials? AIP Advances 13(9) (2023)\n[29] Zhu,R.,Nong,W.,Yamazaki,S.,Hippalgaonkar,K.:Wycryst:Wyckoffinorganic\ncrystal generator framework. Matter 7(10), 3469–3488 (2024)\n[30] Kazeev, N., Nong, W., Romanov, I., Zhu, R., Ustyuzhanin, A., Yamazaki, S.,\nHippalgaonkar,K.:Wyckofftransformer:Generationofsymmetriccrystals.arXiv\npreprint arXiv:2503.02407 (2025)\n[31] Levy, D., Panigrahi, S.S., Kaba, S.-O., Zhu, Q., Lee, K.L.K., Galkin, M., Miret,\nS., Ravanbakhsh, S.: Symmcd: Symmetry-preserving crystal generation with\ndiffusion models. arXiv preprint arXiv:2502.03638 (2025)\n[32] Magar,R.,Wang,Y.,BaratiFarimani,A.:Crystaltwins:self-supervisedlearning\nfor crystalline material property prediction. npj Computational Materials 8(1),\n231 (2022)\n[33] Wang, Y., Cho, K.: Non-convolutional graph neural networks. arXiv preprint\narXiv:2408.00165 (2024)\n[34] Argonne Leadership Computing Facility: Aurora. https://www.alcf.anl.gov/\naurora\n[35] Yan, K., Liu, Y., Lin, Y., Ji, S.: Periodic graph transformers for crystal mate-\nrialpropertyprediction.AdvancesinNeuralInformationProcessingSystems35,\n15066–15080 (2022)\n[36] Chen, C., Ong, S.P.: A universal graph deep learning interatomic potential for\nthe periodic table. Nature Computational Science 2(11), 718–728 (2022)\n[37] Simeon, G., De Fabritiis, G.: Tensornet: Cartesian tensor representations for effi-\ncientlearningofmolecularpotentials.AdvancesinNeuralInformationProcessing\nSystems 36, 37334–37353 (2023)\n[38] Schu¨tt,K.T.,Hessmann,S.S.,Gebauer,N.W.,Lederer,J.,Gastegger,M.:Schnet-\npack 2.0: A neural network toolbox for atomistic machine learning. The Journal\nof Chemical Physics 158(14) (2023)\n[39] Choudhary, K., Garrity, K.F., Reid, A.C., DeCost, B., Biacchi, A.J.,\nHightWalker,A.R.,Trautt,Z.,Hattrick-Simpers,J.,Kusne,A.G.,Centrone,A.,\net al.: The joint automated repository for various integrated simulations (jarvis)\nfor data-driven materials design. npj computational materials 6(1), 173 (2020)\n[40] Ko, T.W., Deng, B., Nassar, M., Barroso-Luque, L., Liu, R., Qi, J., Thakur,\nA.C.,Mishra,A.R.,Liu,E.,Ceder,G.,et al.:Materialsgraphlibrary(matgl),an\nopen-source graph deep learning library for materials science and chemistry. npj\nComputational Materials 11(1), 253 (2025)\n[41] Feng,H.,Tian,H.:Improvingcrystalpropertypredictionfromamultiplexgraph\nperspective. Journal of Chemical Information and Modeling 64(19), 7376–7385\n(2024)\n[42] Gilmer, J., Schoenholz, S.S., Riley, P.F., Vinyals, O., Dahl, G.E.: Neural mes-\nsage passing for quantum chemistry. In: International Conference on Machine\nLearning, pp. 1263–1272 (2017). PMLR\n[43] Goodall, R.E., Lee, A.A.: Predicting materials properties without crystal struc-\nture: Deep representation learning from stoichiometry. Nature communications\n11(1), 6280 (2020)\n[44] Ruff, R., Reiser, P., Stu¨hmer, J., Friederich, P.: Connectivity optimized nested\ngraph networks for crystal structures. arXiv preprint arXiv:2302.14102 (2023)\n[45] New, A., Le, N.Q., Pekala, M.J., Stiles, C.D.: Self-supervised learning for crystal\nproperty prediction via denoising. arXiv preprint arXiv:2408.17255 (2024)\n38"}, {"round": 2, "reviewer_comments": {"reviewer_1": "The authors have addressed my comments properly, and I suggest acceptance of this paper to be published in Nature Communications.", "reviewer_2": "The authors did a good job in addressing my comments. I don't have futher questions", "reviewer_3": "I thank the authors for including additional model comparisons and validation analyses, which help support the model's performance and are central to this study's motivation. That said, my earlier concern remains: the proposed CLOUD method does not consistently outperform PMCGNN across key metrics. According to Table 6, only Ehull (reported in eV) shows a strong result, with other measures lacking clear improvement. To strengthen the case for CLOUD, I encourage the authors to highlight its distinct advantages, whether in terms of scalability, robustness, or applicability in specific scenarios, to better demonstrate its unique value beyond the current benchmarks."}, "author_response": "Response to Reviews on Article “CLOUD: A\nScalable and Physics-Informed Foundation Model\nfor Crystal Representation Learning”\nChangwen Xu1, ShangZhu1, VenkatasubramanianViswanathan1,2*\n1Department of Mechanical Engineering, University of Michigan.\n2Department of Aerospace Engineering, University of Michigan.\n*Corresponding author(s). E-mail(s): venkvis@umich.edu;\nContributing authors: changwex@umich.edu; shangzhu@umich.edu;\nManuscript ID: NCOMMS-25-45609\nWe are very grateful for the reviewers’ comments, suggestions, and valuable input on\nourmanuscript.Wearepleasedthatthepreviousroundofrevisionshasaddressedmost\noftheconcerns.Inwhatfollows,wewouldliketoanswertheremainingquestions.We\nstronglybelievethattherevisedmanuscriptshouldmeetthestandardsforpublication\nin Nature Communications.\nComment: The authors have addressed my comments properly, and I suggest\nacceptance of this paper to be published in Nature Communications.\nAuthor Reply: We appreciate that our previous response has addressed the\nreviewer’s comments. We also sincerely thank the reviewer for the constructive feed-\nback, which has helped improve the quality and clarity of our manuscript. We further\nthank the reviewer for recommending acceptance of our paper for publication in\nNature Communications.\nReviewer #2 Comments\nComment: The authors did a good job in addressing my comments. I don’t have\nfuther questions.\nAuthor Reply: We are glad that our previous response has satisfactorily addressed\nthe reviewer’s concerns. We are grateful to the reviewer for the constructive sugges-\ntions, which have strengthened the clarity and overall quality of our manuscript.\nReviewer #3 Comments\nComment: I thank the authors for including additional model comparisons and val-\nidation analyses, which help support the model’s performance and are central to this\nstudy’s motivation. That said, my earlier concern remains: the proposed CLOUD\nmethod does not consistently outperform PMCGNN across key metrics. According to\nTable 6, only Ehull (reported in eV) shows a strong result, with other measures lack-\ningclearimprovement.TostrengthenthecaseforCLOUD,Iencouragetheauthorsto\n2\nhighlight its distinct advantages, whether in terms of scalability, robustness, or appli-\ncability in specific scenarios, to better demonstrate its unique value beyond the current\nbenchmarks.\nAuthor Reply:Wesincerelythankthereviewerforthethoughtfulcommentandfor\nrecognizing that the additional analyses we provided strengthen the paper’s motiva-\ntion. We agree that CLOUD does not consistently surpass PMCGNN across all five\nJARVIS tasks. PMCGNN is a strong state-of-the-art GNN baseline, and our results\nshowthatwhileCLOUDdoesnotdominateituniformly,CLOUDconsistentlyoutper-\nforms CGCNN and achieves performance comparable to other widely used baselines\nsuch as ALIGNN across these benchmarks.\nAt the same time, we appreciate the reviewer’s suggestion to highlight CLOUD’s\ndistinct advantages beyond benchmark accuracy. We would like to clarify that the\ncore goal of this work is not merely to introduce another benchmark competitor, but\nrather to demonstrate how a scalable scientific foundation model can bridge machine\nlearninganddomainphysicsformaterialsdiscovery.Therefore,therevisedmanuscript\nalready emphasizes several unique strengths of CLOUD:\n• Scalability. In the section “Transferable and scalable crystal representation learn-\ning with CLOUD” in Results, we report scaling analyses showing that CLOUD\nmaintains strong performance with respect to dataset size and model size. In our\nprior response letter for the last round of revision, we also provided scaling plots\nfrom HPC cluster runs illustrating CLOUD’s favorable scaling behavior.\n• Inference efficiency.AsshowninTable5ofourpreviousresponseletter,CLOUD\nachieves faster inference than graph neural networks such as ALIGNN, Matformer,\nMEGNet, M3GNet, TensorNet, and SO3Net—even though CLOUD contains more\nparameters. This advantage is central for large-scale screening workflows.\n• Integration with differentiable physics. A major contribution of this work is\nCLOUD-DEBYE,whichintegratesCLOUDwiththeDebyemodelinadifferentiable\nframework.Withtrainingdataofheatcapacityonlyunder300K,CLOUD-DEBYE\npredictstemperature-dependentC (T)accuratelyfordifferenttemperaturesandin\nv\na physics-consistent manner. These results are detailed in the section “CLOUD-\nDEBYE:IntegrationofCLOUDinadifferentiableandphysics-consistentmaterials\nthermodynamics framework” in Results.\n• Practical materials discovery potential. In the last revision, we expanded the\nevaluation to WBM dataset from MatBench Discovery, adding convex-hull–based\nscreening metrics and a case study of Al Fe. These results highlight CLOUD’s\n2\napplicability to real discovery workflows beyond standard supervised benchmarks.\nbility, inference speed, physics integration, and real-world applicability, even if it does\nnot uniformly outperform PMCGNN on all supervised metrics. We appreciate the\nreviewer’s guidance in clarifying CLOUD’s unique values.\n4"}]

Decision letter