Extended Lagrangian & Shadow QEq Methods

• Extended Lagrangian and Shadow QEq are innovative MD techniques that use auxiliary electronic variables to efficiently model charge transfer and polarization.
• They eliminate costly iterative SCF cycles by employing shadow potentials and harmonic restraints, ensuring strict energy conservation.
• These methods generalize to various flexible charge models, accurately capturing polarization effects in systems like water, alkane chains, and interfaces.

Extended Lagrangian and Shadow Charge Equilibration (QEq) models constitute a class of molecular dynamics (MD) methods designed to efficiently and accurately model electronic polarization and charge transfer phenomena within classical and quantum-based force fields. These approaches circumvent the computational bottlenecks inherent to iterative self-consistent field (SCF) optimization of charge variables at every MD step, while rigorously conserving energy and correctly sampling the associated shadow Born–Oppenheimer (BO) potential. The framework includes generalizations such as the shadow charge–potential equilibration (SChPEq) method based on the Atom-Condensed Kohn-Sham second-order (ACKS2) model, as well as shadow QEq methods for density functional–inspired variable-charge force fields.

1. Foundations and Motivation

Traditional flexible charge models, such as the charge equilibration (QEq) formalism, introduce atomic partial charges as dynamical variables, enabling on-the-fly computation of the BO surface. However, enforcing self-consistency for these electronic variables typically requires iterative solution of large linear systems at every MD step, leading to significant computational expense and, if convergence is incomplete, nonconservative forces and secular energy drift. The extended Lagrangian (XL) strategy resolves this by promoting (approximate) electronic or charge degrees of freedom (DOFs) to dynamical variables within an augmented Lagrangian, allowing for explicit reversible propagation and tightly controlled deviations from ground state (Niklasson et al., 2023, Goff et al., 2023, Stanton et al., 13 Feb 2025).

The concept of a shadow potential formalizes the idea that, for a given integration scheme and approximate electronic state, there exists a modified potential (the shadow BO surface) on which the molecular trajectory is exactly conservative up to a given truncation order in the time step.

2. Extended Lagrangian Formalism

The archetypal extended Lagrangian for shadow QEq or SChPEq methods consists of three essential terms:

• Nuclear kinetic energy: $$T_\text{nuc} = \frac{1}{2} \sum_i M_i |\dot{\mathbf{R}}_i|^2$$ for nuclei at positions $\mathbf{R}_i$.
• Auxiliary electronic/charge DOF kinetic energy: For atomic charges $n$, charge-potential pairs $(n,u)$, or charge-transfer variables $\xi_{ij}$, with fictitious masses $m$, e.g., $$T_\text{aux} = \frac{1}{2} m (\|\dot{n}\|^2 + \|\dot{u}\|^2)$$.
• Shadow potential with harmonic restraint: $U_\text{shadow}(R, n, u)$, plus $$U_\text{harm}(n, u) = \frac{1}{2} m \omega^2 \left[f(R, n, u)\right]^T T \left[f(R, n, u)\right]$$, where the residual $f$ measures deviation from the instantaneous ground-state solution (self-consistent charges and potentials).

The Lagrangian reads

$$L(R, \dot{R}, n, \dot{n}, u, \dot{u}) = T_\text{nuc} + T_\text{aux} - U_\text{shadow}(R, n, u) - U_\text{harm}(n, u)$$

or, for simpler QEq-like cases, omitting $u$, with obvious modifications (Stanton et al., 13 Feb 2025, Goff et al., 2023, Niklasson et al., 2023, Gergs et al., 2020).

The restraint matrix $T = K^T K$ is constructed from the Jacobian or approximated inverse Hessian of the electronic residual mapping, ensuring correct adiabatic separation and controlling dynamical fidelity.

3. Construction of the Shadow Potential

For generic QEq and ACKS2-type models, the shadow potential is derived as follows:

• ACKS2/SChPEq case: The starting energy function,

$$E_{\rm ACKS2}(q, v) = \chi^T q + (q - q^0)^T v + \frac{1}{2} q^T C q + \frac{1}{2} v^T X v$$

is Taylor-expanded with respect to off-diagonal (long-range) components to yield a hybrid shadow function $\mathcal{E}(q, v; n, u)$, then minimized exactly with global constraints $\sum_i q_i = Q_\text{tot}$, $\sum_i v_i = V_\text{tot}$ to provide $U_\text{shadow}(R, n, u)$ in closed form at each step (Stanton et al., 13 Feb 2025).

• QEq/XL-Shadow case: For QEq or second-order QEq–like charge models, the energy

$$E_\text{QEq}(R, q) = V_S(R) + \sum_i \chi_i q_i + \frac{1}{2} \sum_{i \neq j} \gamma_{ij} q_i q_j + \frac{1}{2} \sum_{i} U_i q_i^2$$

leads to a shadow energy $\mathcal{E}(R, \eta, n)$ with off-diagonal blocks linearized around $n$, and the shadow BO potential is $$U_\text{shadow}(R, n) = \min_\eta \{ \mathcal{E}(R, \eta, n)\ |\ \sum_i \eta_i = 0 \}$$, again solvable in closed form (Goff et al., 2023).

The minimization can be efficiently performed using block matrix inversion (e.g., Woodbury formula) with $O(N)$ complexity, as the only nontrivial couplings are of low rank and on-site blocks are strictly diagonal (Stanton et al., 13 Feb 2025). Iterative SCF cycles are thus entirely eliminated.

