High-order Equivariant Flow Matching for Density Functional Theory Hamiltonian Prediction
Authors: Seongsu Kim
Nayoung Kim
Dongwoo Kim
Sungsoo Ahn
Venue: NeurIPS 2025 spotlight ()
Keywords: AI for sciences
Machine learning for science
Hamiltonian prediction
Density functional theory
Quantum physics
Physics
Computational chemsitry
Flow matching
Generative model
Primary area: machine_learning_for_sciences
TL;DR: We propose QHFlow, an SE(3)-equivariant flow model for generating Hamiltonians.

Abstract: Density functional theory (DFT) is a fundamental method for simulating quantum chemical properties, but it remains expensive due to the iterative self-consistent field (SCF) process required to solve the Kohn–Sham equations. 
Recently, deep learning methods are gaining attention as a way to bypass this step by directly predicting the Hamiltonian.
However, they rely on deterministic regression and do not consider the highly structured nature of Hamiltonians.
In this work, we propose QHFlow, a high-order equivariant flow matching framework that generates Hamiltonian matrices conditioned on molecular geometry. Flow matching models continuous-time trajectories between simple priors and complex targets, learning the structured distributions over Hamiltonians instead of direct regression. To further incorporate symmetry, we use a neural architecture that predicts SE(3)-equivariant vector fields, improving accuracy and generalization across diverse geometries. To further enhance physical fidelity, we additionally introduce a fine-tuning scheme to align predicted orbital energies with the target. QHFlow achieves state-of-the-art performance, reducing Hamiltonian error by 71% on MD17 and 53% on QH9. Moreover, we further show that QHFlow accelerates the DFT process without trading off the solution quality when initializing SCF iterations with the predicted Hamiltonian, significantly reducing the number of iterations and runtime.

High-order Equivariant Flow Matching for
Density Functional Theory Hamiltonian Prediction
Seongsu Kim!
Nayoung Kim! KAIST nayoungkim, sungsoo.. ahn) ac.kr kr
Dongwoo Kim²
Sungsoo Ahn
POOSTECH
{seongsu. kim,
dongwoo.. kin@postech ac kr
Abstract
Density functional theory (DFT) is a fundamental method for simulating quantum chemical properties, but remains expensive due the iterative self-consistent field (SCF) process required to solve the Kohn-Sham equations. Recently, deep learning methods are gaining attention as way to bypass this step by directly predicting the Hamiltonian. However, they rely on deterministic regression and do consider the highly structured nature of Hamiltonians. In this work, we propose QHFLOW, high-order equivariant flow matching framework that gen- erates Hamiltonian matrices conditioned on molecular geometry. Flow matching models continuous-tin time trajectories between simple priors and complex targets, learning the structured distributions over Hamiltonians instead of direct regres- sion. To further incorporate symmetry, we use a neural architecture that predicts SB(3)-equivariant vector fields, improving accuracy and generalization across di- verse geometries. To further enhance physical fidelity, we additionally introduce a fine-tuning scheme to align predicted orbital energies with the target. QHFLOW achieves satte-of-the-att performance, reducing Hamiltonian eπror by 73% on MD17 and 53% on QH9 compared the previous best model. Moreover, we further show that QHFLOW accelerates the DFT process without trading off the solution quality when initializing SCF iterations with the predicted Hamiltonian, significantly reducing the number of iterations and runtime.
Introduction
Density functional theory (DFT) [1, 2] is a cornerstone of modern computational physics [3, 4], chemistry [5, 6), and materials science [7, 8). It offers a powerful trade-off between accuracy and computational efficiency in predicting the electronic properties of atoms, molecules, and solids [9- 11]. The practical implementation of DFT relies mainly on solving the Kohn-Shann Shan- equations [12]. which describe system of non non-interacting electrons that yield the same ground-state density as the interacting system [12]. These equations are typically solved using the self-consistent field (SCF) method [13-15]. iterative procedure that updates the electron density until convergence. To this end, at each iteration, the Kohn-Shane Hamiltonian must be reconstructed from the kinetic operator, external potential, Coulomb (Hartree) potential, and exchange sichange-correlation potential [16, 17]. However, this reconstruction of the Hamiltonian becomes computationally demanding for larger systems, limiting the scalability of conventional DFT approaches [18]. To address this challenge, recent work has explored machine learning models to predict the Hamil- tonian directly from atomic configurations, with the aim of bypassing accelerating the SCF loop [19-26]. A key challenge in this task is ensuring that the model preserves the physical symme- tries of molecular systems, particularly SE(3) symmetry, which governs how Hamiltonians transform
39th Conference on Neural Information Processing Systems (NcurIPS 2025).

---

(8) (6)
Hg Milt
H1 Pi
Density PM
M
(h)
CLUMO CHOMO
HC
SCe
SCF / DFT Acceleration
Figure 1: Overview of QHFLOW. (a) Initial Hamiltonian sampled from the tensor expansion-based SE(3)-invariant prior (TE). (b) Initial Hamiltonian sampled from the Gaussian orthogonal ensemble SE(3)-invariant prior (GOE). (c) QHFLOW transforms the initial Hamiltonian Ho into the target Hamiltonian H1 using flow matching, guided by an -quiivaqura vector field uo (t, from an invaraint prior Po- The predicted Hamiltonian H₃ defines learned target distribution p and is used to compute the electronic density P.M as well as the ÉLUMO and еномо. When H₃ is used to initialize SCF, QHFLOW can accelerate conventional DFT and SCF procedures.
under rotations and translations due their construction in spherical harmonics basis. SE(3)- cquivariant networks such as PhisNet [20], QHNet [23], WANet [25], and SHNet [26] have been proposed address this, using tensor representation and Clebsch-Gordan products enforce rota- tional equivariance. These approaches rely on pointwise regression, which may lead to high predictive uncertainty and often struggle capture the structural correlations present in the Hamiltonian matrix.
Contribution. this work, we introduce QHFLOW, high-order SE(3)- for generating Kohn-Sham Hamiltonians. Our work uses flow matching learn neural ordinary equation (ODE) that transforms a prior distribution to the target distribution of Hamiltonians. We also parameterize the vector fields using high-order -Euui-equi neural networks ensure that symmetry is preserved not only in the output, but throughout the entire ODE trajectory [27, 28).
flow matching
To be specific, QHFLOW generates Hamiltonians conditioned on molecular geometries using equiv- ariant architectures and trained match the distribution of Hamiltonians obtained from DFT. To this end, we design two types of (E33)-invaiiant priors to support equivariant flow-based learning: a Gaussian orthogonal ensemble (GOE) and a tensor expansion-ber basonnd (TE) prior that includes a introduce
group theorrtic structure. To further enhance the accuracy of energy-related properties, a fine-tuning stage that aligns the predicted orbital energies with those from the DFT, inspired by the weight alignment loss in WANet [25]. We provide overview of QHFLOW in Figure 1.
Our main contributions are as follows:
Flow matching for Hamiltonian prediction: We are the first to formulate Hamiltonian pre- diction as a generative problem, learning trajectory to generate Kohn-Sham Hamiltonians conditioned on molecular geometry using SE(3)-equivariant vector fields.
Symmetry-aware prior distributions: We introduce two SE(3)-invariant priors, based on the Gaussian orthogonal ensemble (GOE) and the tensor expansion (TE), ensuring equivariant flow matching and symmetrical initialization of Hamiltonian priors.
Energy alignment fine-tuning: We introduce fine-tuning strategy that aligns the predicted orbital energies with target orbital energies. This step encourages the model match the target orbital energies, resulting in physically consistent property predictions.
Empirical performance: Extensive experiments on MD17 and QH9 demonstrate that QHFLOw outperforms recent methods in terms of predicting the Hamiltonian MAE, orbital energy prediction, HOMO, LUMO and gap energy accuracy. We also show that using QH- FLOw's predictions initialization for the SCF procedure leads to significant acceleration of the DFT process, reducing the number of iterations and total computation time.
2 Related work
Hamiltonian matrix prediction. Learning-based approaches for Hamiltonian matrix prediction have emerged promising alternatives the computationally demanding traditional DFT method.

---

Early work such as SchNOrb [19] extends SchNet [29] to enable accurate prediction of molecular orbitals. PhiSNet [20] introduces an architecture that explicitly respects Sym3)-yymmetr using tensor representations and tensor product operations. QHNer [23] improves computational efficiency by reducing the number of high-order tensor products, and WANet [25] further extends this direction to larger molecular systems by incorporating physically grounded loss function. Recently, SPHNet [26] improves scalability through adaptive path selection in the tensor product. Despite these advances, all preceding methods treat Hamiltonian prediction approach leverages flow matching to model the structural correlation in the Hamiltonian.
pointwise regression task. In contrast, our
Flow matching for prediction. Flow matching [30, 31] has gained significant attention in generative modeling for ability to leverage an ordinary differential equation (ODE)-based formulation, enabling faster inference than diffusion models [27]. Flow matching has been increasingly applied to machine learning problems chemistry and drug discovery [32]. particular, recent work has shown its effectiveness in regression-style applications such as protein structure prediction [33], crystal structure prediction [34, 35], and protein binding site identification [36], where modeling complex target distributions is crucial. In this work, we further extend flow matching Hamiltonian matrix prediction. By framing the task as a distributional matching problem, our model learns structured outputs with greater accuracy and physical fidelity than conventional regression-based methods, illustrating the broader utility of flow matching in scientific prediction tasks.
3 Preliminary
3.1 Kohn-Sham density functional theory and Roothaan-Hall equation
Kohn-Sham density functional theory (KS-DFT). The Kohn-Sham (KS) equations [1, 12] are sim- plified analogs of the many-body Schrödinger equation [37] that retain a one-to-one correspondence with the ground-state density of the interacting electron system. This framework has revolutionized quantum mechanical simulations of electronic structures in atoms, molecules, and solids [9-11], offering a practical balance between computational efficiency and predictive accuracy.
Consider N N-particle system associated with the wavefunetion V 3 -> C. In principle, the system is governed by the stationary Schrödinger equation HV = cV, which describes the quantum behavior of a system using the Hamiltonian operator H and the total energy A direct solution is generally intractable, as its computational cost scales exponentially with the number of particles N.
John-Sham density functional theory (KS-DFT) provides a computationally feasible alternative. Instead of solving the complex many-body wavefunction V KS-DFT uses the ground-state electron density p(r) its fundamental variable. This is achieved by introducing a fictitious system of non- interacting electrons that share the same ground-state density with the real, interacting wavefunction V. The behavior of non-interacting system is described by a set of N single-particle orbitals, {tha J.V which are the solutions to the KS equations:
N p(r) = Σ (A))))
H[p 1b E.bss
where p: 3 -> IR is the electron density, H[p is the effective single-particle KS Hamiltonian constructed from p, and e E RR is the orbital energy. The H[p] contains the external, Hartree, and exchange-correlation terms. We provide detailed explanation in Appendix A.
Roothann- (RH) equation. practice, the KS-equations are projected into a finite basis set, converting differential equations into the matrix eigenvalue problem known as RH-equation [38, 39]:
H[C]C = SCe.
(2)
Here, C E IR BxN denotes the set of coefficients where each row {Cm) )b-1 express an orbital th as a linear combination of B basis functions ( ()) = Σ Os Next, H[C) is the Hamiltonian matrix H E IR Bx B constructed by the coefficient C, and e E RN* is the diagonal matrix of orbital energies i.e., diag(e1, eN). Finally, S E BxB the matrix of overlap between basis functions Saa/ f d (r) (r)dr where of denotes the complex conjugate of de. The basis functions (r)) is selected to match the system geometry and commonly chosen between Slater-type orbitals (STOs) [40, 15]. Gaussian-type -ype orbitals (GTOs) [41], plane waves [42, 10].
=
3

---

