Centroid Molecular Dynamics (CMD) Overview

• Centroid Molecular Dynamics (CMD) is a path-integral method that maps quantum particles onto classical ring polymers to simulate quantum dynamical properties.
• It calculates the centroid potential of mean force by integrating out non-centroid modes, enabling the evaluation of real-time correlation functions and vibrational spectra.
• Modern adaptations like f-CMD, f-QCMD, and h-CMD improve computational efficiency and correct curvature artifacts, extending CMD’s applicability to complex condensed-phase systems.

Centroid Molecular Dynamics (CMD) is a path-integral-based methodology for simulating quantum dynamical properties of condensed-phase and molecular systems, especially those manifesting significant nuclear quantum effects (NQEs) such as zero-point energy and tunneling. By mapping quantum nuclei onto classical ring polymers and projecting dynamics onto the centroid of the polymer, CMD yields a classical-like evolution that approximately preserves quantum Boltzmann statistics and enables the calculation of real-time correlation functions and vibrational spectra. CMD serves as the foundational limit for a hierarchy of centroid-based and path-integral quantum dynamics methods.

1. Theoretical Foundations and Formal Structure

CMD arises from the path-integral representation of a quantum partition function. Each quantum particle is mapped to an $N$-bead ring polymer, with beads coupled harmonically to simulate quantum delocalization. The centroid (mean position) of the beads for each atom serves as the classical coordinate of interest: $$\mathbf{Q}_i = \frac{1}{N} \sum_{\alpha=1}^N \mathbf{q}_{i,\alpha}$$ where $\mathbf{q}_{i,\alpha}$ is the position of the $\alpha$-th bead of the $i$-th atom.

The central object is the centroid potential of mean force (PMF) $W(\{\mathbf{Q}_i\})$, defined by integrating out all non-centroid internal (normal) modes in the ring polymer: $Z = \int d\mathbf{Q}\ e^{-\beta W(\mathbf{Q})}$

$$W(\mathbf{Q}) = -\frac{1}{\beta} \ln \int \prod_{\alpha=1}^N d\mathbf{q}_\alpha\ e^{-\beta H_N(\{\mathbf{q}_\alpha\})} \delta\left(\mathbf{Q} - \mathbf{Q}(\{\mathbf{q}_\alpha\})\right)$$

where $H_N$ is the ring-polymer Hamiltonian. The effective CMD Hamiltonian then reads: $$H_{\rm CMD} = \sum_i \frac{\mathbf{P}_i^2}{2m_i} + W(\{\mathbf{Q}_i\})$$ CMD evolves only the centroid variables under classical Hamiltonian equations of motion generated by $W$, with all quantum effects formally encoded in the shape of the centroid PMF. The bulk of computation is in sampling the PMF by path-integral molecular dynamics (PIMD) and extracting its mean force.

Correlation functions of observables linear in position/momentum can then be computed as classical time-correlation functions sampled from the centroid ensemble, providing Kubo-transformed quantum time correlations (Jang, 2013).

2. Algorithmic Variants and Acceleration Strategies

Original CMD (“on-the-fly” or adiabatic) computes the centroid mean force dynamically by tightly thermostatting the high-frequency internal modes and evolving only the centroid in microcanonical (or weakly thermostatted) fashion, enforcing adiabatic separation. Several modern approaches accelerate this workflow or circumvent its cost.

Fast Centroid Molecular Dynamics (f-CMD):

• Reference PIMD simulations are used to tabulate centroid-averaged forces.
• A force-matching procedure fits a classical centroid PMF, $V_{\rm cm}(\mathbf{r})$, ensuring

$$-\nabla V_{\rm cm}(\mathbf{r}) \approx \langle -\nabla V(\mathbf{q}) \rangle_{\rm RP | \mathbf{r}_{\rm cent}}$$

• Production dynamics are then run as classical MD with the fitted $V_{\rm cm}$, eliminating on-the-fly bead propagation (London et al., 2024).

Fast Quasi-Centroid Molecular Dynamics (f-QCMD):