4. Equations of Motion and Integration Schemes

Applying Euler–Lagrange equations in the zero mass, high frequency ($m \to 0$, $\omega \to \infty$, with $m\omega^2=$~const) adiabatic limit gives:

• Nuclear motion:

$$M_i \ddot{R}_i = -\partial_{R_i} U_\text{shadow}(R, n, u)$$

or its QEq/XL analog, with the auxiliary variables $(n, u)$ treated as constant arguments when computing nuclear gradients.

• Auxiliary charge/potential DOF motion:

$$\begin{bmatrix} \ddot{n} \ \ddot{u} \end{bmatrix} = -\omega^2 K ([q(R; n, u) - n; v(R; n, u) - u])$$

or, for simpler charge-only models, $\ddot{n} = -\omega^2 K (q_{\min}(n) - n)$.

Integration is performed using a combination of velocity-Verlet or leapfrog schemes for nuclei and a modified Verlet/dissipative integrator for electronic variables, optionally accelerated by low-rank Krylov subspace updates for the preconditioned kernel $K$ (Stanton et al., 13 Feb 2025, Goff et al., 2023). The method propagates all variables in a symplectic, time-reversible manner.

For the QTE (charge transfer equilibration) model, the auxiliary variables represent pairwise charge transfer ($\xi_{ij}$), and the nuclei and $\xi_{ij}$ are integrated on different time scales using multiple time step techniques to ensure stability while maximizing computational efficiency (Gergs et al., 2020).

5. Avoidance of Iterative SCF and Energy Conservation

By design, the shadow-BO and extended-Lagrangian framework entirely removes the need for iterative self-consistency at each MD step. The auxiliary DOFs, initialized from previous step solutions, are harmonically restrained to remain close to the ground-state solution, ensuring that residuals are $O(\Delta t^2)$. As a result, numerical energy drift is reduced to below $10^{-6}$ eV/atom-ps, and the trajectory remains strictly microcanonical for the corresponding shadow Hamiltonian (Niklasson et al., 2023, Stanton et al., 13 Feb 2025).

Benchmark calculations on representative systems (e.g., water, alkane chains, UO$_2$) show total energy conservation over $10^5$ steps, with charge-dependent shadow energies matching direct self-consistent solutions to within $10^{-9}$ Ha and physical observables (IR spectra, polarizability scaling) reproduced to high fidelity (Stanton et al., 13 Feb 2025, Goff et al., 2023).

6. Model Generalizations and Applications

The XL-shadow paradigm is not specific to QEq, but encompasses a range of flexible charge models including:

• ACKS2/SChPEq: Incorporates both atomic charges and potential fluctuation (softness) terms, enabling physically correct charge fragmentation and improved polarizability scaling ($\alpha \sim N^{1.76}$ for alkane chains, compared to $N^{2.45}$ for QEq) (Stanton et al., 13 Feb 2025).
• QTE (pairwise transfer): Models both metallic and insulating behavior within a single framework and supports correct charge localization and finite-range transfer, critical for surfaces and interfaces. Applicability extends to reactive MD, surface phenomena (e.g., with mirror boundary conditions), and dissociation/recombination (Gergs et al., 2020).
• Flexible coarse-grained multipoles: Generalizations include higher-order charge fluctuations, machine-learned short-range descriptors (ACE, etc.), and fully nonlocal charge interaction kernels (Goff et al., 2023).

Key distinctions among related models are summarized below:

Model	DOF type	Charge transfer regime
QEq	atomic $q_i$	metallic limit, global SCF
QTPIE	atomic $q_i$	insulating limit via $\chi^\text{eff}$
QTE	pairwise $\xi_{ij}$	interpolates, finite-range
SChPEq	charge $q$, pot. $v$	correct fragmentation, softness scaling

7. Accuracy, Stability, and Practical Considerations

The error in forces and energies under XL-shadow integration scales as higher powers of the MD time step ($\Delta t$), specifically: $$\Delta E \lesssim C_E (\Delta t)^4,\quad \Delta F \lesssim C_F (\Delta t)^2,\quad \Delta q, \Delta v \lesssim C_q (\Delta t)^2$$ as demonstrated in explicit benchmarks for molecular clusters and condensed-phase systems (Stanton et al., 13 Feb 2025, Goff et al., 2023). The methodologies support practical time steps ($\Delta t =$ 0.5–2.0~fs) and operate at $O(N)$ per-step computational complexity for the electronic part, making them as efficient as or faster than conventional force-field MD with full polarization.

The use of extended Lagrangian propagation notably suppresses spurious resonances or artificial peaks in dynamical spectra, even with the presence of fast auxiliary DOF oscillations. Weak dissipative terms ensure numerical stability of the auxiliary sector without breaking microcanonicity for the physical subsystem.

In surface and interface simulations, mirror boundary conditions can be used to prevent nonphysical charge transfer into unconstrained or frozen degrees of freedom (Gergs et al., 2020).

---

References:

• "Shadow Molecular Dynamics for a Charge-Potential Equilibration Model" (Stanton et al., 13 Feb 2025)
• "Shadow molecular dynamics and atomic cluster expansions for flexible charge models" (Goff et al., 2023)
• "Shadow Energy Functionals and Potentials in Born-Oppenheimer Molecular Dynamics" (Niklasson et al., 2023)
• "Generalized Method for Charge Transfer Equilibration in Reactive Molecular Dynamics" (Gergs et al., 2020)