3.2 -quuevauiiaee in RH-DFT
-qqiivanaae and irreducible representations. Formally, map f is said to be SO(3)- equivariant if R f(z) (R x) for all R SO(3) where denotes the SO(3) group action. This equivariance can be described in terms of the representation theory [43]. A representation of SO(3) is a matrix-valued function that describes the transforms under rotations, and it is called irreducible if cannot be decomposed into smaller invariant subspaces. These irreducible representations (irreps) form the fundamental components for constructing eqqiivaraaa function, and irrep vector is the vector that transforms under a specific irrep. Each SO(3) irrep is indexed by non-negative integer 1, and represented by a Wigner D-matrix D (R) C(24 (2x(2)+)) These matrices characterize how scalar (e 0), vector (E = 1), and higher-order tensor features (E 2 2) transform under rotation. A rank-& irrep vector w(t) с (2(++1) transforms under rotation R as: w R R w(() D( (R.w w
where R implies the transformation under R. We provide additional background on group theory and equivariance in Appendices A.2 and A.3.
High-order SO(3)--quivariane of RH-DFT. The Hamiltonian matrix in RH-DFT exhibits highly structured form of equivariance. As shown in Figure 2, the matrix can be organized into blocks, where each block corresponds to interactions between two sets of orbital basis functions indexed by pairs of quantum prin- cipal numbers and angular momentum num- bers, i.e., (71, (1) and (m2, f2). For example, block labeled 1s-2p denotes interactions be- tween (n1, (1) = (1 and (n2, = (2, 1). These orbital groupings induce higher-order structure, where the Hamiltonian must remain equivariant not only at the global matrix level but also within each individual block. This makes Hamiltonian modeling more challeng- ing than conventional equivariant tasks such as force or energy prediction, as it requires explicit handling of block-level symmetry constraints.
1s
Zp
15
25
8
3p
8
Figure 2: Illustration of the Hamiltonian matrix structure of H,O. Rows and columns are indexed by quantum and angular momentum numbers (i.e., 1s, 2s) ordered by atom (O, H, H). The full matrix (right) partitioned into atomic blocks with gray dashed line; the top-left 9 x 9 sub-matrix corre- and is shown detail (left), grouped by quantum and angular momentum. Colors indi- cate the sign and magnitude matrix elements.
sponds
This symmetry arises naturally in RH-DFT when the basis set is constructed from spherical harmonics, which serve as irrep vectors of SO(3). In particular, standard basis functions such as STO or GTO are defined using spherical harmonics Y,m where é and denote the angular momentum and magnetic quantum numbers, respectively [44]. When using such spherical-type basis, the full orbital basis set {dr [ can be partitioned into K groups, with the k-th group {on 1 m=-f, sharing the same principal quantum number nt and angular momentum quantum number lg. Each group spans a (2fx 1) functions indexed by m, forming a rank-l, irrep vector. We assume the columns of the coefficient matrix C are ordered according to these groupings. Under a rotation R SO(3), the matrix C transforms as follows:
C R R C DDDRCC
D(R) = (R) (((x) (R),
(3)
where each ((+x) E (x244)) is the rank-&; Wigner D-matrix and e denotes the direct sum over irreps. The Hamiltonian matrix H. being a function of C, also satisfies the squiaeequi
B
R H
[D[R]C] = D(R)H (C)
(4)
We provide further details about the SO(3)-cquivariance structure of Hamiltonian submatrix and the direct sum of Wigner D-matrices in Appendix A.4.
4 Quantum Hamiltonian flow matching (QHFLOW)
4.1 Flow matching for Hamiltonian matrix
Problem formulation. Our goal is to train a model that predicts the Hamiltonian matrix H from molecular configuration M h). Here, for molecule M with M atoms, x E Mx3 and

---

h E RMxD denotes the atom-wise coordinates and D-dimensional atom-wise features, respectively. We train the model using a atomic dataset with a sample (M, H) corresponding to a molecular configuration M and its Hamiltonian matrix H obtained from conventional RH-DFT solvers.
Conditional flow matching (CFM). We first propose the flow matching for Hamiltonian matrix, by parameterizing our model as continuous normalizing flow (CNF) model. Our CNF aims to learn a time-dependent vector field vy (-IMM that outputs the Hamiltonian matrix H from solving the ODE d H, v, (H, (M) for time / (0, 1]. The ODE starts from H, sampled from the prior distribution Po defined over the space BBxB of Hamiltonian matrices and outputs a Dirac distribution centered at the Hamiltonian matrix H. To this end, we learn the vector field of U/ to match the distribution induced from the empirical distribution (M, H1) AA associated with conditional probability paths
Ps (H,HHH Importantly, the conditional path Pa (H [H) is expressed using closed-form conditional
vector field H-): dt H = u, (H, M) for Ho Po-
Following previous works [45, 35), we parameterize the conditional probability path using linear interpolation between Ho and Hj, resulting in the conditional vector field ue [H,[H, H- H. To train the CNF, we use the conditional flow matching objective defined as follows:
||| (H.)) (H.) u, (H, H
LCFM E.H.M))-A---00.... (.
(5)
where U/ denotes an uniform distribution and || practice, we parameterize the vector field by of results in the final loss function:
113 denotes the Frobeneous norm of matrix. In H (H. Hr using a neural network (-). This 1-t
LCrM = B.I.M))-
(H, HIM
(6)
(1 t1-
which one could interpret training H approximate the true Hamiltonian matrix H from its noisy versions H. We provide a full description of the training and sampling algorithms in Appendix B.1.
Fine-tuning for energy alignment. In practice, many important quantum properties derived from DFT rely directly on orbital energies. To improve the accuracy of such downstream quantities including orbital energies, highest occupied molecular orbital (HOMO) energy, loweset unoccupied molecular orbital (LUMO) energy, and the HOMO-LUME gap, we introduce a fine-tuning strategy that explicitly aligns the predicted orbital energies with their ground-truth DFT references. This approach is motivated by the weighted alignment loss (WALoss) introduced in WANet [25], which improves prediction quality by aligning the energy of predicted and reference matrix. Rather than incorporating this loss during full training, we adopt it as a post hoc fine-tuning objective.
Specifically, we define a fine-tuning objective that aligns the approximate orbital energies é derived from our predicted Hamiltonian matrix Hς from CNF with the ground-truth orbital energies e computed from RH-DFT as defined in Equation (2):
((H ((H |
LTT (...M))-A
(7)
Here, the ground-truth orbital energies e are computed via the identity e 6= HCC, which follows from the orthonormality condition C SCC I, where C is the ground-truth coefficient matrix and S is the overlap matrix. Following WALoss, we estimate é using the predicted Hamiltonian: C HfC s (Co) TH where C is the coefficient matrix obtained by solving the RH-DFT for the predicted Hamiltonian Hº. We provide the implementation details in Appendix B.2.
4.2 SE(3)) flow matching with equivariant priors.
Here, we introduce our approach to design our QHFLOW to satisfy -Euivvaquve We follow the framework of Song et al. [27] to learn an equivariant vector field under an invariant prior. To handle translational symmetry, our networks operate in a mean-free system by centering the molecular N
1 11 X given vector 1 E R Nx of ones. To make QHFLOW satisfy equivariance, we (1) parameterize the vector field of using a -quiaaqqivariaat neural network [23] and (2) propose Gaussian orthogonal ensemble (GOE) and tensor expansion- based (TE) priors for the Hamiltonian. The main challenge is to design new priors that are invariant to the rotation of the Hamiltonian matrix as described in Equation (4).
coordinates X at the origin, ie., set X + X
5

---

Gaussian orthogonal ensemble (GOE) prior. We here introduce simple rotational invariant prior, the GOE prior. A GOE is defined over matrices E R"x*, with each element drawn from Gaussian distribution: Hij (O0 o²). The log-density is proportional to the Frobenius norm of the matrix, which invariant under any transformations associated with SO(3) irreps, i.e, the 3ID-rotation matrix (f = 1) and even the high-rank Wigner D-matrices (E 2 2). Therefore, the prior satisfies the desired invariance, i.e., (H) = ((D(R) for all R SO(3). We provide the full proof for the 50(3)-invariance in Appendix C.I.
Tensor expansion-based (TE) prior. We also introduce TE prior, an -nvaa-invar distribution constructed by first sampling an irrep vector w(O), then applying tensor expansion to produce matrix ((4,62) E (2(6 +1)x((++) that transforms covariantly under rotation R E SO(3):
H (1. R, ( (R) (H ()))) (R)-1.
(8)
To this end, we first design an (0(3)-invarian) distribution of irrep vector p(w(fs) by decomposing the vector into radial and rotational parts: w(() rD(()) (R) where r is scalar norm and w(f)) is a fixed unit-norm irrep vector. i.e., Y for fixed direction fq. By sampling independently of
R. and choosing R uniformly on SO(3), the resulting prior p( (/) is SO(3)-invariant. In practice, we take r ~ Normal(p,o* o²) to ensure positivity, and R ~ (aiiform(SO33)))
Next, we apply tensor expansion, which maps rank-l irrep vector w(() (2(++1)) into a matrix of shape (2l1+) +1 x (262+11) +), indexed by the (m1, my) satisfying -61 mm S l and -b2 S my S Eg.
respectively. Each entry of the resulting matrix is defined as follows:
zww()
(B1,ta)
Σ
cle,m) ((2.m)))(( w
where C' denotes Clebsch Corsdan coefficients, @ denotes the tensor expansion operation. Note that the superscript (41, E2) indicates the output irrep structure and can be omitted when clear from context. The rotation of this expanded matrix is represented by the Wigner D-matrix as follows:
(ow RS, ((4)) (R) (xw(()) ( ((2)) (R)-1.
(10)
See Appendix C.2 for proof. This property allows us to construct SO(3)- inverriiva distributions over Hamiltonians. We formalize this in the following theorem:
Theorem 1. Let (ow(1) (Er,f)) where an irrep vector w(() p(w(()) is drawn from SO(3)-
invariant distribution. Then the induced distribution over H is invariant under SO(3) transformation:
p(H)= ((D(t)) (R)) (R)
VR E SO(3).
(11)
We provide proof in Appendix C.3. Using the expansion, we can model the invariant distributions for each sub-matrix of Hamiltonian H. We further extend this construction to build global invariant distributions over full Hamiltonian matrices H by composing multiple high-order irrep vector via tensor expansion. We provide a detailed construction of this multi-bloc block formulation in Appendix C.4.
4.3 Model implementation
We construct QHFLOW by extending QHNet [23], an SE(3)- GNN for Hamiltonian predic- tion. Our architecture introduces time conditioning and two additional inputs: the current Hamiltonian H, and the overlap matrix S. The model takes (M, H, S, t) and predicts the target Hamiltonian H as described in Equation (6). We provided the implementation details in Appendix D.I.
Feature intialization. We first construct the atom-wise features h; S; based on (1) extracting atom- wise submatrices H,, S, from H,, S and (2) mapping the submatrices into SO(3) irrep vectors based on Wigner-Eckarr projections [46], i.e., set h, = ProjBlock (H,). S; ProjBlock (S:) where ProjBlock flattens and applies change of basis transformations to the input matrix. We further encode time / through tensor field network (TFN) [47] layer and sinusoidal encoding e [i)m i.c., set m, TFN (x,, ez)
Message-passinne layers. To propagate information between neighbors, our model updates the irrep vectors h; and s; using -quuivaquv attention-based layers adapted from Equiformer [48]:
h +- EquiBlock( hi: {(bj, diy)) ) Et). <- EquiBlock {(S), dry se.)).)

---

where N2 denotes the neighbors of atom i and dı is distance between atom i and j. We then apply channel mixing module. Mix, which linearly combines the of irrep vectors by flattening and reordering them according to their corresponding irrep indices while preserving -quiivaaaie This module update the m; the residual manner:
f TFN m. + Mix(h, 8,, m.), , {myy diy se).
(12)
This update process is repeated for fixed number of iterations, allowing the model to integrate geometric and physical information across the molecular graph.
Hamiltonian reconstruction. The atom-wise irrep vectors m, are passed into two TFN modules to predict Hamiltonian components. The first TFN, computes atom-wist wise irreps wa using messages from neighboring atoms, while the next TFNp predicts pairwise irreps Wij between atom i and j:
W5 TFN, (m; {(m,, keN)), Wij TFM (m, m.j des).
(13)
Each output consists of set of irrep vectors, Wii {w He,, and Wij {w ij JEEL where Lj is the set of angular momentum l corresponding to the atom i. Then Hamiltonian blocks are constructed via tensor expansion @ and learnable weights F:
(e)
H ((,(2)) = -E (2. (Ow!) (8) (6.
(14)
EEL.
where E L3 and l E Ly. and the case i = j corresponds to diagonal block. These blocks are then assembled into the full Hamiltonian matrix H. which is structured into K x K sub-matrices corresponding to all pairwise combinations of the K groups of basis functions. All operations in this reconstruction preserve -quuvvaeuve ensuring that both H and the learned vector field uf remain equivariant under any SB(3) transformation of the molecular geometry.
5 Experiments
In this section, we present experimental results that demonstrate the effectiveness of QHFLOW. We evaluate the performance on public benchmark and compare it with competitive baselines [19, 20, 23, 26] when metrics are available. We provide training settings and hyperparameters in Appendix
Datasets. We evaluate our model on two datasets: MD17 [49, 19] and QH9 [50, 51, 24J. MD17 contains DFT Hamiltonians along molecular dynamics trajectories for four small molecules, water (4,900 structures), ethanol (30,000), malondialdehyde (26,978), and uracil (30,000). The PBE exchange-correlation [52, 53] functional with a GTO basis set used for DFT.
QH9 consists of two subsets: stable and dynamic-3000 300k, each supporting in nn-distribution (id) and out- of of-distribution (ood) splits. The stable includes 130,831 molecules with up to 29 atoms. The id split randomly samples molecules, while the ood split groups them by size. The dynamic-300k contains 2,998 molecules, each with 100 perturbed geometries. The geo split assigns geometries randomly across splits for each molecule; the mol split partitions at the molecule level to test on unseen molecular identities. All QH9 Hamiltonians are computed using B3LYP [54] with the def2-SVP basis set. We provide the details about the dataset in Appendix E.I.
Evaluation metrics. Our evaluation metrics include the mean absolute error (MAD) of the Hamilto- nian H. the MAE of occupied orbital energies (Coce), and the similarity score of the orbital coefficient matrix (S.). To assess physical fidelity QH9, we also report the energy of LUMO (ELUMO). HOMO (CHOMO), and HOMO-LUME gap (EA). We provide details of the metrics in Appendix E.2.
5.1 Performance on MD17 dataset
We extensively evaluate QHFLOW on the widely used MD17 benchmark dataset to assess its accuracy in predicting Hamiltonian matrices in diverse molecular geometries. As shown in Table 1, we compare QHFlow with baseline models in four representative molecules, and it consistently outperformed all existing methods. With up to 73% reduction Hamiltonian prediction error, QHFlow demonstrates strong predictive performance. This improvement highlights its ability to capture fine-grained variations in quantum Hamiltonians resulting from subtle geometric changes, an essential capability for accurately modeling physical and chemical behavior.
7

---

Table MD17 benchmark. Results shown in bold denote the best result in each column, whereas those that are underlined indicate the second best.
Water (3 atoms)
Fihanol (9 atoms)
Malondialdehydee (9 atoms)
Uracil (12 atoms)
Model
[pEN.]
H
Cacc [pEx]]
S.t
[nEs]]
H
Coer [μEA]
S."
[pEx]
H 4
Fuse [pEa]]
[%]
S²1
[μEn]]
H L
feack [μEn.]
S."
[%)
SchNOrb
165.4
279.3
100.00
187.4
334.4
100.00
191.1
400.6
99.00
227.8
1760.
90.00
PhiSNet
15.67
85.53
100.00
20.09
102.04
99.81
21.31
100.6
99.89
18.65
143.36
99.86
QHNeet'
11.70
26.06
100,00
27.99
99.33
99.99
29.60
100. 16
99.92
26.80
127.93
99.87
SPHNet
23.18
182.29
100.00
21.02
82.30
100.00
20.67
95.77
99.99
19.36
118.21
99,99
Ours
4.93
100.00
5.33
29.03
100.00
3.80
22.68
99.99
3.68
30.54
Table QH9 benchmark. WA-FT implies finetuned model with WAloss. Results shown in bold denote the best result in each column, whereas those that are underlined indicate the second best.
Dataset
Model
H 4 [aEa]
[pEAa]
S. f [%)
4L1MO
4 [1Eba]
CHOMO
[p.Ell
6a [nEx]
QHNer*
77.72
963.45
94.80
18257.34
1546.27
17822.62
WANet
80.00
833.62
96.86
stt-le (id)
SPHNet
45.48
334.28
97.75
Ours
22.95
119.67
99.51
437.96
179.48 48
553.87
Ours (WA-FT)
23.85
101.92
187.48
92.22
QHNer*
69.69
884.97
93.01
25848.83
1045,99
25370.10
Qtt-se
SPHNet
43.33
186.40
98.16
(ood)
Ours
20.01
84.54
99.04
321.20
130.74
395.83
Ours (WA-FT)
20.55
72.64
99.16
171.24
77.96
179.57
QHNer
88.36
1170.50
93.65
23269.41
2040.06
22407.96
WANet
74.74
416.57
99.68
QH9-dynamic
SPHNet
52.18
100.88
99.12
(300k-gea)
Ours
25.94
103.11
99.59
425
175.18
547.33
Ours (WA-FT)
27.12
89.03
99.65
136.63
84.17
154.68
QHNec*
121.39
5554.36
86.02
$3505.09
4352.76
50424.86
QH9-dynamic
SPHNet
108.19
1724.10
91.49
(300k-mo1)
Ours
45.91
442.56
98
1344.68
479.71
1605.03
Ours (WA-FT)
46.60
424.75
98.74
912.10
403.51
1047.88
5.2 Performance on QH9 dataset
We evaluate QHFLOW on the QH9 dataset, a more challenging benchmark that includes various molecular sizes and compositions. This setting requires the model to generalize effectively across chemical and geometric variations. We trained on four different data splits. As shown in Table QHFLOW consistently outperforms all prior methods in all metrics. In particular, it shows significant improvements properties sensitive to eigenvalues such as éLUMO, éHоMO, and εΛ, showing the benefits of the flow-based model in capturing physically meaningful structure in Hamiltonians. These results demonstrate QHFLOW's 's strong potential as a surrogate model for DFT.
We also apply our energy fine-tuning stage to pretrained QHFLOW using the weighted alignment loss (WAloss) [25] for an additional 60,000 steps, which we denote WA-FT, As shown Table this improves predictions of orbital energies and their downstream task, confirming that energy alignment fine-tuning serves an effective inductive bias. We compare training from scratch with the WALoss and flow matching objectives in Appendix F.I.
5.3 DFT acceleration performance
Here, we show that our QHFLOW can accelerate DFT through initializing the SCF iterations [13-15] using the predicted Hamiltonian matrix, replacing the conventional SCF initialization methods.
Preliminaries on SCF method. Starting from an initial Hamiltonian H each SCF iteration solves
the generalized eigenvalue problem ( (()) (C()) (A) to obtain the coefficients ((A) The coef- ficients are then used to construct the density ((t) and subsequently update the Hamiltonian to H (k+1)
through the mapping (A) 2, (kk) H. ((+1),, where & and H denote the density construction and Hamiltonian update operators, respectively. The iteration continues until convergence.
WWe use re- re-trained version of QHNet, based the publicly available implementation.

---

019-stable (55
0119 stable (ood)
52.4
52.8
40.0
31.7 31.3
30.9
0,5 0.8
9,4 0.7
SCT Iter Rarlo
SC7 Retio
In Retio
Total T Ratio
SCT Her Ratio
SCT Retio
Irr Rebio
Total Retio
Qr9-dynamic (300k-gao)
QH0-dynamic (300kkmmll)
52.7. 65.3 45.7
54. 46.1 40.6
60
54.5 45.4 46.0
25.4
47.1 46.7
51.5
30.9 31.2
0.5 0.8 0.9
0.4 0.7 In T Racio
SCF Ber Ratio
SCP Ratio
In Rabio
Tdtal T Ratio
SCF ftar Raiio
SCF T Ratio
Tocal Rallo
OHNect
Ours
Ours (WAAFFT)
Figure 3: DFT acceleration performance on 300 samples from the QH9 dataset. All metrics are reported percentages relative conventional DFT (initialized with minao), which serves as the 100% baseline. The SCF Iter Ratio measures the ratio of SCF iterations required, while Inf T Ratio, SCF T Ratio, and Total T Ratio measure time. Lower SCF Iter Ratio and Total Ratio values indicate faster convergence. For example, a Total T Ratio of 46% means QHFLOW converges in 46% of the conventional DFT time, including the negligible model inference time.
DFT acceleration via SCF initialization. We use QHFLOW's prediction H to replace the con- ventional Hamiltonian guess by beginning the SCF loop closer to convergence. Performance is summarized with four relative metrics normalized by the conventional DFT baseline metric: SCF Iter Ratio: ratio of SCF iteration count, SCF Ratio: ratio of SCF wall time, Inf Ratio: ratio of inference time as a fraction of baseline SCF wall time, Total Ratio: Ratio of net runtime after adding inference overhead. We provide the definition of metric Appendix E.2.
We evaluate QHFLOW on 300 test molecules from the QH9 dataset, selected via fixed index slice using a 3-step ODE-based sampling. As shown in Figure 3, QHFLOW significantly accelerates the convergence of SCF. For instance, on the QH9 id benchmark, QHFLOW achieves convergence in 46% of the runtime required by conventional DFT, corresponding 54% relative speedup. Compared to previous ML-based initializations methods like QHNet which requires 53% of the conventional runtime, QHFLOW provides an additional reduction of 7% points. These results demonstrate the practical utility of QHFlow for accelerating real-world DFT simulations.
5.4 Ablation studies
Here, ablate design choices, with more results (e.g., orbital energy metrics) in Appendix F.3.
GOE, TE-based priors. We explore the effects of the flow matching prior for Hamiltonian prediction by comparing the
Table 3: Effect of prior distribution.
two priors introduced Section 4: GOE and TE prior. The GOE prior treats each symmetric matrix structure-free noise, whereas the TE prior injects group-theoretic bias that matches the blockwise symmetry of Hamiltonians.
Data
Pršor
[p
faxe [rF))
S.1 [%6]
GOE
25.93
151.65
99.39
22.95
119.61
99.51
GOE
20.4
87.32
98.55
add
TE
20.01
84.54
99.04
Table compares GOE and TE priors, where the TE prior results are identical to those of QHFlow reported in Table 2. We found that the TE prior consistently yields lower errors than the GOE prior across all splits. This highlights the importance of designing task-aligned and symmetry-avee aware priors for accurate Hamiltonian prediction.
geo
GOE
29.39
122.14
99.49
GOE
16.78
419.68
98.65
mol
TE
45.91
442.56
98.65
Predictive variance. To evaluate the robustness of QHFLOW stochastic initialization, we repeat infer- ence 5 times per molecule using the TE prior, chang- ing only the seed that generates the initial noise sam- ple. Results are shown in Table 4, with all metrics reported tet over five runs. The predicted Hamiltonians exhibit mean absolute deviation of just 0.03% and a standard deviation of predicted energy below 0, much smaller than the prevalent chem- ical accuracy thresholds. These results confirm that QHFLOW produces stable and consistent predictions despite stochastic initialization.
Table 4: Predictive variance of QHFLOW.
Data
Model
[pESs]
fwed [aEx]
S.t [%]
id
WA-FT
Ours
23.85=
22,95 noci
101.92:
119.67 szu
99.56=-3az
Ours
20.01 noc
84.54; son
99. 99:04::.00:
aod WA-FT
20.55g
72.64:
99.161
Ours
25.94 nucs
103 26 aes
99.59.-a00
geo
WA-FT
27.12
89.03
99.65
mol
Ours
46.60
45.91
443. 443.56G
424.75
98.74
98. 05.655

---

6 Conclusion
In this work, we introduce QHFLOW, high-order equivariant flow-based generative model designed to predict Hamiltonian matrices with high accuracy and physical consistency. By designing SE(3)- invariant priors GOE and TE, QHFLOW improves the accuracy and generalization across diverse geometries. Extensive experiments on the MD17 and QH9 datasets demonstrate that QHFLOW not only achieves state-of-the-aas performance across multiple evaluation metries but also significantly improves efficiency in downstream SCF calculations. Our results highlight the power of dlw-bbseed learning for scientific tasks and establish QHFLOW scalable and physically grounded framework for accelerating quantum chemistry simulations.
Limitation
We did not fully verify the generalizability of approach, e.g., extension to new architectures (WANet [25], SPHNet [26]) or datasets (PubChemQH [25]), due to the lack of resources, public codes, and public datasets. Our flow-based formulation introduces computational overhead for solving ODEs. In the future, one could consider removing the overhead with the one-step generative models [55]. Our experiments are limited to gas-phase molecules with up to 29 atoms and two commonly used XC functionals (PBE and B3LYP). The algorithms are left to be verified for more general settings, such as periodic solids or higher-accuracy functionals (e.g., hybrid or double-hybridy We also acknowledge that our datasets consist of relatively small molecules. One could further test the scalability of our approach once new datasets on larger systems are available. We also note that our work does not give a theoretical explanation of how Hamiltonian prediction error translates to SCF acceleration. Developing a theoretically grounded algorithm specifically for SCF acceleration would be interesting future research.
Broader Impact
This work focuses on advancing machine learning methods for quantum chemistry applications and does not involve human subjects, personal data, or sensitive information. As such, we believe that there are no direct ethical concerns associated with this research. Nevertheless, we acknowledge that improved computational tools for molecular modeling could have downstream applications in sensitive areas such as pharmaceuticals or materials development. We encourage responsible and ethical use of our methods in accordance with applicable laws and scientific guidelines.
Acknowledgements and Funding Disclosure
This work was partly supported by Institute for Information & communications Technology Planning & (valuation(IITPP) grant funded by the Korea (oVerN (RS- R0S-2011-11190005. Artificial In- telligence Graduate School Support (oogram(KAIST)) National Research Foundation of Korea(NRF) grant funded by the Ministry of Science and ICT(MSIT) (No. RS-2022-NR0721 Research Development Center) Cooperative Hub Program through the National Research Foundation of Korea(NRF) grant funded by the Ministry of Science and ICT(MSIT) (No. and the Institute of Information & Communications Technology Planning & Evaluation(IITP) grant funded by the Korea government(MSIT) thank Minsu Kim, Hyomin Kim, and Hyosoon Jang for helpful discussions and feedback.
GRDC(Global
Al Star Fellowship(KAISS We
References
[1] Pierre Hohenberg and Walter Kohn. Inhomogeneous electron gas. Physical review, 136(3B): B864, 1964. 1, 22
[2] Robert G Parr. Density functional theory of atoms and molecules. In Horizons of Quantum Chemistry: Proceedings the Third International Congress of Quantum Chemistry Held at Kyoto, Japan, October Noveemer 3, 1979, pages 5-15. Springer, 1989.
[3] Nathan Argaman and Guy Makov. Density functional theory: An introduction. American Journal of Physics, 68 68(1))69-77 2000.
10

---

[4] Axel D Becke. Perspective: Fifty years of -fsstty-funtiioall theory in chemical physics. The Journal of chemical physics, 140(18). 2014.
[5] Kieron Burke. Perspective on density functional theory, The Journal of chemical physics, 136 (15), 2012.
[6] Narbe Mardirossian and Martin Head-Gordon. Thirty years of density functional theory in computational chemistry: an overview and extensive assessment of 200 density functionals. Molecular physics, 115(19)2315-23322 2017. 1
[7] Jörg Neugebauer and Tilmann Hickel. Density functional theory in materials science. Wiley Interdisciplinary Reviews: Computational Molecular Science, 3(5):438-44) 2013.
[8] Andrew M Teale, Trygve Helgaker, Andreas Savin, Carlo Adamo, Bálint Aradi, Alexei V Arbuznikov, Paul W Ayers, Evert Jan Baerends, Vincenzo Barone, Patrizia Calaminici, et al. Dft exchange: sharing perspectives the workhorse of quantum chemistry and materials science. Physical chemistry chemical physics, 24(47):28700-287811 2022. 1
[9] Michael te Vrugt, Hartmut Löwen, and Raphael Wittkowski. Classical dynamical density functional theory: from fundamentals to applications. Advances Physics, 69(2):121-247, 2020. 1,3
[10] Richard M Martin. Electronic structure: basic theory and practical methods. Cambridge university press, 2020.
[11] Priyanka Makkar and Narendra Nath Ghosh. A review on the use of dft for the prediction of the properties nanomaterials. RSC advances, 1145):27877-27944, 2021. 3.3
[12] Walter Kohn and Lu Jeu Sham. Self-consistent equations including exchange and correlation effects. Physical review, 140(4A):A1133, 1965. 1, 3, 22
[13] JC Slater. A generalized self-consistent field method. Physical Review, 91(3):528, 1953. 1,
[14] Ehrenreich and Morrel H Cohen. Self-consistent field approach to the many-electron problem. Physical Review, 1 115(4))786, 1959.
[15] Warren J Hehre, Robert F Stewart, and John A Pople. Self-consistent molecular-orbite methods.
i. use gaussian expansions of slater-type atomic orbitals. The Journal of Chemical Physics, 1969. 1, 3, 8, 25
5 11(6):2657-266
[16] C David Sherrill. An introduction hartree-fock molecular orbital theory. School of Chemistry and Biochemistry Georgia Institute of Technology, 2000.
[17] Eva Pavarini, Erik Koch, Jeroen van den Brink, and George Sawatzky. Quantum materials: Experiments and theory. Modeling and Simulation, 6:420, 2016.
[18] Beatriz G del Rio, Brandon Phan, and Rampi Ramprasad. A deep learning framework to emulate density functional theory. npj Computational Materials, 911):18,, 2023.
[19] KT Schütt, M Gastegger, A Tkatchenko, K-R Müller, and RJ Maurer. Unifying machine learning and quantum chemistry with deep neural network for molecular wavefunctions. Nature Communications, 10(1), 2019. 1, 7,35
[20] Oliver Unke, Mihail Bogojeski, Michael Gastegger, Mario Geiger, Tess Smidt, and Klaus-Robon Müller. Se (3)-equivariant prediction of molecular wavefunctions and electronic densities. Advances in Neural Information Processing Systems, 34: 14434-14447, 2021. 2, 3, 3,
[21] He Li, Zun Wang, Nianlong Zou, Meng Ye, Runzhang Xu, Xiaoxun Gong, Wenhui Duan, and Yong Xu. Deep-learning density functional theory hamiltonian for efficient ab initio electromic-structuree calculation. Nature Computational Science, 2(6):3673377, 2022.
[22] Xiaoxun Gong, He Li, Nianlong Zou, Runzhang Xu, Wenhui Duan, and Yong Xu. General framework for (u)eequivar neural network representation of density functional theory hamiltonian. Nature Communications, 14(1):2848, 2023.

---

[23] Haiyang Yu, Zhao Xu, Xiaofeng Qian, Xiaoning Qian, and Shuiwang Ji. Efficient and equiv- ariant graph networks for predicting quantum hamiltonian. In International Conference on Machine Learning, pages 40412-40422 PMIR, 2023. 3, 5, 6, 7, 33, 35,
[24] Haiyang Yu, Meng Liu, Youzhi Luo, Alex Strasser, Xiaofeng Qian, Xiaoning Qian, and Shuiwang Ji. Qh9: A quantum hamiltonian prediction benchmark for qm9 molecules. Advances in Neural Information Processing Systems, 36:40487-40500 2023. 7, 35,36
[25] Yunyang Li, Zaishuo Xia, Lin Huang, Xinran Wei, Samuel Harshe, Han Yang, Erpai Luo, Zun Wang, Jia Zhang, Chang Liu, al. Enhancing the scalability and applicabilityy of kohn-sham hamiltonians for molecular systems. In The Thirteenth International Conference on Learning Representations, 2025. 2, 3, 5, 8, 10, 28, 33, 39
[26] Erpai Luo, Xinran Wei, Lin Huang, Yunyang Li, Han Yang, Zun Wang. Chang Liu, Zaishuo Xia, Jia Zhang, and Bin Shao. Efficient and scalable density functional theory hamiltonian prediction through adaptive sparsity. arXiv preprint arXiv:2502.0117., 2025. 1, 3, 7, 10, 33
[27] Yuxuan Song, Jingjing Gong, Minkai Xu, Ziyao Cao, Yanyan Lan, Stefano Ermon, Hao Zhou, and Wei-Ying Ma. Equivariant flow matching with bybrid probability transport for molecule generation. Advances in Neural Information Processing Systems, 36:549-568, 2023. 2, 3,5
[28] Leon Klein, Andreas Krämer, and Frank Noé. Equivariant flow matching. Advances in Neural Information Processing Systems, 36:59886 59910, 2023. 2
[29] Kristof T Schütt, Huziel E Sauceda, P-J Kindermans, Alexandre Tkatchenko, and K-R Müller. Schnet-a deep learning architecture for molecules and materials. The Journal of Chemical Physics, 148(,4), 2018. 3
[30] Yaron Lipman, Ricky TQ Chen, Heli Ben-Hamu Maximilian Nickel, and Matthew Le. Flow matching for generative modeling. In The Eleventh International Conference on Learning Representations, 2023.
[31] Alexander Tong, Kilian FATRAS, Nikolay Malkin, Guillaume Huguet, Yanlei Zhang, Jarrid Rector-Brooks, Guy Wolf, and Yoshua Bengio. Improving and generalizing flow-based genera- tive models with minibatch optimal transport. Transactions on Machine Learning Research, 2023.7
[32] Majdi Hassan, Nikhil Shenoy, Jungyoon Lee, Hannes Stark, Stephan Thaler, and Dominique Beaini. Et-flow: Equivariant flow-mater www-mathhing for molecular conformer generation. In The eight Annual Conference on Neural Information Processing Systems, 2025.
[33] Bowen Jing, Bonnie Berger, and Tommi Jaakkola. Alphafold meets flow matching for generating protein ensembles. arXiv preprint arXiv:2402.04845 2024. 3
[34] Benjamin Kurt Miller, Ricky TQ Chen, Anuroop Sriram, and Brandon M Wood. Flowmm: Generating materials with riemannian flow matching. In Forty-first International Conference on Machine Learning, 2024. 3
[35] Nayoung Kim, Seongsu Kim, Minsu Kim, Jinkyoo Park, and Sungsoo Ahn. Mofflow: Flow matching for structure prediction of metal-organic frameworks. arXiv preprint 14I::2210.11270, 2024. 3,5
[36] Hannes Stark, Bowen Jing, Regina Barzilay, and Tommi Jaakkola. Harmonic self-conditioned flow matching for joint multi-ligand docking and binding site design. In Forty-first International Conference on Machine Learning, 2024.
[37] Erwin Schrödinger. An undulatory theory of the mechanies of atoms and molecules. Physical review, 28(6):1049, 1926. 3, 22
[38] C. C. J. Roothaan. New developments in molecular orbital theory. Reviews of Modern Physics, 23(2):69 89, Apr 1951. doi: https://doi.org/10. http:///wooog//1111111eeeogggggp 103/revmodphys.23.69..
[39] G G Hall. The molecular orbital theory of chemical valency viii. method of calculating ionization potentials. Proceedings of the Royal Society of London, 205(1083):541-552 Mar 1951. doi: https://doi.or//// 10. 1090/..//............... 1951 9951990011....
12

---

[40] John C Slater. Atomic shielding constants. Physical review, 36(1):57, 1930. 3,25
[41] Peter MW Gill. Molecular integrals over gaussian basis functions. In Advances in quantum chemistry, volume 25, pages 141 141-205. Elsevier, 1994. 3, 25
[42] Milman, MH Lee, and MC Payne. Ground-state properties of cosi determined by total-energy pseudopotential method. Physical Review B, 49(23):16600,, 1994. 3
[43] Yvette Kosmann-Schwar/onn Schwarzbach et al. Groups and symmetries. Universitext, Springer, New York, 2010. 4. 23
[44] David J Griffiths and Darrell Schroeter. Introduction to quantum mechanics. Cambridge university press, 2018. 4, 25
[45] Jason Yim, Andrew Campbell, Emile Mathicu, Andrew YK Foong, Michael Gastegger, Jose Jimenez-Luna, Sarah Lewis, Victor Garcia Satorras. Bastiaan S Veeling, Frank Noe, al. Improved motif-scafion soaff-scafooldig with se (3) flow matching. Transactions on Machine Learning Research, 2024. 5
[46] Eugene Wigner. On the matrices which reduce the kronecker products of representations of groups. In The Collected Works of Eugene Paul Wigner: Part A: The Scientific Pupers, pages 608-654 Springer, 1993. 6, 33
[47] Nathaniel Thomas, Tess Smidt, Steven Kearnes, Lusann Yang. Li Li, Kai Kohlhoff, and Patrick Riley. Tensor field networks: Rotation-ane translation-equivaane equiivariint neural networks for 3d point clouds. arXiv preprint urXiv:/802.088) 19, 2018. 6, 24
[48] Yi-Lun Liao and Tess Smidt. Equiformer: Equivariant graph attention transformer for 3d atomistic graphs. arXiv preprint 2tti:::200.119990 2022. 6
[49] Stefan Chmiela, Alexandre Tkatchenko, Huziel E Sauceda, Igor Poltavsky, Kristof T Schütt, and Klaus-Robert Müller. Machine learning of accurate energy-consere molecular force fields. Science advances, 3(5):e6603015, 2017. 3,35
[50] Lars Ruddigkeit, Ruud Van Deursen, Lorenz C Blum, and Jean-Louis Reymond. Enumeration of 166 billion organic small molecules in the chemical universe database gdb-17. Journal of chemical information and modeling, 52(1 12224-282755 2012. 7,35
[51] Raghunathan Ramakrishnan, Mia Hartmann, Enrico Tapavicza, and O Anatole Von Lilienfeld. Electronic spectra from tddft and machine learning in chemical space. The Journal of chemical physics, 143(8), 2015. 7, 35
[52] John P Perdew, Kieron Burke, and Matthias Ernzrrhof. Generalized gradient approximation made simple. Physical review letters, 77(18):3865, 1996. 7, 35
[53] Florian Weigend and Reinhart Ahlrichs. Balanced basis sets of split valence, triple zeta valence and quadruple zeta valence quality for h to m.: Design and assessment of accuracy. Physical Chemistry Chemical Physics, 721833329733355 2005. 7,
[54] Philip J Stephens, Frank Devlin, Cary Chabalowski, and Michael J Frisch. Ab initio calculation of vibrational absorption and circular dichroism spectra using density functional force fields. The Journal of physical chemistry, 16645):11623-11122 1994. 7, 35
[55] Nikita Kornilov, Petr Mokrov, Alexander Gasnikov, and Aleksandr Korotin. Optimal flow matching: Learning straight trajectories in just one step. Advances in Neural Information Processing Systems,37, 37:104180-114204,, 2024. 10
[56] John F Cornwell. Group theory in physics: An introduction, volume Academic press, 1997.
23
[57] Eugene Wigner. Group theory: and its application the quantum mechanics of atomic spectra, volume 5. Elsevier, 2012. 24
[58] Wolfram Koch and Max C Holthausen. A chemist's guide to density functional theory. John Wiley & Sons, 2015. 25
13

---

[59] Zun Wang, Chang Liu, Nianlong Zou, He Zhang. Xinran Wei, Lin Huang, Lijun Wu, and Bin Shao. Infusing self-consistency into density functional theory hamiltonian prediction via deep equilibrium models. In The Thirty-eighth Annual Conference on Neural Information Processing
Systems, 2024. 33
[60] Zhuoran Qiao, Matthew Welborn, Animashree Anandkumar, Frederick Manby, and Thomas Miller. Orbnet: Deep learning for quantum chemistry using symmetry-adapted atomic-orbital features. The Journal of chemical physics, 153(12), 2020.
[61] Frank Neese, Frank Wennmohs, Ute Becker, and Christoph Riplinger. The orca quantum chemistry program package. The Journal of chemical physics, 152(22), 2020. 35
[62] Qiming Sun, Timothy C Berkelbach,, Nick S Blunt, George Booth, Sheng Guo, Zhendong Li, Junzi Liu, James D McClain, Elvira R Sayfutyarova, Sandeep Sharma, et al. Pysef: the python- based simulations of chemistry framework. Wiley Interdisciplinary Reviews: Computational Molecular Science, 8(1),e1340, 2018. 35,
[63] Mario Geiger and Tess Smidt. arXiv:2207.09453, 2022.
e3nn: Euclidean neural networks.
arXiv preprint
[64] Ask Hjorth Larsen, Jens Jørgen Mortensen, Jakob Blomqvist, Ivano E Castelli, Rune Christensen, Marcin Dulak, Jesper Friis, Michael N Groves, Bjørk Hammer, Cory Hargus, al. The atomic simulation environment-4 a python library for working with atoms. Journal of Physics: Condensed Matter, 29(27):273002, 2017. 37
14

---

NeurIPS Paper Checklist
1. Claims
Question:
Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
Justification: Discussions the methods and empirical results can be found in Sections and
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made in the paper and important assumptions and limitations. A No or
NA answer to this question will not be perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalize to other settings.
is fine to include aspirational goals as motivation as long as is clear that these goals not attained by the paper.
2. Limitations
Question: Does the paper discuss the limitations the work performed by the authors?
Answer: [Yes]
Justification: Discussion on the experiments results are found in the Section
Guidelines:
The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those are not discussed in the paper.
The authors are encouraged to create separate "Limitations" section their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect how these assumptions might be violated in practice and what the implications would be.
The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on few datasets or with few runs. In general, empirical results often depend on implicit assumptions, which should be articulated.
The authors should reflect on the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution is low images are taken in low lighting. Or speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon.
The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness.
While the authors might fear that complete honesty about limitations might used by reviewers as grounds for rejection, worse outcome might be that reviewers discover limitations that aren't 't acknowledged judgment and recognize that individual actions in favor of transparency play impor- tant role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.
the paper. The authors should use their best
3. Theory assumptions and proofs
Question: For each theoretical result, does the paper provide the full set of assumptions and complete (and correct) proof?
Answer: [Yes]
15

---

We provide the full set of assumptions and proof of the theretical statements in Appendix B.
Justification:
Guidelines:
The answer NA means that the paper does not include theoretical results.
All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated or referenced the statement of any theorems.
The proofs can either appear in the main paper or the supplemental material, but if they appear in the supplemental material, the authors are encouraged to provide short proof sketch provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix or supplemental material.
Theorems and Lemmas that the proof relies upon should be properly referenced.
4. Experimental result reproducibility
Does the paper fully disclose all the information needed to reproduce the main ex- perimental results of the paper to the extent that it affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided or not)?
Question:
Answer: [Yes]
Justification: We provide the main objective function of our work with the description the implementation and hyperparamters.
Guidelines:
The answer NA means that the paper does not include experiments.
If the paper includes experiments, a No answer to this question will not perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
If the contribution is dataset and/or model, the authors should describe the steps taken make their results reproducible or verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is novel architecture, describing the architecture fully might suffice, or if the contribution a specific model and empirical evaluation, it may be necessary either make it possible for others replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way accomplish this, but reproducibility can also be provided via detailed instructions for how to replicate the results, access to a hosted model (e.g., in the case of large language model), releasing of a model checkpoint, or other means that appropriate to the research performed.
While NeurIPS does not require releasing code, the conference does require all submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) the contribution is primarily a new algorithm, the paper should make clear how reproduce that algorithm. (b) If the contribution is primarily new model architecture, the paper should describe the architecture clearly and fully.
(c) the contribution is a new model (e.g., a large language model), then there should either be way to access this model for reproducing the results or way to reproduce the model (e.g., with an open-source dataset instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome deseribe the particular way they provide for reproducibility. In the case of closed-source models, it may be that access to the model is limited in some way (e.g., to registered users), but it should possible for other researchers to have some path to reproducing or verifying the results.
5. Open access to data and code
16

---

Does the paper provide open access to the data and code, with sufficient instruc- tions faithfully reproduce the main experimental results, as described in supplemental material?
Question:
Answer: [Yes]
Justification:
We used the publicly available data. Code will be on GitHub after publication.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurIPS code and data submission guidelines (https://nips... cc/
for more details.
While we encourage the release of code and data, we understand that this might not be possible, "No" is an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (c.g., for new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurIPS code and data submission guidelines (https: //niss. wwtpcccppppppeeeeeeeccecccccceee gguidds/Codeefee
ssisssiiiccc for more details.
The authors should provide instructions on data access and preparation, including how to access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. If only subset of experiments are reproducible, they should state which ones are omitted from the script and why.
At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable).
Providing as much information as possible in supplemental material (appended the paper) is recommended, but including URLs to data and code is permitted.
6. Experimental setting/details
Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary understand the results?
Question:
Answer: [Yes]
Justification: Details are in Appendix D.
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core the paper to a level of detail that is necessary appreciate the results and make sense of them.
The full details can be provided either with the code, appendix, as supplemental material.
7. Experiment statistical significance
Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Question:
Answer: [Yes]
Justification: Improvement percentages (> 50%) are discussed in Section 5. Also, provide the variance of our model prediciton which is related with the nature of the generative model.
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by error bars, confi- dence intervals, statistical significance tests, least for the experiments that support the main claims of the paper.
The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall with given experimental conditions).
17

---

The method for calculating the error bars should be explained (closed form formula, call a library function, bootstrap, etc.)
The assumptions made should be given (e.g., Normally distributed errors).
It should be clear whether the error bar is the standard deviation the standard error of the mean.
It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have 96% CI, if the hypothesis of Normality of errors is not verified.
For asymmetric distributions, the authors should be careful figures symmetric error bars that would yield results that are out of range (e.g. negative rates).
show in tables or
If error bars are reported in tables or plots, The authors should explain the text how they were calculated and reference the corresponding figures tables in the text.
8. Experiments compute resources
For each experiment, does the paper provide sufficient information on the puter resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?
Question:
Answer: [Yes]
Justification:
Discussions are in Appendix D.
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage.
The paper should provide the amount of compute required for each of the individual experimental runs as well estimate the total compute.
The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary failed experiments that didn't 't make it into the paper).
9. Code of ethics
Question: Does the research conducted in the paper conform, in every respect, with the NeurlPS Code of Ethics /tttps//neeuriie wwtpppublic/thicstiitiiiiie
Answer:
[Yes]
Justification:
The research conforms the NeurIPS Code of Ethics.
Guidelines:
The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.
If the authors answer No, they should explain the special cireumstances that require deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g., if there is special consid- cration due to laws or regulations in their jurisdiction).
10. Broader impacts
Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?
Answer: [Yes]
Justification: The work has positive impact on the drug discovery and material science, and this is discussed in Section and 2.
Guidelines:
The answer NA means that there is no societal impact of the work performed.
If the authors answer NA or No, they should explain why their work has no societal impact why the paper does not address societal impact.
18

---

Examples of negative societal impacts include potential malicious unintended uses (e.g., disinformation, generating fake profiles, surveillance), faireess considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.
The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster.
The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional unintentional) misuse of the technology.
If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how system learns from feedback over time, improving the efficiency and accessibility of ML).
11. Safeguards
Does the paper describe safeguards that have been put in place for responsible release of data models that have high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [NA]
Justification:
The work does present issues of high-risk misuse.
Guidelines:
The answer NA means that the paper poses no such risks.
Released models that have high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere usage guidelines restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make best faith effort.
12. Licenses for existing assets
Are the creators or original owners of assets (e.g.. code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?
Question:
Answer: [Yes]
Justification: All data sources and baseline models open-sourced. They have been properly credited and mentioned, as outlined in Section 4 and Appendix D.
Guidelines:
The answer NA means that the paper does use existing assets.
The authors should cite the original paper that produced the code package or dataset.
The authors should state which version of the asset is used and, if possible, include URL.
The name of the license (e.g., CC-BY 4.0) should be included for each asset.
For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should be provided.
19

---

If assets are released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, paperswi thcode. com/c datasttss has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
For existing datasets that are re-packaged, both the original license and the license of the derived asset (if has changed) should be provided.
If this information not available online, the authors are encouraged to reach out to the asset's s creators.
13. New assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [NA]
Justification:
No new assets.
Guidelines:
The answer NA means that the paper does not release new assets.
Researchers should communicate the details of the dataset/code/moode as part of their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent obtained from people whose asset used.
At submission time, remember to anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file.
14. Crowdsoureing and research with human subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well details about compensation (if any)?
Answer: [NA]
Justification:
No crowdsourcing experiments were conducted.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Including this information the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then as much detail as possible should be included in the main paper.
According to the NeurIPS Code of Ethics, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector.
15. Institutional review board (IRB) approvals or equivalent for research with human subjects
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or equivalent approval/review based the requirements of your country institution) were obtained?
Answer: [NA]
Justification: The paper does not involve crowdsoureing nor research with human subjects.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing research with human subjects.
Depending the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
20

---

We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors adhere the NeurIPS Code of Ethics and the guidelines for their institution.
For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review.
16. Declaration of LLM usage
Does the paper describe the usage of LLMs if it is an important, original, or non-standard component of the core methods in this research? Note that if the LLM is used only for writing, editing, formatting purposes and does not impact the core methodology. scientific rigorousness, or originality of the research, declaration is not required.
Question:
Answer: [NA]
Justification: Our methodology is not involved with the LLM.
Guidelines:
The answer NA means that the method development in this research does not involve LLMs as any important, original, or non-standard components.
Please refer to our LLM policy (https: //ttrr//s.u...pee / for what should or should be described.
21

---

A Additional Preliminary
A.1 Kohn-Sham density functional theory
Density functional theory (DFT) is cornerstone of quantum chemistry and materials science, offering a computationally tractable framework for approximating the electronic structure many- body systems. As directly solving the many-electron Schrödinger equation [37] is prohibitively expensive for systems with many particles, DFT reformulates the problem by expressing ground-stan -tate properties functional of the electron density p(r). Since p(r) depends only three spatial coordinates, this significantly reduces the computational cost regardless of the number of electrons.
In practical applications, DFT is most commonly implemented via the Kohn-Sham formulation [1, 12], which introduces a fictitious system of non-interacting electrons designed to reproduce true ground-state density of the real interacting system. The wavefunction of this non-interacting system, often referred to as the Kohn-Shan wavefunction V/ Ks. is defined as follows. Since electrons are fermions, V/ Ks must obey the Pauli exclusion principle. This requires the wavefunction to be antisymmetric with respect to the exchange of any two electrons, a property that is achieved by constructing V ks Slater determinant built from the N single-particle Kohn-Sham orbitals [1 I,M
4/1 (r1)
(12(
ty (r1)
t/ (12)
(12))
un (12)
(Ks( FN) = VN!
det det(ub(())) (Ty))
(15)
VN!
(TN)
(2(
(n((N)
Originally, the electron density is calculated from the wavefunction v::
p(r) N S.VV..VV...
Fw)]drg ... drv-
(16)
However, assuming orthonormal orbitals, the density expression simplifies to:
p(r) Σ (.(r)13
(17)
Instead of explicitly constructing the full KS wavefunction, which scales exponentially with system size due to the determinant, we can directly obtain the electron density from the single-particle orbitals dhg The Langle-particle orbitals (r)((r) and their corresponding energies {& are found by solving the Kohn-Shand (KS) equations:
1 v2 Vea(r) Vril)(() Viele((r) (b(r) (v)a()), H(pl)r(r)) [
(18)
Here, the entire term in brackets is the effective KS Hamiltonian Híp]. which functional of the electron density p. It consists of the Laplacian operator (V²) for the kinetic energy of the non- interacting electrons, the external potential (Ver) from the atomic nuclei, the Hartree potential (Vr) representing classical electrostatic repulsion, and the crucial exchange-correlation potential (V,)), which encapsulates all complex many-body quantum effects.
The ground-state electron density p(r) is then constructed from the N occupied KS orbitals obtained from solving these equations:
N
p(r) = Σ (T)
(19)
It is important to note that the KS orbital dyy is not, by itself, physical wavefunction. It is mathematical construct introduced solely to yield the same ground-state electron density as the real, interacting system. Additionally, the Hamiltonian H(p) depends on the density p. and p itself is constructed from the orbitals the that solve the KS equations, therefore forming self-consistent problem.
22

---

A.2 Group theory
We briefly introduce key concepts from group theory, providing the necessary background under- stand the motivation behind equivariant designs in Hamiltonian prediction.
Groups. A group G' is non-empty set equipped with binary operation o GG G. satisfying the following axioms:
Associativity: For all a, b, c e G, (a 6) oc a (bo
Identity: There exists identity element G such that for all G. G.eoar aac Inverse: For each E G, there exists an inverse a E G such that a a a oa c.
Closure: For all a, G. the result a is also in G.
We typically simplify notation by omitting o, writing ab instead of b. Groups can be finite or infinite, discrete or continuous, compact or non-compact. Group theory [43] provides a rigorous mathematical foundation for describing symmetries in physical systems [56], an essential perspective for designing physically informed machine learning models.
Group actions. To make group theory applicable to structures in mathematics and physics, we define how groups act on sets, A group G is said to act on a set X' there exists a function xx Gx + X such that for all g, G and x: E X, the following conditions hold: (1) x 2; where is the identity element of G. and (2) (gh) xx' x = g 'x (h 'x I). These conditions ensure that the group operation is compatible with its action to geometric transformations, such as translations, rotations, or reflections, on points or functions defined over Euclidean space.
the set. In practice, such group actions often correspond
When working with structured data such as molecular geometries or fields, is often desirable that the operations f X + y is said to be equivariant with respect to group actions xx on the input space X and -y on the output space ) if it satisfies
Equivariance.
apply preserve the symmetries of the underlying space. A function
flg x = y f(x), Ag G, E.
(20)
This means that applying a group transformation before or after the function yields consistent results. In physical systems, such as those governed by rotational symmetry (e.g. SO(3) in molecules), equivariance ensures that rotated inputs yield correspondingly rotated outputs. This property is
critical in designing models that generalize well across symmetrically equivalent configurations.
Group representations. Group actions on vector spaces are formalized through the notion of representations. A representation of group G on a vector space V is homomorphism p GL.(V,, where GL(VV denotes the group invertible linear transformations on V. This mapping associates each group element with linear operator that describes how the group acts the space. The vector space V is referred to as the representation space, and its dimension defines the representation's 's size.
real spaces = Rd, these representations are often expressed as invertible matrices D(g) E Rd For instance, the group SO(3) of 3D rotations is represented by orthogonal matrices E SO(3) and its 3D rotation matrix R E p3x3 acting on vectors x e IR via matrix multiplication D(g)x Rx.
In
Two representations D(g) and (g) are said to equivalent if they are related by a change of basis:
D' D(g)
(21)
for some invertible matrix Q Such equivalence indicates that the two representations describe the same group action under different coordinate systems.
A.3 Irreducible representations
Irreducible representations (Irreps). A representation p: G GL(V) of group G a vector space V is called irreducible if it does not contain any nontrivial invariant subspace. More formally,
'The subscript of group action operator for indicating the input and output space is often omitted when clear from context, iee. 9 x g h.
23

---

subspace W _ V is said to be invariant under the representation p if for all group elements g E G and all vectors w W, the action satisfies:
p(g)w E W.
(22)
Then, the representation p is irreducible the only invariant subspaces admits are the trivial subspace {0} and the entire representation space V.
Irreducible representations play a fundamental role because any finite dintessnmnnii unitary representa- tion of a compact group can be uniquely decomposed into direct sum of irreducible representations, up to isomorphism [57]. This decomposition foundational in analyzing the symmetry properties and constructing equivariant neural network architectures.
For the rotation group SO(3), the irreps completely characterized by non-negative integers l No- Each irreducible representation D SO(3) ((2++1)x(2(+)) has dimension 2/++1 and is explicitly realized by the action on spherical harmonics Y,, (0, Φ). For a function f : S - can be expanded in the basis of spherical harmonics
f(0,0)
Σ Σ
fim (0,
(23)
1=0 m=-t
where fim C is the complex coefficients of the spherical harmonics. The action of rotation R E SO(3) on ff is defined as:
(R ()0 = S(R (0, =
D(C mame (R. Sem T (0, Ø),
(24)
1=0 m=
where (R) is the Wigner D-matrix corresponding the irreducible representation labeled by C(24+1),
acting linearly on the vector of expansion coefficients f( (f-6f-411...... Explicitly, under rotation:
(()) R ( (R)
(25)
Thus, irreps serve the fundamental building blocks in constructing functions and features that transform equivariantly under the action of SO(3).
Cledsch-Gordan tensor product. Given irreps D(A), D(4) of SO(3), their tensor product D(4) & ((fa) acts the tensor product space and is generally reducible. It decomposes into direct sum of irreps:
fit
D(a)) EcG ((ba)
D
(26)
t= 11
E2l
where the decomposition is governed by MobschhGGrrdan (CG) coefficients. CG coefficients enable a basis transformation from the product space to irreps, playing crucial roles in quantum angular momentum coupling and equivariant neural networks.
Explicitly, for vectors u (f) E ((2h+1) and (f2) E C() their CG tensor product produces an irreducible feature (f) E ((2+1);
w(()
Σ
clem)) u (
(21)
(27)
(Ga,)))
This construction ensures that (() transforms correctly under R E SO(3):
(R)(()))(())) Sea (D((*) (R)v((s)) = D( (R)
(28)
Thus, CG tensor products ensure equivariant transformations, fundamental to maintaining consistency and symmetry in equivariant neural architectures.
Tensor field networks (TFN). The TFN [47] is widely used architecture that achieves SE(3)- equivariance by explicitly encoding features as irreps of SO(3) and using spherical harmonics to process directional information. Each feature vector node is decomposed into spherical tensor
components (() C(2++)), where & denotes the order of the irreducible representation. The key
24

---

component of the TFN layer is equivariant convolution and self interaction. For the equivariant convolution, the filter function is generated by the spherical harmonic functions Y (Ft)) are applied to the direction between atoms i and i.e, Paj = (r: r/)//|r; rj |l and are modulated by learnable radial function R([ra,]])) Il), implemented as an MI.P to produce the filter:
Fllner) (ray) (( (I|r,; lI) (f1)).
(29)
The equivariant convolution aggregates information from neighboring nodes set / using the tensor product between the filters and neighboring features:
Vilea) = =E
P (elinet) ( (mi) @ v!(w)) (fau)
(30)
JEN
where four E {lem - .,...... (enn + Ey)}. cquivariantly under rotation.
This convolution ensures that the result transforms
The self-interaction step applies a linear transformation over the channel dimension to each irreps independently:
(o(0) ECITI
Σ
Weee
V(() ie m
(31)
preserving equivariance due to the linearity and diagonal structure in the irrep indices. By construc- tion, all operations in TENs are equivariant to SE(3), allowing them to model tensorial and directional quantities while respecting physical symmetries. This makes them especially well-suited for applica- tions in molecular modeling, quantum chemistry, and materials science, where equivariance is crucial for generalization and physical accuracy.
A.4 Symmetry and equivariance of RH-DFT
Spherical harmonics and atomic orbital basis. To numerically solve the electronic wavefunctions are typically expanded in basis of atomic orbitals. These orbitals are constructed using spherical harmonics Y (0, σ)), which form complete orthonormal basis on the sphere S and transform under the irreducible representations ( of the rotation group SO(3). An atomic orbital basis function takes the form [58]:
equations,
pntm (r) R, (II=I)) (0,
(32)
where (ne([1 ([[r]]) is radial function (e.g., Gaussian [41] Slater-type [40, 15]), and Y captures the angular component. Here, n is the principal quantum number, l is the orbital angular momentum quantum number, and E {-l,.... (} is the magnetic quantum number [44].
These basis functions are central to Kohn-Sham DFT calculations, as molecular orbitals are expressed as linear combinations of atomic orbitals:
(r)
Σ
Cos Pnem (r-rrr))
(33)
oo (nfm)
where C E BBO is the orbital coefficient matrix, B is the number of basis functions, and O is the number of occupied orbitals.
Block Hamiltonian matrix. The Hamiltonian matrix in RH-DFT can be decomposed according to the angular momentum quantum numbers grouping associated with the atomic orbitals. Specifically, let (e,e)) ((2(-1)x(26' +1) denote the submatrix representing the coupling between two angular momentum numbers and f. Under rotation R E SO(3), each such sub-matrix H ((,e") transforms equivariantly according to the SO(3) irreps as follows:
H [D(R)C] ( (R) (H (())(()(()))
(34)
where atomic orbitals are grouped by their and e", D(R) is block diagonal Wigner-D matrix, and De and DC are Wigner D-matrices for and f', respectively. These submatrices H(e,") called orbital-wise Hamiltonian blocks, are the minimal building blocks of the full Hamiltonian that respect sO(3)-equivariance and have (2f + 1) (26 + elements which corresponds to the combination of the magnetic quantum number m and m'.
25

---

Full Hamiltonian matrix. In practice, basis functions are grouped not only by their angular momentum, but also by the atom they are centered on. For atom i, we denote its angular momentum support as I.j = (fi, (2.,2...... where each Eg corresponds to its atomic shell (e.g., 8, p. d) represented in the basis set.
From this atomic partitioning, the Hamiltonian can be reorganized into block form where each submatrix H!!ee) corresponds to the interaction between angular momentum of atom and E L, of atom j. The full Hamiltonian H is thus a structured matrix composed of K K blocks, where K = E Lil is the total number of irreps (angular momentum groups) across all atoms, where |L; denotes the cardinality of the set.
To construct the full Hamiltonian matrix, we first form the atom-pair matrix Hij for every atom-pair (i, j). This atom-pair matrix gathers all interactions between the angular momentum orbital groupings on the two atoms:
(G.44
(.)
H
E E Li, l E Ln,
(35)
where H (())) i and El E Tu of atom j.
E (22++)) denotes the interaction between angular momentum & E L; of atom
Next, place these atom-pair blocks into global block matrix obtain the full RH Hamilto- nian matrix: The full Hamiltonian matrix H E Cuxn, where Ee + 1), is assembled by concatenating the atom-pair blocks:
H
H.
H
(28 + 1), EEL
Σ Σ
(36)
for every atom-pair (i,j) in the molecule M.
For better intuition, we illustrate schematic structure of the Hamiltonian matrix in Figure 4a. Here, blocks corresponding to s, p, and d orbitals are arranged according to their angular momentum, highlighting the blockwise symmetry pattern induced by the underlying spherical harmonics basis.
Rotational equivariance full Hamiltonian. The global action of a rotation R on the basis is encoded by the block-diagonal matrix follows Equation (3):
(R)
D(R) (R))(R)) B ... (R)
(
(37)
where each block (RR) corresponds the Wigner D-matrix for angular momentum l. The Hamiltonian matrix as a whole then transforms under rotation as:
H[C R H(R CJ R H[C D(R) H(C)
(38)
This global equivariance condition guarantees that the Hamiltonian transforms consistently with molecular rotations, preserving the physical structure of the system. Leveraging this symmetry is critical for constructing SE(3)-equivariant machine learning models, which not only respect physical laws but also improve model generalization and interpretability in quantum chemistry tasks. For better intuition, we illustrate schematic structure of the Wigner D-matrix in Figure 4b.
26

---

(8)
1s2=3s
2p
3p
3d
1s
2p
1s2s
(b)
1s2s3s
2p
3p
3d
1s
2p
1s2s
2p
2s
3s
38
2p
3g
3d
3d
2s
H
2p
2p
Figure 4: (a) Schematic illustration of the full Hamiltonian matrix H for water molecule (H₃O). Color intensity indicates the magnitude of matrix elements, with red representing larger values and blue representing smaller values. (b) Schematic illustration of the full Wigner D-matrix D corresponding H, where green denotes larger values and purple denotes smaller values. Gray solid and dashed lines separate molecular blocks and orbital blocks, respectively, corresponding to submatrices defined by atomic orbital pairs.
27

---

B Training and sampling algorithms, and fine-tuning objective
B.1 Training and sampling algorithms
We investigate two types of prior distributions for the initial Hamiltonian Ho: Gaussian orthogonal ensemble (GOE) and tensor expansion-based (TE) priors. For the GOE prior, we sample each matrix entry independently as Maj ~ ~N(00, o²), where we set o² = 1.0 for MD17 and x² = 0.1 for QH9.
For TE prior, we sample each irreducible component w/ by drawing its radial norm from LogNormal(1, o² 0.1) and applying uniform SO(3) rotation to construct the full equivari- ant basis, which is then mapped to the Hamiltonian space via tensor expansion.
These choices provide flexibility in modeling different types of initial distributions, depending on the target dataset and symmetry constraints. Formal definitions and proofs of equivariance for both priors are provided in Appendix C.I and Appendix C.4.
The training and sampling procedures based the above rotationally invariant priors summarized in Algorithm and Algorithm 2, respectively.
Algorithm QHFlow training procedure
Require: Dataset of molecular configurations {Mi}, target Hamiltonians HI,M, overlap matrix S.M. model fo
1: for each training step do
2: Sample minibatch B {M
3: for each M, in B do
4: Sample initial Hamiltonian Ho Po and time ~ (((0,))
5: Compute interpolated Hamiltonian: HAM = - t)Ho (HMM
6: Predict (())) (H,,,M S.M,M
7: Compute flow matching loss: L.
[ 1,M HIMM
8: end for
9: Compute C = / Ei=1 Li, update 0 + 0 VoL
10: end for
Algorithm 2 QHFlow sampling procedure
Require: Molecular configuration M, initial Hamiltonian Ho ~ Po, model fe. time discretization (to. (1o.t,...... (x) where to and tx
1: Initialize H H,
2: for
K do 3: Predict target Hamiltonian: H T.M + fo fo S,M,tt)
(9)
4: Compute conditional vector field: ut + (H H )/(1 tx)
5: Update Hamiltonian:
Heg +- H + (le41 - (x) vt,.0
6: end for
7: Output H, -H,x
These algorithms ensure that training aligns the predicted vector field with the true conditional flow, while inference produces final Hamiltonian via integration through the learned ODE trajectory. This design ensures SE(3)--quivariance and physical consistency throughout the training and prediction processes.
B.2 Fine-tuning objective
After pre-training with the flow matching objective, we further fine-tune QHFlow to enhance its spectral accuracy by introducing energy alignment objectives. To implement the energy alignment objective, adopt the weighted alignment loss (WALoss) from WANet [25].
28

---

Property of the RH-equation. The relationship between the Hamiltonian matrix H, the coefficient matrix C, and the orbital energy matrix governed by the Roothaan Hall equation:
HC SCc,
(39)
where S is overlap matrix. Using orthonormality condition CC = I, this relation simplifies to:
C' HC
(40)
This identity forms the foundation for designing spectral alignment objectives.
Weighted alignment loss (WALoss). WALoss encourages the predicted Hamiltonian H to match the spectral structure of the converged SCF Hamiltonian H Let C* and E* denote the eigenvectors and eigenvalues of H", respectively. We define WALoss as:
diag "HO*)"
LWALoss
(41)
where diag(-) extracts the diagonal elements and || |/2 denotes the L.2 norm.
To emphasize physically important states, we place larger weights eigenvalues up to the Lowest unoccupied molecular orbital (LUMO) i.e., the k + lowest-energy orbitals. By projecting H onto the fixed cigenbasis C WALoss promotes alignment with the ground-trutt spectrum. However, it is important to note that (C*) THC does not exactly diagonalize H unless H = H* Nevertheless, minimizing this discrepancy helps guide H toward the correct cigenspectrum, improving orbital energy prediction and SCF behavior.
Overall fine-tuning objective. Our final finetuning objective combines the flow matching loss with WALoss:
Lital LCFM + AwA LWA,
(42)
where Awa is the weighting coefficients controlling the contributions of the energy alignment terms, and we chose AWA 2.0 for the finetuning.
29

---

C Proofs and additional details
C.1 Proof of -euueequivariince of GOE
Let M E R'x* be a symmetric matrix drawn from the Gaussian orthogonal ensemble (GOE), characterized by:
M[M = 0,
Var(M,) = =0; with My = Mai. and independent entries for i S j. (43)
Let D D() denote orthogonal matrix representing Wigner D D-transformation corresponding
to rotation R SO(3), satisfying orthogonality i.e. DD = I. We define the rotated matrix as:
M' DMD
with Mi = E
DuDDD
(44)
We aim to show that the distribution of M' is the same as that of M, i.e., that GOE invariant under conjugation by orthogonal transformations:
(M) (DMD
(45)
To establish this, it suffices to verify that the first and second moments of the entries in M' match those of GOE. Mean. Since [[[Ma] for all k 1, we have:
B[M
E = 0.
(46)
R
Covariance. We now compute the covariance of the entries Mab and M
Σ DaeeaaM
E[M Madl =E
Σ
Dos Dog
(47)
i,j
Ril
Σ
Dby Dek Dai B[Mi Mal].
(48)
AJ.kd
Using the independence and zero mean of GOE entries, we know:
Mpl o 2.
(49)
Thus,
E[M Mc = o² ab
Σ
Das Dbj Des Dej
(50)
0²
Dai
Σ D Dos
(DD Jac (DD )od
Dacd
(51)
This matches the covariance structure of the original GOE matrix. Since both M and M' are jointly Gaussian with zero mean and identical covariances, conclude:
DMD
M, i.e., GOE is invariant under orthogonal conjugation.
(52)
This proves that the GOE prior is iov3)-invariant, even for the high-rank Wigner D-matrices & 2 2.
30

---

C.2
Proof of the SEuia-equi of tensor expansion operation
Since the w() is the irrep vector, the m-th entry w transforms under rotation R E SO(3) as:
(f m
R R
w(1)
=
Σ Y
D mmn (R)
((())
(53)
where (R) (2x is the Wigner D-matrix.
By substituting above Equation (53) in Equation (9), the entry rotated tensor expansion is
(()(t)
(f1.62)
R
2 R w())
=
(h.m) ((amm)
Σ Do mm (R)
(54)
(m)))ag)
and exchange the order of summation:
(Sw(!)
(61,)))
R Σ
Σ
(l.m) (m.mm) (222
D mmm (R)
wmm
(55)
(m).mz)
Now crucial identity from representation theory:
Theorem 2. The Clebsch-Gorlar coefficients provide basis change between the product represen- tation ((G)) (R) & (((2)()) and the irreducible representation ()0) (R):
Σ
D mam! (R) = Σ
D(li) ming
(R) (R) mamb
(56)
((,,m))())
mi ml
We substitute the identity, then we have:
(!)
(E1.2) (mimmm))
((G)) (R) mimi
(R)
Σ
(emmm
w
(57)
mams
mi.mn
m
and the last sum implies the
(f,mn (ilm w (Sw(()) (22,f2)
(58)
(m/).))))
Therefore,
(w)())
(file2) (m)))
R
(())
D
(62)
(R) (3)) (6)
(59)
mimi
mams
(mf.m)))
in the matrix form.
(w)())
R
(R) (Sw!2) (
(60)
C.3 Proof of the invariance of the tensor expansion based (TE) prior
To construct an SO(3)-invariant distribution over matrices via tensor expansion, we begin with a rotationally invariant distribution over irreducible features (f) E ((++1), We define this distribution by factorizing it into two independent components:
A magnitude
||( drawn from radial distribution, e.g., ~ Normal(1, o²).
A direction sampled uniformly the sphere 8?, induced by a Haar- rotation
R. SO(3), applied to canonical unit-norm feature vector w (6) such that w(!) TD(0 (R) (()
31

---

This construction ensures that w(0) is distributed isotropically, i.., p(w()) =p(D(( (R)w(() for any R SO(3). Now, from the tensor expansion defined as:
(f (6)w(!) (41,62)
(61)
we know from the previous subsection that this construction transforms under SO(3) as:
(f.) = D(6)) (R) (6,,2) ((ea)
(62)
To prove invariance, consider the pushforward distribution:
(())) p(w = D( (R) VR E SO(3).
(63)
Because (a) is deterministie, equivariant function of (((), and ((w(9) is rotation-invariant the
induced distribution p(H (())) must also be invariant under conjugation by ( (R) and ( (R).
Thus:
((H( p(D( (R) H (f.62) D ( (R)
(64)
This proves that the prior distribution constructed via tensor expansion from rotationally invariant (C)) yields an SO(3) -Quui-equivarian (conjugation-invariin) distribution over matrix-valued outputs.
C.4 Proof of the invariance of the multiple tensor expansion prior
We now generalize the single-irrep vector expansion to construct a full Hamiltonian matrix E R2xx composed of multiple irrep vectors components. Let the atomic orbital basis consist of a set of irrep vectors indexed by subspace of dimension (2l; + 1).
{e1,
{B}, where each angular momentum & corresponds to
For each irrep vector rank 4, we define random irrep vector w(4) E C(24++1, sampled indepen- dently from a rotationally invariant distribution:
(6)
((E)) wR.) (0:)
r; Log,orraa(1,o"), R, Uniform(SO3))
(65)
We then define a factorized prior:
B
p(w (1))
w((l)) = p.(wwl.)).
(66)
with each (:(w(f being SO(3))invarian)
Using the tensor expansion, each pair (Ls, Ej) generates a block of the Hamiltonian matrix via:
(4.6) &e.,e wilk), for some la fle, f,,,..,, eg}.
(67)
The full Hamiltonian then assembled as:
H H iJ
(68)
From the single-component proof, we know that each block transforms as:
H ()))) R D(t.) (R)H((H) D(G) (())(R)
(69)
Therefore, under global rotation R SO(3), the full matrix H transforms
H R, D(R) HDR))
(70)
where D(R) the block block-diagonal rotation operator acting on the full basis.
Since each (6)) is sampled independently from an SO(3)-invariant distribution, the joint distribution of irreducible components p(H) remains invariant under this global conjugation:
p H) = p (D(R) HD(() YR E SO(3).
(71)
This completes the proof that the full Hamiltonian matrix constructed via multiple tensor expansions, using independently sampled spherical features from invariant priors, defines an -O(3)-invariant distribution over matrices.
32

---

D Implementation details
D.1 Model implementation
We build our model upon the official QHNet [23] and DEQHNet [59] codebases, which implement -qquvaqqivaarian GNN for Hamiltonian matrix prediction. QHNet is selected as our backbone for its balance architectural simplicity and strong predictive performance. Although more recent models such WANet [25] and SPHNet [26] have demonstrated improvements, their codebases not publicly available and thus are not considered here. Notably, our method is model-agnostic and could be paired with more expressive backbones if desired.
We extend QHNet by incorporating additional physically meaningful inputs: the current Hamiltonian matrix H,, the overlap matrix S, the molecular configuration M, and the time l. Unlike the original QHNet, which processes only M, our model is designed to handle H,, and as well, following ceemetty-rreservig strategies inspired by DEQHNet and OrbNet-Equi [60].
Project Block. To map the equivariant matrix features and S into the SO(3) irrep vectors. we used diagonal reduction in OrbNet [60] based on the Wigner-Eckarl theorem [46]. For atom i, we first extract the atom-wise matrices, H; == H (see Equation (35)). Each element of H; carries two orbital index, (f.m) and (E', m²). We project these matrices onto rank-&, irrep vectors
(P.)
Ea)
h1 (forfo))
] E (((66+1) by
n (o))
Σ (E,6" (mLm)
Σ
H
(.mst m') ( mm m²) m.
(72)
where H. (C.mall m()
is an element of submatrix H. (e,e))
that corresponds to the m and m' index, and (..(e..m.)) which plays the role of Clebsch-Gordann
(e,e")
set of atom-wise projection weights (.e.m,at' m²)
projection weights. Because the map Equation (72) an SO(3) tensor contraction, the resulting feature vectors remain equivariant. The same procedure is applied to overlap matrix S, giving two parallel streams of SE(3)-equivarian atomic features that are subsequently fused in the network.
Construction of projection weights. To construct the coefficients atom-wise projection coefficients, compute three-center integrals involving orbital basis function o (f,m) and auxiliary spherical-tyy type basis functions क (l.m))
Qf (Ronst' .m)) (fo.
- (1/ (o/
(aamm) dr.
(73)
In practice, the angular parts of these integrals are related with the Clebsch-Gordan coefficients due to their relation to spherical harmonics:
(l m²) mltommm
/
2e (
(+) (f)) dr.
(74)
gg
D.2 Training objective
Our training objective adopts residual learning strategy, following prior works such as WANet and SHNet. Rather than directly predicting the converged Hamiltonian, the model learns to approximate the residual between the initial and converged Hamiltonians using time-dependent vector field within the conditional flow matching framework. We define the residual target as done in prior works:
HIM := HM - (
(75)
(0) (0) where M is the initial guess of Hamiltonian matrix2, and H is the SCF solution. The final prediction of Hamiltonian obtained by summing the predicted residual and the initial Hamiltonian guess.
The conditional flow matching loss defined as:
L - (.H.
Ur (Hr,)) ue (HrM Hi,M)
(76)
H (0) 1.M HI,
(77)
We use the minao initialization from the PySCF package, and this guess does not contains SCF steps.
33

---

where the second line follows from the analytical form of the conditional velocity field ug. showing that the model is penalized based on the squared residual error at each time step. This formulation encourages smooth convergence from the initial guess to the target solution. To simplify the objective, we omit the 1/(1 1)2 time-scaling for our implementation.
34

---

E Experimental study settings
E.1 Dataset preparation
To demonstrate the effectiivnness of now-matching-bee bassed training, we conduct experiments on two molecular datasets: MD17 and QH9. The MD17 represents a relatively simple task compared to the QH9, focusing solely on small systems and their conformational space. The PubChemQH9 is not considered since their dataset and codebase are not publicly released.
MD17. The MD17 [49, 19] dataset consists of quantum chemical simulations for four small or- ganic molccules: water (H2O), ethanol OH), malondialdehyde (CH,(CHO)p), and uracil (C²,,N,02). It provides comprehensive set of molecular properties, including geometries, total energies, forces, Kohn-Sham Hamiltonian matrices, and overlap matrices. All reference com- putations were implemented via the ORCA electronic structure package [61] using the PBE ex- change-correlation functional [52, 53] and the def2-SVP Gaussian-type orbital (GTO) basis set. We follow the standard data split protocol used in prior work [19, 20, 23] to divide each molecule's con- formational data into training, validation, and test sets. The detailed dataset statistics are summarized in Table and MOs in the table imply molecular orbitals (i.e.. s.p.d.f)
Table 5: The statistics of MD17 dataset [19].
Dataset
# of structures
Train
Val
Test
f of atoms
# of orbitals
ff of occupied MOs
Water
4,900
500
500
3,900
24
5
Ethanol
30,000
25,000
500
4,500
72
10
Mulondialdehyde
26.978
25.000
500
1,478
90
19
Uracil
30,000
25,000
500
4,500
12
132
26
QH9. The QH9 dataset [24] is large-scale quantum chemistry benchmark designed support the training and evaluation of machine learning models for Hamiltonian matrix prediction across diverse chemical structures. It is based the QM9 [50, 51] molecular dataset and includes 130.831 Hamiltonian matrices from stable molecular geometries, as well as 2,698 molecular dynamics trajectories. The dataset covers small organic molecules composed of up nine heavy atoms (C, N. O. and F). All Hamiltonians were computed using the PySCF [62] quantum chemistry package with the B3LYP [54] exchange-correlatio functional and the def2-SVP Gaussian-type orbital (GTO) basis set. We provide the detailed dataset statistics in Table 6.
QH9 is organized into two main subsets, @H9-stable and 1199-dymamic33000 and provides four standard evaluation splits: stable-id, stable-ood, dynamic-300k-gen and dynamic-300k-mon The id split randomly partitions the QH9-stable subset into training, validation, and test sets, while the ood split is constructed based on molecular size. Specifically, it groups molecules with 3-20 atoms in the training set, molecules with 21-22 atoms in the validation set, and molecules with 23-29 atoms in the test set. This ood setup allows for a rigorous evaluation of out-of-distrriutiin (OOD) generalization in terms of molecular complexity.
The geo and mol splits are based on molecular dynamics (MD) trajectories, where each molecule is associated with 100 geometry snapshots. In the geo split, geometries randomly assigned within each molecule: 80 for training, 10 for validation, and 10 for testing. Thus, while all molecules are present across the splits, the specific conformations are disjoint, enabling assessment of geometric generalization. In contrast, the mol split divides the 2,698 molecules themselves into disjoint sets with an 80/10/10 vrannvaliiation//estt ratio. All 100 geometries of a molecule are assigned to the same subset. This split presents a more challenging task than the geo split, as it requires generalization to entirely unseen molecular identities rather than just new conformations.
E.2 Evaluation Metrics
To comprehensively evaluate the performance of Hamiltonian prediction models, we adopt several metrics that measure accuracy, physical fidelity, and downstream impact on quantum chemical properties. Below are detailed descriptions of each metric used in our experiments.
35

---

Table 6: The statistics of QH9 dataset [24].
Dataset
# of structures
# of Molecules
Train
Val
Test
Stable-id
130,831
130,831
104,664
13,083
13,084
Stable-ood
130,831
130,831
104,001
17,495
9,335
Dynamic-300k-2
299,800 8000
2,998
239,840
29,980
29,980
Dynamic-300k
299,800
2,998
239,800
29,900
30,100
Hamiltonian MAE. This metric measures the mean absolute error (MAE) between the predicted Hamiltonian matrix H and the ground-truth matrix H
MAE(H)
Σ JJ=1
Hij
H
(78)
T2
where n is the length of the Hamiltonian matrix (H E R_xx), and || 113 is the Frobenius norm. This metric directly reflects the quality of the predicted Hamiltonian in element-wise terms.
Occupied orbital energy MAE (Coce). This metric measures the MAE between the predicted and true orbital energies for only the occupied orbitals. Let é and et be predicted and the ground-trr truth generalized eigenvalues of the Hamiltonian, respectively, and let Ioce be the set of indices of the occupied orbitals:
MAE((uuc))
Σ TE
12- 1
(79)
card
This is particularly important for capturing the energy spectrum relevant to ground-state electronic properties. We identify the occupied orbitals by selecting the [N/2] lowest eigenvalues, where N is the number of electrons in the system.
Orbital coefficient similarity score (Sa). To compare the predicted and ground-truth molecular orbital coefficient matrices and C we use a cosine similarity- based score averaged over columns:
C.C.C (o(c,c)-((c.c) -CCICC
(80)
where C; is the i-th column vector of the coefficient matrix and (.--) denotes the inner product. This metric evaluates how well the predicted orbitals align with the ground truth.
HOMO, LUMO, and energy gap MAE (CHOM, LUMO, CA). These metrics quantify the predictive accuracy of specific frontier orbital energies, which are critical for determining electronic properties such as reactivity and charge transport. The highest occupied molecular orbital (HOMO) and lowest unoccupied molecular orbital (LUMO) energies, along with their difference (the HOMO-LUMO gap), are particularly sensitive model generalization, especially in challenging settings like QH9, where predictions are evaluated on unseen molecular structures or conformations, These metrics are computed as follows:
<HOMO
(81)
4LUMO
fLuuO
<iuMo
(82)
EA (ÊLUMO êHоMо) (Eiumo {hoMo)]
(83)
Here, (HомO corresponds to the [N/2|-tt lowest eigenvalue, and CLUMO to the (1N/2) + 1)-th, where N is the number of electrons in the system.
SCF Iter Ratio. This metric measures the reduction the number SCF iterations required when using a predicted Hamiltonian as the initial guess, relative to standard reference initialization (e.g..
from a minimal basis or atomic superposition guess). defined as
SCF Iter Ratio Iter
(84)
Iteræ
where pedd is the number of SCF iterations with the predicted initialization, and Iter is the number iterations under the reference setup. In our setting, we used the reference value, which measured from the conventional DFT metric to compare the acceleration of our method.
36

---

SCF T Ratio. This represents the reduction in total wall-clock time spent during SCF convergence using the predicted Hamiltonian. It reflects actual runtime savings in quantum chemical simulations:
SCF Ratio peed-
(85)
TSCF-set
where Tece-peed and TsCr-set denote the SCF runtimes under predicted and reference initializations, respectively.
Inf T Ratio. This metric captures the inference overhead of the predictive model itself. It the fraction of time required to generate the predicted Hamiltonian relative to the baseline SCF time:
Inf T Ratio Tinterncce
(86)
TSCF-ref
where inference is the runtime of the Hamiltonian prediction model.
Total T Ratio. This quantifies the net cost of using the predictive model, combining both inference and SCF convergence times. It provides holistic view of overall computational efficiency:
Total T Ratio Tinforence + TsCP-pred
(87)
TSCF-eef
Note. All energy values are reported in units of Hartree (1E,, = 27. 27.211eV), and similarity scores are unitless with values in [0, 1]. where higher is better. While the SCF Iter Ratio deterministic metric reflecting algorithmic convergence behavior, the wall-clock timing ratios (SCF T, Inf T, and Total T) are subject to variability due to system-level factors. In our experiments, all timing evaluations were conducted on A A100-80GB GPUs. However, due to resource allocation via SLURN and shared CPU environments, precise time measurements were not always reproducible. As such, reported time ratios should be interpreted as indicative trends rather than absolute benchmarks.
E.3 Experimental setup
Environment. Experiments were conducted using single GPU per model. For the MD17 dataset, we used NVIDIA RTX 3090 and A5000 GPUs, while for the QH9 dataset, experiments were performed on NVIDIA 100 80GB GPUs. Our implementation is based on Py'Torch 2.1.2 and PyG 2.3.0, both compiled with CUDA 12.1. Detailed package versions [62-64] and environment specifications will be released upon publication for full reproducibility.
Training was performed on a SLURM--anagee cluster, where slight variability in runtime was observed due system-wide GPU and CPU resource contention. Random seeds were fixed where possible to ensure training stability; however, minor non-determinism remains due to the inherent variability of distributed computing environments.
Training QHFlow took approximately 2-4 days on MD17 and 7-10 days on QH9, which is slightly longer than QHNet due to the added complexity of incorporating the current Hamiltonian and orbital information as inputs. Finetuning on QH9 required an additional 3-4 days. While the total training time could be further reduced with optimized infrastructure, experiments SLURM-managed environment with periodic reinitialization every 72 hours. Importantly, QHFlow is model-agnostic framework; with more scalable architectures, training efficiency can be significantly improved.
conducted under
Hyperparameters. For fair comparison, we adopted the same hyperparameters used by the bascline model QHNet [23] whenever possible. Our QHFlow shares the majority of its architecture and training settings with QHNet to ensure that performance improvements are attributable to the proposed flow matching design rather than hyperparameter tuning. A summary of the key hyperparameters used across different datasets is provided in Table 7.
For evaluating DFT acceleration performance, we used the PySCF [62] framework with the B3LYP exchange-correlation functional. The numerical integration grid level was set to 3, and the maximum number SCF cycles was 50. All other parameters followed the default settings of the PySCF RKS implementation.
37

---

Table 7: Training hyperparameters of QHFlow used across datasets.
Hyperparameter
Description
QH9
MD17 (Water)
MD17 (Others)
Leamingg Rate
Initial learning rate
Se-4
5e-4
Minimum Learning Rute
Minimum learning rate
le-7
lc-9
Ic-9
Batch Size
Number of molecules per hatch
32
10
Scheduler
Leaming rute scheduler
Polynomial
Polynomial
Polynomial
LR Warmup Steps
Warmup steps linearly increase learning rate
1,000
1,000
1,000
Max Steps
Maximum number of training steps
260,000
200,000
200,000
Fine tenintunng
Initial learning rate fine tuning stage
Fine-tuning Minimum LR
Minimum learning rate of fine-tuning Maximum number of Finetuning steps Prior distribution for flow matching
lc-7
Finetuning Steps
60,000
Prior Distribution
GOE/
GOE Thue False
GOE
Using If, Block
Use time-dee dee-eependent Hamiltonian IT, input
True
Thue
Using Block
Use overlap matrix as input
True
False
Model Order
Maximum degree of spherical harmonies
Embedding Dimension
Node feature embedding dimension
128
128
128
Bottle Hidden Size
Hidden size bottlencck layer
32
Number of GNN Layers
Number of graph neural network layers
5
5
Max Radius
Cutoff radius for neighbor search
15
15
15
Sphere Chunnels
Number of chunnels spherical busis
Edge Channels
Number of channels for edge features
32
32
32
38

---

F Additional experimental results and limitations
F.1
Ablation study of the finetuning objective
We conduct additional experiments to evaluate various fine-tuning strategies to improve the Hamilto- nian prediction performance QHFlow after standard flow matching pretraining. In this section, we describe each strategy in detail and provide empirical observations.
For the WALoss implementation, we follow the coefficient settings proposed in WANet [25] (AwA 2.0). During fine-tuning, we train for an additional 60,000 steps starting from the pretrained model with TE prior, using an initial learning rate of 10 and a polynomial learning rate scheduler.
Full Flow Matching with WA Loss training from scratch (WA-Full), In this strategy. we train the model from scratch by jointly optimizing the flow matching objective and WALoss for 260,000 steps, matching the default training steps as done in WANet. However, this approach sigmificantly degraded the Hamiltonian error. We hypothesize that imposing strong spectral constraints too early disrupts the learning of the global Hamiltonian structure, leading to unstable convergence and poor generalization.
WA Loss finetuning (WA-FT), After completing standard flow matching pre-training, we fine-tune the model solely using WALoss, which encourages alignment between the predicted Hamiltonian's eigenstructure and the ground-truth SCF solutions. This approach improves orbital energy prediction and SCF convergence compared to pure flow matching training, confirming the value of adding spectral supervision during fine-tuning.
Summary. Table 8 summarizes the quantitative results of all fine-tuning strategies. Here, w/o-FT refers to QHFlow trained solely with the flow matching objective using the TE prior. Among the approaches, WA-FT achieves the best balance between Hamiltonian prediction accuracy and practical DFT acceleration, demonstrating the effectiveness of spectral alignment in the fine-tuming theettnnng stage.
F22 Full results of DFT acceleration via SCF initialization
For completeness, include in Figure 5 the figure with the inset from Figure 3. As shown in Figure 5, the inference time of QHFLOW is longer than that of QHNet because QHFLOW performs multi-step inference for ODE integration. However, this additional cost negligible compared the overall DFT computation pipeline, as reflected in the SCF T ratio.
F.3
Full results
of ablation study
For completeness, we provide the full ablation study results in Table 9, and Table 10, including additional metrics beyond those presented in the main text. While the main manuscript focused primarily on Hamiltonian prediction error (H), occupied orbital energy error (Enez), and electronic density accuracy (S.) for clarity, here we also report the LUMO (eLuMo)) HOMO (eHOMo), and HOMO-LUMO gap (EA) energy errors
Overall, the trends observed in the full table are consistent with those discussed in the main text. In particular, the improvements in Hamiltonian accuracy and density prediction are accompanied by consistent reductions in orbital energy errors (CLUMO, (HMMO) and gap errors (cA). further demonstrat- ing the broad effectiveness of our proposed methods. These results reinforce the conclusion that flow matching enhance both Hamiltonian structure prediction and downstream electronic properties.
Also, better visualize the role of the prior distributions, we provide interpolation trajectories from the GOE and TE prior toward target Hamiltonians. Figure 6 shows example interpolations for both settings, highlighting how the initial distribution affects the learned flow dynamics. As seen, the TE prior results in smoother and more structured trajectories compared to the GOE prior, which starts from unstructured noise.
F.4 Ablation studies on model design
In this section, we present additional ablation studies on the architectural design choices of QH- FL.w, as summarized in Table 11, to clarify their individual roles and contributions the overall performance.
39

---

53.0
22.4
31.4
0,5 8
0.7 87
SCF nor Ranio
SCFT Rab
Inf Roco
Tothl Rat
SCF fer Raco
SCF Rato
Inf Rato
Total Rato
dal-ddyyrnia (300k)gge))
Oiti-dynamic (100k-mall
35.7
14.3
di11 46.6
14.9
064 den
25.4
4.1 467
SCF Bar Babio
SCF Rabio
Rako
Tntal Ratio
SC1 kar fakeo
SCFT Rabo
In T Ralio
Total Rabio
OHNict"
Ours
Ours (WA-FT)
Figure 5: DFT acceleration performance on 300 samples from the QH9 dataset. All metrics are reported percentages relative conventional DFT (initialized with minao), which serves as the 100% baseline. The SCF Iter Ratio measures the ratio of SCF iterations required, while Inf Ratio, SCF T Ratio, and Total T Ratio measure time.
Regression. QHFLOW (Regression) is baseline that directly predicts the Hamiltonian correction matrix (AH H' H)) using the same objective described above, without incorporating any physical priors or structural constraints. This variant isolates the effect of the flow-matching mechanism itself and evaluates how well the model performs under purely regression-basse setting. Although this approach achieves improved Hamiltonian prediction compared QHNet, its overall performance remains inferior to the flow-matching version.
Overlap matrix S. We further examine the impact of incorporating the overlap matrix embedding. which encodes pairwise correlations among atomic orbitals. Comparing models with and without embedding shows that leveraging overlap information significantly enhances model stability and accuracy across both Hamiltonian prediction and SCF convergence.
F.5 Analysis of error distribution of the predictions
To better understand the properties of QHFLOw, visualise the entire error distribution for six key quantities Hamiltonian MAE (H), occupied-orbital MAE (Eoce), ÉHOMO, ELUMO EA. and Se on all four QH9 splits (id, ood, geo, mol). Figure shows log-scale violin plots for QHNet (grey), QHFLOW (orange), and its WA-finetuned variant (cyan). (only the similarity score is kept on linear axis.) Each violin displays the full empirical distribution; the thick dashed line marks the median, and the two thinner lines divide the inter-quartile ranges.
beyond single-number summaries and
Across every split and metric, the QHFLOW violins are both lower and narrower. Median Hamiltonian error is reduced by roughly one order of magnitude on id and geo, and by factor of 3-5 on the tougher ood and mol sets. In addition, the long error tails of QHNer* largely vanish, demonstrating that our flow objective mitigates the occasional catastrophic failure modes that plague point-wise regression models.
Figure 8 plots log H against log face for every test geometry. QHFLOW occupics the bottom-lefi corner of each panel, while QHNet* points fan out towards larger joint errors, especially in ood and mol, where chemical diversity is highest. These results show the QHFLOW superior performance on both the geometric and unseen molecules.
Finally, Figure compares Hamiltonian error versus atom count. QHFLOW remains flat up to the largest molecules in QH9 (29atoms). This suggests that modeling a distribution over Hamiltonians confers form size transferability.
40

---

Table 8: Ablation experimental results on QH9 dataset Finetuning, Results shown in bold denote the best result in each column, whereas those that are underlined indicate the second best.
Dataset
Type
Gooc [pball J
S 1%]
GLMMO (pLal J
GROMO IpEal
<a [xE
Q11--stable (id)
WA-Full
62.69
153,82
98.96
240.07
117.26
251.81
W/o-FT
22.95 wozs
119.67 0211
99.51 wacs
437.96 7482
79.48-cume
553.87 7456
WA-FT
23.85 wozs
101.92 ozn
99.56
187.48 5454
92.22 0234
206.15 aww
QH9-stable (ood)
WA-Full
57.11
110.34
98.12
416.08
113.55
387.83
WTo-FT
20.01 1oozs
84.54 0
wocs
321.20
130.74 d.MS
395.83
WA-FT
20.55
72.64
99.16
171.24
77.96_
179.57
QH9-dynamic
WA-Full
75.38
108.96
98.77
231.63
114.14
219.20
(300k-geo)
wo-FFT
25.94 n/0
103.11+
99.59 00001
425.18
175.18=0025
547,33
WA-FT
27.12 1o022
89.03 + 0.15
99.65 woos
136.63
02.17-0021
154.68 4.49
QH9-dynamic
WA-Full
103.47
692.17
97.02
1046.21
760.85
1408 36
w/o-FT
45.91 noo
om
98.65
1344.68 7.88
4
1605.03 2.706
(300k-mo1)
WA-FT
46.60a 0033
424.75 0.504
98.74 0000
912.10= 2.01
403.51
1047.88 2,043
Table Full experimental results on QH9 dataset along the initial distribution. Results shown in bold denote the best result in each column.
H
Focs 4 [10
f [10
CLUMO | [10 [10-
<HOMO & [10 [x.]
CA [10- [x,]
Dataset
Prior
(10-9E,)
QH9-stable
GOB
25.93 0001
154.65 1997
99.39 noon
638.03
220.49
764.4627M8 7.15
(id)
22.95 0001
119.61 0211
99,51 0.001
437.96 7.82
179.48 = 0.00
553.87
QH9-stable (ood)
GOE
21.93 0001
87.32+
98.95
0.002
382.09 13:87
120.16 0.000
432.93+
TE
20.01 0001
84.54 0.0on
99.04 0.00
321.20
1.497
130.74 0.433
395.83
QH9-dynamic
GOE
29.39 2 0001
122. 14
99.49 0.001
618.75 3,099
215.41
756.45
(300k-geo)
TE
25.94 0.001
103.11 0003
99.59 0,001
425.18= 1.119
175.18 0255
547.33
QH9-dynamic
GOE
46.78 2 on01
419.68 0.001
98.65 n/e0n
1409.07=
478.12 12 0.14
1718.80
(300k-ma1)
TE
45.91 + 0.01
0.171
98.65 0.001
1344.68 2:08
479.71 0.150
1605.03 2.216
Table 10: Full predictive variance results QH9 dataset. WA-FT implies the finetune the model
with WA loss. Metrics for Ours and Ours (WA-FT) are meansistd over five random seeds. Results shown in bold denote the best result in each column, whereas those that are underlined indicate the second best.
Dataset
Model
H [1ff,]]
feee [p 4
S, [%]]
'LMO [e 4
HOOOO [pEiJ 4
<a (xPN
QHNet*
77.72
963.45
94.80
18257.34
1546.27
17822.62
QH9-stable (id)
WANet
80.00
833.62
96.86
SPHNet
45.48
334.28
97.75
Ours
22.95=0amm
119.67 ezu
+0
437.96+ vAs
179.48 aas
$53.87= CAs
Ours (WA-FT)
23.85
101.92
99,56:
187.48
92.22
206.15
QHNet"
69.69
884.97
93,01
25848.83
1045.99
25370.10
H99-tablle (ood)
SPHNet
43.33
186.40
98. 16
20.01=00001
84.54+
99
321,20
130.74
395,83-
Ours (WA-FT)
20.55 0.202
72.64
99.16: oxes
171.24 0273
77.96: 0205
179.57 9271
QHNet"
88.36
1170.50
93.65
23269.41
2040.06
22407.96
O99-dyaamic
WANet
74.74
416.57
99.68 68
SPHNet
$2.18
100.88
99.12
(300k-gao)
Ours
103.11
09,59 0001
425. 1.1
18
$47.33
Ours (WA-FT)
27.12
89.03
99.65
136.63 4061
84.17
154.68
QHNet
121.39
$554,36
86.02
53505.09
4352.76
50424.86
QH9-dynamic (300k-no1)
SPHNet
108.19
Ours
442,56
98 65
1344.68 2.se
479.71
1605,03- 2.09
Ours (WA-FT)
46.60=
424.75 0204
98.74 0o01
912.10 2041
403.51
1047.88 2003
Table 1I: Ablation study on the model design of QHFLOW. The models are trained on the QH9-Stable (id) split. Bold indicate the best.
Model
H [aEb]]
Guv J [eEil
S. [%)
GLUMO IpEid
снпмо [(EEs)]
ca 4 [|bEa]]
QHNet*
77.72
963.45
94.80
18257.34
1546.27
17822.62
QHFLOW (Regression)
35.85
213.45
98.89
1548 6448.64
330.77
1783.99
QHFLOW (Wo embed)
25.24
159.41
99.38
734.40
239.18
851.33
QHFLOW (W/ emhed)
22.95
119.67
99.51
437.96
179.48
553.87
41

---

(a)
t= 0
t
t 2/3
t=1
(b) GOE
H
H
t=0
t= 1/3
t= 2/3
t=1
Figure 6: Visualization of Hamiltonian interpolation trajectories from different prior distributions of H2O: (a) tensor expansion-based (TE) prior and (b) Gaussian orthogonal ensemble (GOE) prior. Color intensity indicates the magnitude of Hamiltonian elements across flow time
9+3
for (4)
Fammo Fil
Kunna FidS
FAI
WMT
NM-1
WM-FT
Gg
10T
Imco
Imcll
UTAMM 0.0
WM/7
AM77
110 0.5 1677
WKLT
WA-7
Figure Violin plots compare QHNet* (gray), QHFLoW(orange), and its weighted-alignmeet fine-tune (WA-FT, cyan) on six metrics: Hamiltonian MAE (H), occupied-orbital MAE (Eoce). Fhомо- eLUMO. (A), and coefficient similarity (S.). Results are shown for the four evaluation splits id, ood, geo, mol. All axes are logarithmic except S.. Thick dashed lines mark the median; thin dashed lines indicate the quartiles. Across every split the QHFLowviolins are both lower and narrower than those of QHNet", and WA-FT delivers an additional (though smaller) improvement on the energy- related metrics.
42

---

H vs. Eocc (id)
H vs. éocc (ood)
10a
10
QHNet*
Ex: [(E,,)]
Ours
Eog [2F
10
10
WA-FT
10a
10
10
10a
10
10
10a
10a
10
10a
10a
log H [[ES]]
log H [nE,]]
H (gea)
H v5. (mol)
10a
10a
log Eace [uEn]
10
log Fooe. [uEn]
10
10a
10a
10a
10a
10
140 10
10
10
10a
10
10
log H [WEa]
log H [uEn]
Figure 8: Scatter plots of log versus log foce for each QH9 split. Points for QHNer* (gray) extend to the upper right, revealing frequent large joint errors. QHFLOW (orange) and WA-FT (cyan) cluster tightly in the lower-left region, indicating that the flow-matching model achieves simultaneously low Hamiltonian and occupied -obbital errors, even on the challenging ood and mol splits.
Number of atoms vs. (id)
Number atoms vs. H (ood)
100
500
QHNet Ours WA-FT
log [nEn]
log (NE,)
100
40
20
20
10
10
3
5
15
20
25
30
23
24
25
26
27
28
29
Number of atoms
Number atoms
Number of atoms vs. (geo)
Number atoms vs. H (mol)
100
log [uEn]
log H (uEn)
100
40
40
20
10
10
25
Number of atoms
Number atoms
Figure 9: Mean log Hamiltonian error (markers) and one-standard-deviation bands (vertical bars) versus atom count for the four QH9 splits. QHNet* (gray) errors rise with molecular size in ood, whereas QHFLOW (orange) and WA-FT (cyan) remain nearly flat up to 29 atoms. This suggests that modelling the full Hamiltonian distribution with symmetry-aware priors confers robust transferability to larger, more complex molecules.
43


Reviews:
[{"review_id": "V4C37adaks", "reviewer": "Reviewer_YiLn", "summary": "DFT is an approximation to solve the Schroedingers Equation. It has received nobel prizes. DFT is computationally demanding and takes hours of compute for a single molecule. In the past, starting with Snyder et al PRL 2012 ML approaches have been used to estimate DFT densities. Other work followed. Also a direct prediction of the Hamiltonian has been proposed (Schütt et al 2019). In the meantime a number of approaches have been conceived. Also some that attempt to directly improve and accelerate the DFT loop. \n\nThe present ms proposes to use higher-order flow matching for DFT level Hamiltonian prediction in an equivariant fashion and presents impressive improvements over SOTA. Also SCF iterations in the DFT loop can be accelerated using the proposed approach.", "questions": "One area that could benefit from further clarification is the chemical meaning of the observed progress. While the results are promising, it would strengthen the paper to elaborate more explicitly on what the improvements represent in terms of underlying chemical principles or phenomena. Providing a clearer interpretation from a chemical standpoint would help bridge the gap between the computational findings and their practical or theoretical implications in chemistry.\n\nAdditionally, the paper would be improved by a more detailed discussion of the limitations of the proposed method. While the results are strong, every approach has its constraints, and a candid evaluation of where the method may fall short—whether in terms of scalability, generalizability, computational cost, or domain-specific applicability—would offer a more balanced perspective. Such a discussion would not only increase the transparency of the work but also guide future research directions.", "limitations": "yes", "rating": 5, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 4, "originality": 3, "strengths_and_weaknesses": "strength\n\nThe paper is clearly written, with well-structured sections and precise language that make it easy to follow the authors' arguments and methodology. The clarity of presentation contributes significantly to the overall readability and accessibility of the work.\n\nThe experimental section is impressive. The authors conduct a thorough and well-designed set of experiments that convincingly support their claims. The results are clearly presented, and the comparison with relevant baselines demonstrates the effectiveness of the proposed approach.\n\nThe paper is also well-situated within the existing body of literature. The authors provide a comprehensive and thoughtful review of related work, and they effectively position their contribution in relation to previous studies. This contextual grounding enhances the reader’s understanding of the paper’s novelty and relevance.\n\n\nweaknesses \n- it remains unclear to reviewer from the ms what the differences in the tables mean chemically, i.e. are the improvements insightful or chemically meaningless?\n(To give context to this question: many newer ML approaches in MLFFs provide improvements that are chemically meaningless and yield unstable MDs)\n- I haven't found limits discussed, please discuss them carefully", "was_revised": true, "final_justification": "I stay with my assessment in the light of the rebuttal and discussion. I would accept the ms.", "rebuttal": "Dear Reviewer YiLn,\n\nThank you for your thoughtful and encouraging review. We appreciate your recognition of the clarity of our presentation, the strength of our experiments, and the relevance of our contribution within the existing literature. Below, we address your comments in detail.\n\n---\n\n**W1, Q1: Are the improvements insightful or chemically meaningless?**\n\nOur improvements are chemically meaningful since they (1) reliably accelerate DFT calculations and (2) maintain chemical accuracy even for out-of-distributions. \n\nFirst, since QHFlow accelerates DFT (as shown in Fig. 3), it reduces the overall computational cost in cases where DFT was previously required. Even for out-of-distribution data, the additional SCF iterations ensure that the solver reliably achieves DFT-level accuracy. This is particularly impactful for chemical simulations, such as *ab initio* molecular dynamics (AIMD), which heavily rely on DFT calculations and often encounter out-of-distribution inputs.\n\nFurthermore, even without additional SCF iterations, QHFlow achieves the level of chemical accuracy (error < 1 kcal/mol ≈ 1,600 μEₕ [1,2]) for energy-related properties such as HOMO, LUMO, and the HOMO–LUMO gap, where baseline models fall short. Notably, as shown in Table 2, this level of accuracy is maintained even on unseen molecular configurations from molecular dynamics trajectories, specifically the QH9-dynamic (300k-mol) set.\n\n---\n\n**W2, Q2: Lack of discussion on limitations in terms of scalability, generalizability, computational cost, or domain-specific applicability.**\n\nThank you for the suggestion! We originally included the discussion in Appendix F.4, but now realize that its visibility is low despite its importance. To alleviate your concern, we will move the limitation to the main text and add more content to incorporate your comments. \n\nTo be specific, we will include the following text:\n> **(Generalizability)** We did not fully verify the generalizability of our approach, e.g., extension to new architectures (WANet, SPHNet) or datasets (PubChemQH), due to the lack of resources, public codes, and public datasets. **(Computational cost)** Our flow-based formulation introduces computational overhead for solving ODEs. In the future, one could consider removing the overhead with the one-step generative models [3]. **(Domain-specific applicability)** Our experiments are limited to gas-phase molecules with up to 29 atoms and two commonly used XC functionals (PBE and B3LYP). The algorithms are left to be verified for more general settings, such as periodic solids or higher-accuracy functionals (e.g., hybrid or double-hybrid). **(Scalability)** We also acknowledge that our datasets consist of relatively small molecules. One could further test the scalability of our approach once new datasets on larger systems are available. We also note that our work does not give a theoretical explanation of how Hamiltonian prediction error translates to SCF acceleration. Developing a theoretically grounded algorithm specifically for SCF acceleration would be an interesting future research.\n\nWe mark the keywords, e.g., generalizability and computational cost, as bold to highlight their relevance to your concerns. They will be removed from the final manuscript.\n\n---\n\nWe hope our responses have addressed your questions and concerns. Please let us know if there is anything we may have missed or if any points require further clarification.\n\n**Reference**\n\n[1] Pople, John A. \"Nobel lecture: Quantum chemical models.\" Reviews of Modern Physics 71.5 (1999): 1267.\n\n[2] Kuang, Yang, Yedan Shen, and Guanghui Hu. \"Towards chemical accuracy using a multi-mesh adaptive finite element method in all-electron density functional theory.\" Journal of Computational Physics 518 (2024): 113312.\n\n[3] Kornilov, Nikita, et al. \"Optimal flow matching: Learning straight trajectories in just one step.\" 2024 NIPS."}, {"review_id": "lotVuAvyi8", "reviewer": "Reviewer_NZgB", "summary": "The authors introduce QHFlow, a flow matching method for predicting converged Kohn-Sham density functional theory (KS-DFT) Hamiltonians. The authors introduce two novel SO(3)-invariant base densities in the same vector space of KS-DFT Hamiltonians. This enables the authors to train a flow matching model conditioned on the molecular structure mapping from one of the base densities to the Dirac distribution of the converged Hamiltonian. In their experiments, the authors find that their QHFlow reduces MAE errors significantly compared to direct Hamiltonian predictors.", "questions": "* Figure 3 confused me as the caption and the appendix read that the SCF T Ratio is the ratio of converging a DFT calculation starting from QHFlow vs minao. Though, then a factor of 46 would mean that the QHFlow initialization takes 46 times as long to converge as minao? However, the text reads it's 13-15% faster? An explanation on how the figure should be read would be greatly appreciated.\n* How sensitive is your approach to the choice of hyperparameters? Is there a good way to estimate these a priori?", "limitations": "No, the authors could discuss limitations such as scalability. It appears that the method does not scale to large systems, as the largest ones included here are from QM9.", "rating": 5, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 2, "originality": 2, "strengths_and_weaknesses": "Strengths:\n* The paper is well written, easy to follow, and provides the necessary background for understanding their work.\n* The authors propose an interesting application for flow matching-based regression.\n\nWeaknesses:\n* While the empirical performance motivates the work very well, the qualitative reasoning provided is weak. For instance, in line 39, it reads \"[...] approaches [...] often fall short of capturing the full distributional structure and uncertainty inherent in the Hamiltonian.\" With the exception of degenerate states, there's neither a distribution of Hamiltonian nor an uncertainty for a given structure. As you correctly note yourself in l.169, where you describe it as a Dirac distribution.\n* The technical novelty of this work is constrained to the two new base densities.\n* Hyperparameters vary quite a bit between datasets, and suggestions such as the WA Loss provide inconsistent results.", "was_revised": true, "final_justification": "I thank the reviewers for their thorough responses and increased my score accordingly.", "rebuttal": "Dear Reviewer NZgB,\n\nThank you for your thoughtful review. We appreciate your comments on the clarity of the writing with necessary background for understanding our work, and the proposed application of flow matching. Below, we address your points in detail.\n\n---\n\n**W1: While the empirical performance motivates the work very well, the qualitative reasoning provided is weak. There is neither a distribution of Hamiltonian nor an uncertainty for a given structure.**\n\nThank you for pointing this out. We agree that the Hamiltonian is deterministic for a given structure. Our goal is not to model the uncertainty in data, but to better capture the predictive uncertainty arising from model approximation. In what follows, we provide a similar argument to Han et al. [1].\n\nWe first clarify what we meant by distributional structure and uncertainty (in the prediction). Regression models typically minimize element-wise losses (e.g., MSE), implicitly modeling conditional independence between matrix elements: $p(\\mathbf{H}|\\mathcal{M}) = \\prod_{i,j} p(H_{ij}|\\mathcal{M})$, where $\\mathbf{H}$ is the Hamiltonian, $H_{ij}$ is its $(i,j)$-th element, and $\\mathcal{M}$ denotes the input molecule. In contrast, conditional generative models, such as diffusion or flow-based models, introduce a prior $z$ that globally conditions the output (i.e. $p(\\mathbf{H}|\\mathcal{M}) = \\int_z p(z)\\prod_{i,j} p(H_{ij}|\\mathcal{M},z)dz$). This allows the model to better capture structural correlations among Hamiltonian elements. Even in deterministic tasks, this structure-aware generation leads to more coherent predictions, as evidenced in recent advances like AlphaFold3 [2], DiffCSP [3], image segmentation [4], and an object detection model [5]. \n\nWe further demonstrate how better capturing the uncertainty leads to improved performance. When the model is uncertain between two plausible Hamiltonians $\\mathbf{H}_1$ and $\\mathbf{H}_2$, a regression model often predicts their average $(\\mathbf{H}_1 + \\mathbf{H}_2)/2$, which may not correspond to a physically meaningful solution and results in high uncertainty. A generative model, however, can represent both modes explicitly, yielding sharper, more reliable outputs.\n\nEmpirically, our flow-based model achieves consistently lower MAE compared to regression baselines as shown in Tables 1 and 2, supporting its advantage in modeling both the structural dependencies and predictive uncertainty of Hamiltonian prediction.\n\nNevertheless, we fully agree that our earlier phrasing was ambiguous. We will revise the manuscript to clarify this point:\n\n>  (lines 38-39) These approaches rely on pointwise regression, which may lead to high predictive uncertainty and often struggle to capture the structural correlations present in the Hamiltonian matrix.\n\n---\n\n**W2: The technical novelty of this work is constrained to the two new base densities.**\n\nWe resonate with this concern. Nevertheless, we would like to re-emphasize that our core contribution is not to methodological novelty. Instead, we propose to approach the Hamiltonian prediction problem with generative models to incorporate the underlying structure, which demonstrates strong empirical performance. We consider this idea to be novel, since the existing works mainly focused on developing new predictive architectures instead of generative models. \n\n---\n\n**W3-1, Q2: Hyperparameters vary quite a bit between datasets. How sensitive is your approach to the choice of hyperparameters? Is there a good way to estimate these a priori?**\n\nWe appreciate the reviewer’s concern and would like to clarify that our framework is generally robust to hyperparameter choices. Across seven out of eight datasets, we use nearly identical training hyperparameters, except for the minimum learning rate and the use of S-block embedding. The only exception was MD17–water, which contains only 500 training samples; For this case, we adopted a smaller learning rate to ensure stable optimization while keeping the other model hyperparameters fixed.\n\nAs a general guideline, we recommend selecting the largest batch size based on available GPU memory, since increasing batch size or training steps consistently leads to improved performance. For other hyperparameters, we generally recommend following our setting. However, hyperparameter tuning is inevitable for achieving the best performance in machine learning. One can start by tuning the learning rates to further improve performance. \n\n---\n\n**W3-2: WA Loss provides inconsistent results.**\n\nWhile we do not fully understand the reviewer’s concern, we believe it refers to the inconsistent improvements of QHFlow+WA-FT on Hamiltonian MAE in Table 2. For this case, we would like to clarify that this is, in fact, consistent with the motivation of WA-FT loss. \n\nWe first note that the WA-FT loss was mainly designed to enhance energy prediction performance, which is consistently improved across datasets. When minimizing the weighted combination of the original QHFlow loss (for Hamiltonian prediction) and the WA-FT loss (for energy prediction), there exists a natural trade-off between the metrics. Our experiments show that one can find a \"sweet spot\" that yields high accuracies for both metrics, but the trade-off is inevitable.\n\n---\n\n**Q1: Unclear explanation of Figure 3 (SCF accelerating performance)**\n\nThank you for pointing this out. We will revise the explanation of Figure 3 (both the caption and the main text) to improve clarity.\n\nSpecifically, we will update the caption as follows:\n> (Figure 3 caption) The SCF Iter Ratio reflects the relative improvement in the number of SCF iterations compared to conventional DFT. All time-related values (Inf T Ratio, SCF T Ratio, Total T Ratio) are expressed as percentages relative to the total runtime of conventional DFT initialized with minao, which is set to 100%. Lower values indicate faster convergence. For example, an “Total T Ratio” of 46% means that QHFlow converges in 46% of the time required by standard DFT, including the model inference time. \n\nIn the main text, we will revise Section 5.3 as follows:\n> (line 321-323) As shown in Figure 3, QHFlow significantly accelerates SCF convergence. For instance, on QH9-stable (id), QHFlow achieves convergence in just 46% of the runtime compared to conventional DFT, corresponding to a 54% relative speedup. Compared to prior ML-based initializations such as QHNet (53%), QHFlow yields an additional 7%p reduction. These results demonstrate the practical utility of QHFlow for accelerating real-world DFT simulations.\n\nWhile we originally reported relative ratios (i.e., 13\\% $\\approx 1 - \\frac{0.46}{0.53}$), we agree that we should present the improvement more clearly, such as expressing it as a percentage-point difference (e.g., a 7%p gain from 53% to 46%).\n\n---\n\n**Limitations: Lack of discussion on limitations in terms of scalability**\n\nThank you for pointing this out! We originally included the discussion in Appendix F.4, but now realize its low visibility despite the importance. To alleviate your concern, we will move the limitation to the main text and revise the contents. \n\nTo be specific, we will include the following text:\n> We did not fully verify the generalizability of our approach, e.g., extension to new architectures (WANet, SPHNet) or datasets (PubChemQH), due to the lack of resources, public codes, and public datasets. Our flow-based formulation introduces computational overhead for solving ODEs. In the future, one could consider removing the overhead with the one-step generative models [6]. Our experiments are limited to gas-phase molecules with up to 29 atoms and two commonly used XC functionals (PBE and B3LYP). The algorithms are left to be verified for more general settings, such as periodic solids or higher-accuracy functionals (e.g., hybrid or double-hybrid). **We also acknowledge that our datasets consist of relatively small molecules. One could further test the scalability of our approach once new datasets on larger systems are available.** We also note that our work does not give a theoretical explanation of how Hamiltonian prediction error translates to SCF acceleration. Developing a theoretically grounded algorithm specifically for SCF acceleration would be an interesting future research.\n\nWe mark the text in bold to highlight its relevance to your concerns. They will be removed from the final manuscript.\n\n---\n\nWe hope our responses have addressed your questions and concerns. Please let us know if there is anything we may have missed or if any points require further clarification.\n\n**References**\n\n[1] Han, Xizewen, Huangjie Zheng, and Mingyuan Zhou. \"Card: Classification and regression diffusion models.\" NIPS 2022\n\n[2] Abramson, J., Adler, J., Dunger, J. et al. Accurate structure prediction of biomolecular interactions with AlphaFold 3. Nature 2024\n\n[3] Jiao, Rui, et al. \"Crystal structure prediction by joint equivariant diffusion.\" NIPS 2023.\n\n[4] Chen, Shoufa, et al. \"Diffusiondet: Diffusion model for object detection.\" ICCV 2023.\n\n[5] Baranchuk, Dmitry, et al. \"Label-efficient semantic segmentation with diffusion models.\" ICLR 2022.\n\n[6] Kornilov, Nikita, et al. \"Optimal flow matching: Learning straight trajectories in just one step.\" 2024 NIPS."}, {"review_id": "APZfP3I0B9", "reviewer": "Reviewer_upG2", "summary": "The authors propose a flow-matching framework for generating Hamiltonians corresponding to ground-state solutions of DFT. Namely, the method builds on the QHNet architecture while adopting a generative approach in place of direct regression. To ensure SE(3) equivariance during flow matching, the method introduces two SE(3)-invariant priors: a Gaussian orthogonal ensemble and a tensor expansion-based prior. Additionally, the model extends their approach with a post-hoc fine-tuning step using a weighted alignment loss to improve the accuracy of orbital energy predictions. Experimental results on the MD17 and QH9 datasets demonstrate that the method outperforms existing regression-based baselines across multiple metrics.", "questions": "**Questions**\n1. How does the method perform (only inference) on even larger systems e.g. Alanine tripeptide (42 atoms) and DHA (56 atoms) from the MD22 dataset?\n\n2. Is there any intuition why the TE prior generally performs better than GOE prior?\n\n3. Line 38-39: “However, these approaches rely on pointwise regression and often fall short 39 of capturing the full distributional structure and uncertainty inherent in the Hamiltonian.” Could the authors clarify what is meant by \"distributional structure and uncertainty\" in this context? To my understanding, there is generally a one-to-one mapping between structure and the ground-state solution, and it is unclear what uncertainty is referring to here.\n\n**Suggestions**\n1. Perhaps I missed it, but it would be helpful to explicitly state that Table 3 reports results from QHFlow using the TE prior. I had to match the numbers between tables to infer this connection.\n\n**Typos**\n1. Table 1: Malonaldehyde -> Malondialdehyde", "limitations": "yes", "rating": 5, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 4, "clarity": 4, "significance": 4, "originality": 2, "strengths_and_weaknesses": "**Strengths**\n1. The paper is well written and clearly structured.\n\n2. The proposed approach provides SOTA results by leveraging flow matching with equivariance and tailored prior to predict Hamiltonian matrices for the ground-state solution of DFT. The incorporation of SE(3)-invariant priors (GOE and TE) are novel and well-motivated by the symmetry properties of Hamiltonians.\n\n3. The experiments are comprehensive and rigorous, comparing QHFlow with other leading regression-based approaches (e.g QHNet, WANet, SPHNet) across various datasets. The various experiments and metrics strengthen the validity of the findings.\n\n**Weaknesses**\n1. While QHFlow demonstrates strong performance on small to medium-sized systems (e.g., MD17, QH9), including out-of-distribution results within QH9, its scalability to larger and more complex molecules such as those in MD22 is untested. Evaluating larger molecules would help establish applicability to larger chemistry workloads and help strengthen the paper even further.\n\n2. QHFlow primarily extends QHNet with the flow matching and symmetry-aware priors. Although, these proposed components lead to significant empirical improvements, the novelty of the approach is limited as it largely builds on existing methodologies and does not significantly advance methodology.", "was_revised": true, "final_justification": "I will maintain my score, as I believe it appropriately reflects the paper’s strong empirical performance and their rigorous experimental design, as detailed in my review. Although the paper achieves strong results using generative modeling and symmetry-aware priors for Hamiltonian prediction, the work does not introduce a fundamentally novel methodological innovation or conceptual insight. For that reason, I do not see a basis for increasing the score.\n\nThat said, the authors have addressed my questions and most of my concerns during the rebuttal, and I appreciate the clarity they brought to several points.", "rebuttal": "Dear Reviewer upG2,\n\nThank you for your thoughtful and encouraging review. We sincerely appreciate your positive feedback on the clarity of the manuscript, the novelty and motivation of our SE(3)-invariant priors, and the rigor of our experimental evaluations. Below, we address your comments in further detail.\n\n---\n\n**W1, Q1: Lack of evaluation on larger molecules such as MD22.**\n\nWe fully resonate with the reviewer's concern. However, this is infeasible at the moment since there are no public datasets with large molecule sizes to evaluate our method. Note that the MD22 dataset is not a quantum Hamiltonian dataset. We also note the existence of PubChemQH, which consists of molecules with 40–100 atoms, but the dataset is not public.\n\nNevertheless, we strongly believe that our improvements will generalize to larger systems. The improvements may be larger since flow matching (and continuous-time generative models) has proven to be quite scalable for large systems. We note that our computational cost is similar to that of predictive models for each training step. For inference, the cost of multiple ODE integration steps is negligible compared to SCF iterations.\n\n---\n\n**W2: The novelty of the approach is limited as it largely builds on existing methodologies and does not significantly advance methodology.**\n\nWe resonate with this concern. Nevertheless, we would like to re-emphasize that our core contribution is not to methodological novelty. Instead, we propose to approach the Hamiltonian prediction problem with generative models to incorporate the underlying structure, which demonstrates strong empirical performance. We consider this idea to be novel, since the existing works mainly focused on developing new predictive architectures instead of generative models.\n\n---\n\n**Q2: Is there any intuition why the TE prior generally performs better than the GOE prior?**\n\nRegrettably, this is more of an empirical finding rather than an insightful design. We hypothesize that the observed performance difference arises from the structural alignment between the TE prior and the Hamiltonian matrices. Specifically, the TE prior utilizes tensor expansion, which mirrors the structure and operations used to construct the Hamiltonian matrix in our model. Rather, GOE is constructed by the element-wise Gaussian distribution, which lacks this structural inductive bias.\n\n---\n\n**Q3: There is generally a one-to-one mapping between structure and the ground-state Hamiltonian. Clarification on “distributional structure and uncertainty.**\n\nThank you for pointing this out. We agree that the Hamiltonian is deterministic for a given structure. Our goal is not to model the uncertainty in data, but to better capture the predictive uncertainty arising from model approximation. In what follows, we provide a similar argument to Han et al. [1].\n\nWe first clarify what we meant by distributional structure and uncertainty (in the prediction). Regression models typically minimize element-wise losses (e.g., MSE), implicitly modeling conditional independence between matrix elements: $p(\\mathbf{H}|\\mathcal{M}) = \\prod_{i,j} p(H_{ij}|\\mathcal{M})$, where $\\mathbf{H}$ is the Hamiltonian, $H_{ij}$ is its $(i,j)$-th element, and $\\mathcal{M}$ denotes the input molecule. In contrast, conditional generative models, such as diffusion or flow-based models, introduce a prior $z$ that globally conditions the output (i.e. $p(\\mathbf{H}|\\mathcal{M}) = \\int_z p(z)\\prod_{i,j} p(H_{ij}|\\mathcal{M},z)dz$). This allows the model to better capture structural correlations among Hamiltonian elements. Even in deterministic tasks, this structure-aware generation leads to more coherent predictions, as evidenced in recent advances like AlphaFold3 [2], DiffCSP [3], image segmentation [4], and an object detection model [5]. \n\nWe further demonstrate how better capturing the uncertainty leads to improved performance. When the model is uncertain between two plausible Hamiltonians $\\mathbf{H}_1$ and $\\mathbf{H}_2$, a regression model often predicts their average $(\\mathbf{H}_1 + \\mathbf{H}_2)/2$, which may not correspond to a physically meaningful solution and results in high uncertainty. A generative model, however, can represent both modes explicitly, yielding sharper, more reliable outputs.\n\nEmpirically, our flow-based model achieves consistently lower MAE compared to regression baselines as shown in Tables 1 and 2, supporting its advantage in modeling both the structural dependencies and predictive uncertainty of Hamiltonian prediction.\n\nNevertheless, we fully agree that our earlier phrasing was ambiguous. We will revise the manuscript to clarify this point:\n\n>  (lines 38-39) These approaches rely on pointwise regression, which may lead to high predictive uncertainty and often struggle to capture the structural correlations present in the Hamiltonian matrix.\n\n---\n\n**Suggestions and Typos**\n\nThank you for your helpful suggestions. We will revise Section 5.4 to clearly indicate that the results shown in Table 3 correspond to QHFlow with the TE prior, as reported in Table 2.\n\nSpecifically, we will update section 5.4 as follows:\n> (line 333-334) Table 3 compares the performance of GOE and TE priors, **where the TE prior results are identical to those of QHFlow reported in Table 2.** We found that the TE prior consistently yields lower errors than the GOE prior across all QH9 splits. This highlights the importance of designing task-aligned and symmetry-aware priors for accurate Hamiltonian prediction.\n\nWe mark the text in bold to highlight its relevance to your concerns. They will be removed from the final manuscript.\n\nWe also appreciate the correction regarding the molecular name in Table 1 and will correct it in the final version.\n\n---\n\nWe hope our responses have addressed your questions and concerns. Please let us know if there is anything we may have missed or if any points require further clarification.\n\n**References**\n\n[1] Han, Xizewen, Huangjie Zheng, and Mingyuan Zhou. \"Card: Classification and regression diffusion models.\" NIPS 2022\n\n[2] Abramson, J., Adler, J., Dunger, J. et al. Accurate structure prediction of biomolecular interactions with AlphaFold 3. Nature 2024\n\n[3] Jiao, Rui, et al. \"Crystal structure prediction by joint equivariant diffusion.\" NIPS 2023.\n\n[4] Chen, Shoufa, et al. \"Diffusiondet: Diffusion model for object detection.\" ICCV 2023.\n\n[5] Baranchuk, Dmitry, et al. \"Label-efficient semantic segmentation with diffusion models.\" ICLR 2022."}, {"review_id": "nruAPou7Wd", "reviewer": "Reviewer_UyjR", "summary": "Deep learning approaches for prediction of the Hamiltonian matrix offer a way to improve the efficiency of density functional theory (DFT) calculations by forgoing or reducing the number of self-consistent field (SCF) steps needed. The authors propose a SE(3)-equivariant flow matching model for the prediction of this Hamiltonian. This is achieved through the introduction of a symmetry-aware prior distribution and orbital energy finetuning procedure. On several datasets, this method demonstrates improved performance in Hamiltonian error metrics, and reduction in SCF steps compared to prior work that predicts the Hamiltonian through regression.", "questions": "* How much of the performance improvement over QHNet is due to architectural changes (see *Weaknesses*) vs. the introduction of flow matching? It would be helpful to expand the ablation studies in section 5.4 to include either a QHNet baseline with the architectural improvements and energy alignment fine-tuning, or a QHFlow implementation without the architectural improvements, so we can better understand the source of the improved performance.\n* Fig. 3: Given QHFlow is evaluated with \"3-step ODE-based sampling\" wouldn't we expect the inference time ratio for QHFLow to be roughly 3 times that of QHNet? While the bar colors of \"Inf T Ratio\" can be inferred by order, it would be helpful to modify this figure so that this is more obvious.\n* What datapoints are the \"71%\" M17 and \"53%\" QH9 Hamiltonian derived from? The authors should be more precise with these statements in both the abstract and results section. Measuring the mean % improvement over the most recent state-of-the-art method would make the most sense to me.", "limitations": "The authors do not include any limitations of their work. One notable limitation in this line of research is that Hamiltonian prediction error is not always well correlated with SCF step reduction. I.e. a model that has a very low $\\mathbf{H}$ MAE might not necessarily improve DFT runtimes, which could suggest a different task if SCF step reduction is the main goal. This does not seem to be the case in this work (at least in aggregate) but could be worth discussing (along with other things).", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 4, "significance": 3, "originality": 3, "strengths_and_weaknesses": "### Strengths\n\n* *Clarity*: the paper is generally well written and easy to understand.\n* *Scope*: the problem is well defined, and the work includes two distinct benchmark tasks.\n* *Originality*: while SE(3) equivariant flow matching models have been proposed in the past, this work is the first to apply them to the task of Hamiltonian prediction, which has generally been done through direct regression. This is non-trivial due to the higher-order nature of the Hamiltonian matrix.\n* *Significance*: the proposed method offers an improved accuracy of both Hamiltonian prediction and SCF step reduction compared to prior work.\n\n### Weaknesses\n * The authors make several modifications to the QHNet architecture in their method, including the additional input of the overlap matrix $\\mathbf{S}$, attention-based message-passing layers, and a channel mixing module. Furthermore, the introduced energy alignment fine-tuning could also be applied to the regression-based model. Given the focus of the work is on the proposed flow-matching method, it is unclear if the performance gains are due to the generative task or improvements in model architecture.", "was_revised": true, "final_justification": "The authors propose a SE(3)-equivariant flow matching model for the prediction of the molecular Hamiltonian from density functional theory (DFT). While there are existing methods for predicting this quantity through regression, this work is the first to demonstrate the use of a generative framework for the task. With this flow matching framework, they demonstrate an improvement over prior methods both in terms of prediction error, but also self-consistent field (SCF) step reduction in simulation on several benchmark tasks. They demonstrate through ablation that this is due to both architectural improvements and the proposed flow-matching method. \n\n* *Why I am not rating it higher:* While the work demonstrates clear improvements in the Hamiltonian prediction task, I would consider it more of an incremental improvement in modeling than a \"groundbreaking impact\". Furthermore the methods and analysis introduced likely have limited applicability in other areas of AI.\n* *Why I am not rating it lower:* The use of flow matching in the prediction of Hamiltonians is a novel and interesting contribution, and the authors provide adequate evaluation and analysis.", "rebuttal": "Dear Reviewer UyjR,\n\nThank you for your thoughtful and encouraging review. We greatly appreciate your recognition of the clarity of our writing, the originality of applying SE(3)-equivariant flow matching to Hamiltonian prediction, and the significance of our results in both Hamiltonian accuracy and SCF step reduction. Below, we address your comments in detail.\n\n---\n\n**W1, Q1: It is unclear if the performance gains are due to the generative task or improvements in model architecture.**\n\nWe appreciate the reviewer’s question regarding the source of performance improvement. To address this concern, we conducted an ablation study comparing our flow-based framework to a regression-based baseline under identical architectural conditions:\n\n| Model (QH9-id) | H mae $\\downarrow$ [$\\mu E_h$] | $\\epsilon_{occ}$ $\\downarrow$ [$\\mu E_h$] | $\\mathcal S_c$ $\\uparrow$ [%]    |   LUMO $\\downarrow$ [$\\mu E_h$] | HOMO $\\downarrow$ [$\\mu E_h$] | GAP $\\downarrow$ [$\\mu E_h$] |\n| -------- | --------: | --------: |--------: | --------: | --------: |--------: |\n| QHNet     |  77.72     | 963.45     | 94.80  |   18257.34 |  1546.27 |  17822.62 |\n| QHFlow (*Regression*)    | 35.85     | 213.45     | 98.89  |   1548.64 |  330.77 |  1738.99 |\n| QHFlow (GOE)    | 25.93     | 154.65    | 99.39 | 638.03  |  220.49  |       764.46 |\n| QHFlow (TE)   | **22.95**     | 119.67     | 99.51   |      437.96   | 179.48   | 553.87  |\n| QHFlow (TE, WA-FT)    | 23.85     | **101.92**     | **99.56**   |  **187.48**       |   **92.22**    |  **206.15**    |\n\nThese results demonstrate that both architectural modifications and the adoption of flow matching contribute meaningfully to our performance gains. Notably, switching from regression to flow-based training yields consistent improvements across all metrics, even under the same model architecture. Furthermore, the choice of prior distribution (TE vs. GOE) and the use of a fine-tuning strategy further highlight the trade-off between Hamiltonian accuracy and energy prediction.\n\nWe will include the complete ablation results in the final version and publicly release our codebase to support further research in this area.\n\n---\n**Q2-1: (Figure 3) Given QHFlow is evaluated with \"3-step ODE-based sampling,\" wouldn't we expect the inference time ratio for QHFlow to be roughly 3 times that of QHNet?**\n\nThis is due to additional latency that is constant across QHFlow and QHNet. The latency is caused by system-level overheads such as data movement and Python runtime costs. Once we make the (unrealistic) assumption that the latency is not part of the computational cost, QHFlow would cost roughly three times that of QHNet. \n\nMore importantly, we emphasize that the primary computational bottleneck in DFT workflows is the SCF iteration process. Since our method substantially reduces the number of SCF iterations, as shown in Figure 3, the added inference cost becomes negligible in the context of end-to-end simulation time.\n\n---\n\n**Q2-2: While the bar colors of \"Inf T Ratio\" can be inferred by order, it would be helpful to modify this figure so that this is more obvious.**\n\nThank you for the helpful suggestion. In the final version, we will revise the \"Inf T ratio\" part to be more visible and clear. Specifically, we will add a small inset above each “Inf T Ratio” bar group to zoom in on it and better highlight the differences.\n\nAs external link sharing is not allowed, we are unable to upload the revised figure. We will ensure to update it in the final manuscript.\n\n---\n\n**Q3: What datapoints are the \"71%\" M17 and \"53%\" QH9 Hamiltonian derived from?**\n \nThank you for the question. We report the mean relative improvement in Hamiltonian MAE compared to the previous state-of-the-art model, as reducing Hamiltonian error is one of the primary objectives. We measure the improvement by $1-\\frac{A}{B}$ where $A$ is the Hamiltonian MAE of our model and $B$ is that of the previous best model.\n\nFor example, we observed the following improvements on MD17:\n\n- Water: 58% = $1 - \\frac{4.93}{11.7}$\n- Ethanol: 73% = $1 - \\frac{5.33}{20.09}$\n- Malondialdehyde: 82% = $1 - \\frac{3.80}{20.67}$\n- Uracil: 80% = $1 - \\frac{3.68}{18.65}$\n- Average: 71%\n\nAnd we observed the following improvements on QH9 (without WA-FT):\n- Stable-id: 50% = $1 - \\frac{22.95}{45.48}$\n- Stable-ood: 54% = $1 - \\frac{20.01}{43.33}$\n- Dynamic-300k-geo: 50% = $1 - \\frac{25.94}{52.18}$\n- Dynamic-300k-mol: 58% = $1 - \\frac{45.91}{108.19}$\n- Average: 53%\n\nWe will revise the manuscript accordingly to reflect clarify how the improvements were calculated, as follows:\n> (line 15-16) QHFlow achieves state-of-the-art performance, reducing Hamiltonian error by 71% on MD17 and 53% on QH9 **compared to the previous best model**.\n\n> (line 289-290) With up to a 71% reduction in Hamiltonian prediction error **on average across the entire MD17 dataset** compared to the previous state-of-the-art model, QHFlow demonstrates strong predictive performance.\n\n> (line 297-298) As shown in Table 2, QHFLOW consistently outperforms all prior methods in all metrics up to a 53% reduction in Hamiltonian prediction error **on average over all QH9 splits.**\n\nWe mark the text in bold to highlight its relevance to your concerns. They will be removed from the final manuscript.\n\n---\n\n**Limitations: Lack of discussion on limitations; weak correlation between Hamiltonian prediction error and SCF acceleration.**\n\nThank you for pointing this out! We originally included the discussion in Appendix F.4, but now realize that its visibility is low despite its importance. To alleviate your concern, we will move the limitation to the main text and revise the contents. \n\nFurthermore, we agree with your points on the lack of correlation between Hamiltonian prediction error and SCF acceleration. We think this is very insightful and a very interesting venue for future research!\n\nIn the final manuscript, we will include the following text:\n> We did not fully verify the generalizability of our approach, e.g., extension to new architectures (WANet, SPHNet) or datasets (PubChemQH), due to the lack of resources, public codes, and public datasets. Our flow-based formulation introduces computational overhead for solving ODEs. In the future, one could consider removing the overhead with the one-step generative models [1]. Our experiments are limited to gas-phase molecules with up to 29 atoms and two commonly used XC functionals (PBE and B3LYP). The algorithms are left to be verified for more general settings, such as periodic solids or higher-accuracy functionals (e.g., hybrid or double-hybrid). We also acknowledge that our datasets consist of relatively small molecules. One could further test the scalability of our approach once new datasets on larger systems are available. **We also note that our work does not give a theoretical explanation of how Hamiltonian prediction error translates to SCF acceleration. Developing a theoretically grounded algorithm specifically for SCF acceleration would be an interesting future research.**\n\nWe mark the text in bold to highlight its relevance to your concerns. They will be removed from the final manuscript.\n\n---\n\nWe hope our responses have addressed your questions and concerns. Please let us know if there is anything we may have missed or if any points require further clarification.\n\n[1] Kornilov, Nikita, et al. \"Optimal flow matching: Learning straight trajectories in just one step.\" 2024 NIPS."}]

Author rebuttal:

Meta-review: ## **Summary of the Paper**
This paper proposes **QHFlow**, an $\mathrm{SE}(3)$-equivariant generative model for predicting Kohn–Sham Hamiltonians in a localized atomic orbital basis. Unlike prior deterministic regression approaches (e.g., QHNet), QHFlow frames prediction as a *flow matching* problem, mapping from a symmetry-aware prior (Gaussian Orthogonal Ensemble or Tensor Expansion) to the target Hamiltonian. The method also incorporates architectural changes such as attention-based message passing, channel mixing, and use of the overlap matrix $S$ as an additional input. The model is evaluated on MD17 and QH9 datasets, showing improved Hamiltonian MAE, derived property accuracy, and reduced SCF wall-time.

## **Strengths**
- Addresses a relevant and active problem in quantum chemistry with a novel generative framing.  
- Consistent, clear empirical improvements over strong baselines (QHNet, WANet, SPHNet) across multiple datasets.  
- Incorporates physically meaningful priors and evaluates on metrics directly relevant to downstream SCF convergence.  
- Well-written and accessible, with relevant background context.

## **Weaknesses**
- Main novelty (flow matching) is introduced alongside several other changes — notably the use of structured priors so the model is not starting from scratch, inclusion of $S$, and architectural updates— without sufficient ablations to attribute gains specifically to the generative formulation.  
- Scalability to larger systems not demonstrated (currently being addressed by authors with additional experiments on MD22).
- Minor presentation issues (initial ambiguity in “uncertainty” wording, figure clarity) were largely addressed in rebuttal.

## **Reviewer Consensus and Discussion**
Reviewers agree the performance gains are convincing and the generative approach is interesting, but several noted that attribution is unclear. In particular, the introduction of $S$ as an input and the use of structured priors are significant changes from baselines and may account for part of the improvement. While the rebuttal included some ablations, I recommend two **relatively light additional experiments** that could still be conducted during the camera-ready phase to clarify attribution and strengthen the contribution:  
1. QHFlow **without** $S$, to isolate the effect of the overlap matrix.  
2. A QHNet variant trained to predict $\Delta H$ from the same prior initialization (e.g., TE prior), to separate the effect of the structured prior from the flow mechanism.

**Note:** These ablations are my own suggestions as AC to address reviewer and community concerns—they were not explicitly requested by the reviewers. The NeurIPS 2025 camera-ready deadline is **October 23, 2025 AoE**, which may still allow time to run these tests before finalizing the manuscript.

## **Recommendation**
I recommend **Accept (Spotlight)**. This is the highest-scoring paper in my AC batch. The contribution is timely, broadly relevant, and demonstrates strong empirical gains over strong baselines, warranting wider visibility. While some attribution questions remain, they do not diminish the significance or potential impact of the work. Including the suggested ablations before the camera-ready deadline could make the case for the generative formulation even stronger.
Decision: Accept (spotlight)