Part I showed that replacing the exact quantum Liouvillian by a classical one is too blunt for nonlinear systems. Matsubara dynamics keeps a more structured classical limit: only the smooth low-frequency normal modes of the imaginary-time path are propagated.

From beads to modes

A discretized Kubo path can be viewed as a periodic imaginary-time loop. Expanding that loop in normal modes gives

\[ q(\tau)\approx \sum_{k=-K}^{K} Q_k \phi_k(\tau), \qquad M=2K+1. \]

The real basis used in the code is

\[ \phi_0(x)=1,\qquad \phi_k(x)=\sqrt 2\sin(2\pi kx),\qquad \phi_{-k}(x)=\sqrt 2\cos(2\pi kx), \]

where \(x=\tau/\beta\). The low modes change smoothly along imaginary time; the high modes contain rapid bead-to-bead roughness.

Matsubara Hamiltonian

After truncating to the \(M\) smooth modes, the effective Hamiltonian has the form

\[ H_M(P,Q)=\sum_{k=-K}^{K}\frac{P_k^2}{2m}+U_M(Q), \]

with an averaged path potential

\[ U_M(Q)=\frac{1}{\beta}\int_0^\beta V\!\left(\sum_{k=-K}^{K}Q_k\phi_k(\tau)\right)\,d\tau . \]

For the quartic benchmark, the code evaluates this integral by quadrature over \(q(\tau)\) and computes the gradient directly in mode space.

The phase term

The Matsubara distribution is not merely \(e^{-\beta H_M}\). A complex phase appears:

\[ \rho_M(P,Q)\propto e^{-\beta H_M(P,Q)} e^{i\beta\theta_M(P,Q)}, \qquad \theta_M(P,Q)=\sum_{k=-K}^{K}\omega_k Q_{-k}P_k . \]

This phase is the conserved quantity associated with imaginary-time translation symmetry. Keeping it is what lets Matsubara dynamics conserve the quantum Boltzmann distribution. It is also the source of the sign or phase problem that appears in practical calculations.

Workflow

  1. Generate a periodic path from many Fourier components.
  2. Keep only the \(M=11\) lowest Matsubara modes.
  3. Reconstruct the low-mode path and compare it with the full path.
  4. Compute the retained spectral power and RMS high-frequency component.

Result

Low-mode Matsubara filtering of an imaginary-time path
Matsubara low-mode filtering. Left: the horizontal axis is normalized imaginary time \(\tau/\beta\), and the vertical axis is the path value \(q(\tau)\). The blue solid curve uses all \(N=101\) Fourier modes, while the orange solid curve keeps only \(M=11\) low modes. Right: the horizontal axis is Fourier index \(k\), and the vertical log-scale axis is mode power. Dark bars are retained modes and pale bars are discarded modes. In this reproducible path, the retained power fraction is 0.573 and the RMS removed high-frequency component is 0.529, illustrating that the Matsubara subspace is a low-pass projection rather than an arbitrary bead subset.

Connection to CMD and RPMD

When only the centroid mode is retained, the construction points toward centroid molecular dynamics. When the phase is analytically continued and the resulting spring term is retained while the imaginary pieces are discarded, one recovers the logic behind RPMD-like dynamics. These are practical descendants, not exact replacements: they simplify the phase structure that makes full Matsubara dynamics hard.

Code used in this note

  • matsubara_mode_filter.py generates the Fourier path, applies the \(M=11\) low-mode truncation, and prints the retained-power diagnostics.

References

  • Hele, Willatt, Muolo, and Althorpe, J. Chem. Phys. 142, 134103 (2015).
  • Willatt, Matsubara Dynamics and its Practical Implementation, PhD thesis, University of Cambridge (2017).
  • Royer, J. Math. Phys. 25, 2873 (1984).

Matsubara series map

  1. Part I: From Kubo to LSC-IVR - phase-space representation and the harmonic/quartic benchmark.
  2. Part II: Low Modes and the Matsubara Phase - smooth imaginary-time paths, Matsubara truncation, and phase conservation.
  3. Part III: Quartic Benchmark and Phase Reweighting - practical estimator, phase problem, and noisy higher-mode convergence.
  4. Part IV: 绝热与非绝热 Matsubara 实时动力学推导 - parallel derivations with a cyclic electronic trace and explicit nuclear Gaussian integrals.