• For molecules with significant curvilinear (e.g., angular) fluctuations, f-QCMD defines quasi-centroids by averaging bond lengths and angles over beads.
• The quasi-centroid PMF is constructed as a sum of corrections to a classical potential, split into intra- and inter-molecular parts, and computed via regularized iterative Boltzmann inversion (IBI) on PIMD-derived distributions (Lawrence et al., 2023).

Hybrid CMD (h-CMD):

• Heterogeneous systems partitioned so that f-QCMD is applied to moieties with strong curvature artifacts (e.g., water), while more rigid or "classically" behaving subsystems use f-CMD.
• The total effective Hamiltonian is a direct sum of subsystem-specific potentials (London et al., 2024, Limbu et al., 2024).

Neural Network-Accelerated CMD (ML-CMD):

• A machine-learned force field is trained on centroid forces extracted from PIMD (Loose et al., 2022).

Partially-Adiabatic and Elevated-Temperature CMD:

• Internal and centroid modes are thermostatted at different temperatures; using an elevated internal mode temperature corrects flattening of the centroid PMF at low temperatures (“curvature problem”) (Castro et al., 22 Aug 2025).

Bead-Fourier CMD (BF-CMD):

• Imaginary-time trajectories are expanded in both discretized beads and a few Fourier sine modes, enabling efficient, accurate PMF evaluation with fewer degrees of freedom (London et al., 19 May 2025).

3. Path-Integral Curvature Problem and Solutions

A core limitation of traditional CMD is the so-called “curvature problem”: as temperature decreases or along large-amplitude curvilinear motions (e.g., bond stretching), the constraint of a centroid in Cartesian coordinates leads to an artificial flattening of the PMF. This results in spurious red shifts and broadening of vibrational bands, particularly for X–H stretching modes in H-bonded systems (London et al., 2024, Lawrence et al., 2023, Trenins et al., 2019).

Table: Comparison of Approaches to the Curvature Problem

Method	Curvature Correction	Computational Cost
Standard CMD	None (problematic)	High
f-CMD	None (problematic)	Low (after fit)
f-QCMD	Curvilinear coords	Moderate
h-CMD	Hybrid, per region	Moderate
PA-Tₑ-CMD	Elevated $T_e$ PMF	Similar to PIMD

Both f-QCMD and h-CMD, as implemented in DL_POLY Quantum 2.1, accurately suppress the curvature-induced red shift by employing curvilinear constraint spaces for critical subsets of the system, with negligible compromise of computational efficiency versus f-CMD (London et al., 2024, Lawrence et al., 2023, Limbu et al., 2024). For generic systems, PA-Tₑ-CMD (CMD with elevated internal mode temperatures) provides a curvature-free alternative, with the $T_e$ parameter empirically chosen to minimize artifacts (Castro et al., 22 Aug 2025).

4. Practical Computation of Dynamical Properties

The principal application of CMD is to the calculation of quantum real-time correlation functions—infrared (IR) absorption spectra, velocity autocorrelation functions, thermal conductivity, and neutron cross sections—where NQEs are essential.

IR Spectra:

• CMD and its variants use the Kubo-transformed dipole time-derivative autocorrelation function. Calculated as

$$\tilde{I}(\omega) = \frac{1}{2\pi} \int_{-\infty}^{\infty} dt\, e^{-i\omega t} \langle \dot{\boldsymbol{\mu}}(0)\cdot\dot{\boldsymbol{\mu}}(t) \rangle\,f(t)$$

with $f(t)$ a Hann window; no further quantum corrections are applied (London et al., 2024).

• f-QCMD/h-CMD methods yield band positions and widths for water and ice within 5–10 cm$^{-1}$ of experiment, resolving features (e.g., O–D doublets at interfaces) that are inaccessible to T-RPMD or classical MD (London et al., 2024, Limbu et al., 2024, Lawrence et al., 2023).

Transport Properties:

• CMD can be used to compute thermal diffusivity $a$ and thermal conductivity $\lambda$ by extracting time-correlation functions from centroid dynamics, combined with quantum heat capacity from PIMD (Sutherland et al., 2021).

Neutron Cross Sections:

• CMD-derived velocity autocorrelation functions serve as input to Gaussian-approximation formulas for the dynamic structure factor, yielding quantitative agreement with experimental neutron cross sections for quantum liquids (Guarini et al., 2015).

5. Formal Relation to Other Path-Integral Methods

CMD is a mean-field approximation to quantum Matsubara dynamics, exact for linear position or momentum operators in the harmonic limit, but neglecting non-centroid (fluctuational) modes that yield quantum coherence and higher-order effects (Hele et al., 2015, Jang, 2013). Relative to RPMD, CMD discards real-time propagation of the internal ring-polymer modes, thereby eliminating fictitious "ring-spring" artifacts but introducing the curvature problem in anharmonic systems (London et al., 2024). Mean-field Matsubara dynamics recovers the correct centroid PMF and $T$-independent vibrational peaks by partially reintroducing non-centroid fluctuations (Trenins et al., 2018).

Advances such as QCMD, h-CMD, and PA-Tₑ-CMD are motivated by the need to overcome the formal limitations of centroid-only sampling in strongly curved phase-space regions, demonstrating that strategic restoration of curvilinear centroids or Matsubara subspace averaging is essential for achieving quantum fidelity in vibrational spectra (Trenins et al., 2019, Lawrence et al., 2023, Limbu et al., 2024, Castro et al., 22 Aug 2025).

6. Implementation in Software and Computational Considerations

DL_POLY Quantum 2.1 provides modular, scalable, and MPI/OpenMP-parallelized implementations of f-CMD, f-QCMD, and h-CMD (London et al., 2024). Key workflow steps include:

• Reference PIMD trajectory generation with many beads (32–64) to sample centroid/quasi-centroid distributions.
• Fitting of effective centroid potentials by force-matching or IBI (usually 20–30 iterations, 50 ps trajectories).
• Production MD of classical centroids under these potentials with many independent replicas (20–50), yielding trajectories several times cheaper than on-the-fly adiabatic CMD for system sizes up to several hundred molecules.
• No thermostats are applied to centroid degrees of freedom in production runs.

The general cost scaling for fast centroid-based methods (f-CMD, f-QCMD, h-CMD) is only a small multiple of classical MD, enabling quantum-corrected IR spectra collection for complex heterogeneous and interfacial systems previously inaccessible to full path-integral dynamics (London et al., 2024, Lawrence et al., 2023).

7. Applications, Performance, and Scope

CMD and its modern centroid-based successors have established broad impact in the simulation of:

• Water and ice vibrational spectra, accurately capturing NQEs such as red-shifts in O–H stretching bands and temperature-dependent spectral features (London et al., 2024, Lawrence et al., 2023).
• Electrolyte solutions, where h-CMD uniquely reproduces interfacial spectral features such as the O–D doublet in concentrated Li-TFSI water-in-salt systems (London et al., 2024).
• Heterogeneous interfaces (e.g., water in MOFs), where h-CMD outperforms both T-RPMD (which over-broadens and obscures spectral resonance structure) and f-CMD (which misplaces and fails to resolve fine features) (Limbu et al., 2024).
• Bulk quantum fluids (para-hydrogen, helium), providing quantitative quantum corrections to thermal transport and neutron scattering cross sections (Sutherland et al., 2021, Guarini et al., 2015).
• Multi-region and multiscale simulations, leveraging adaptive resolution frameworks to combine quantum fidelity where needed with classical efficiency elsewhere (Agarwal et al., 2016, Kreis et al., 2017).

The methods are applicable to any system with sufficiently converged classical reference potentials; fitting of corrective potentials or neural network models is system-specific. For light-atom or proton-transfer systems, careful partitioning and curvilinear centroid construction are essential.

---

CMD and its fast, hybrid, and quasi-centroid extensions provide a rigorous and scalable approach to simulating quantum-dynamical observables in condensed matter, overcoming limitations of both classical MD and previously available path-integral schemes. These methods achieve a high degree of agreement with experiment for IR spectra, dynamical, and transport properties at quantum cost not much exceeding that of classical simulation, particularly when implemented with force-matching or machine-learning accelerations (London et al., 2024, Lawrence et al., 2023, Limbu et al., 2024, Loose et al., 2022).