Symplectic Optimization on Gaussian States

TL;DR

The paper addresses the constrained nature of optimizing Gaussian ground states for quadratic bosonic Hamiltonians by introducing a globally unconstrained symplectic optimization (Sopt) framework. It parameterizes covariance matrices as positive-definite symplectic objects through unit-triangular factorizations, enforcing the uncertainty principle exactly and transforming the problem into an unconstrained trace-minimization with an explicit symplectic structure. The method yields accurate ground-state energies, covariance matrices, and spectral gaps for dipole-coupled quantum Drude oscillator lattices, with efficient scaling ($\mathcal{O}(d^{3})$ per step) and the ability to warm-start from nearby configurations to accelerate convergence. This unconstrained variational primitive is well suited to large-scale, iterative tasks and can be extended with tensor networks or perturbative non-quadratic terms for broader applications in materials science and quantum many-body physics. Key formulas include $E_0=\frac{1}{4}\mathrm{tr}(\mathcal{L}_3^{T}\mathcal{L}_3 H)$ and the condition $\gamma_0=(\mathcal{L}_3^{0})^{T}\mathcal{L}_3^{0}$, with the symplectic spectrum governed by the eigenvalues $\epsilon_i$ of $H$.

Abstract

Computing Gaussian ground states via variational optimization is challenging because the covariance matrices must satisfy the uncertainty principle, rendering constrained or Riemannian optimization costly, delicate, and thus difficult to scale, particularly in large and inhomogeneous systems. We introduce a symplectic optimization framework that addresses this challenge by parameterizing covariance matrices directly as positive-definite symplectic matrices using unit-triangular factorizations. This approach enforces all physical constraints exactly, yielding a globally unconstrained variational formulation of the bosonic ground-state problem. The unconstrained structure also naturally supports solution reuse across nearby Hamiltonians: warm-starting from previously optimized covariance matrices substantially reduces the number of optimization steps required for convergence in families of related configurations, as encountered in crystal lattices, molecular systems, and fluids. We demonstrate the method on weakly dipole-coupled lattices, recovering ground-state energies, covariance matrices, and spectral gaps accurately. The framework further provides a foundation for large-scale approximate treatments of weakly non-quadratic interactions and offers potential scaling advantages through tensor-network enhancements.

Figures (5)

• Figure 1: (a) Schematic of the Gaussian variational manifold. The ground-state CM $\gamma_{0}$ lies in the intersection of the sets of symmetric positive-definite and symplectic matrices. The symplectic optimization (Sopt) procedure performs unconstrained optimization over this manifold using unit-triangular factorization, starting from an initial CM $\gamma_{\mathrm{in}}$. (b--d) Validation of Gaussian ground-state recovery. (b) Agreement between the symplectic spectrum obtained via Sopt and that from symplectic diagonalization (SD) for a $d = 3(5\times5\times5)$ cubic lattice at $\rho = 1.9$. (c) Element-wise deviation between the CM retrieved via Sopt and that obtained from SD, shown for position correlations (left) and momentum correlations (right), displayed as log-scale heatmaps. (d) Convergence of the Sopt algorithm for cubic lattices of different sizes at $\rho = 1.9$, shown as the ground-state energy deviation $E^{\mathrm{Sopt}}_{0} - E^{\mathrm{SD}}_{0}$ versus the number of optimization steps $s$. Downward-pointing markers correspond to initialization from $\gamma_{\mathrm{in}}=\gamma_{\mathcal{T}}$, while upward-pointing markers indicate a warm-started optimization initialized from the ground-state CM of the same lattice at $\rho = 1.92$. This preconditioning reduces the number of optimization steps required to reach a given accuracy by approximately a factor of two. (e--f) Energy gap $\Delta$ obtained via Sopt for cubic and chain lattices at $\rho = 2$ and $\rho = 1.9$, respectively, compared with SD. (f) Residual gap error after $s = 400$ optimization steps as a function of the number of lattice modes $d$. The inset highlights the low-$d$ regime.
• Figure 2: (a) The graphical representation of a tensor $T^{\sigma_\nu,\sigma_\nu'}_{\alpha_{\nu},\alpha_{\nu+1}}$. (b) The graphical representation of the product of two matrices $A$ and $B$. The connection of two legs represents the contraction of shared matrix dimensions. (c) Tensor network diagram of the cost function $E_0$. (d) networks required to compute the gradient of $E_0$ with respect to a local tensor, in this case, $M_2(2)$. The region wrapped by the dotted line containing two holes (one with a red tensor) is the environment of the local tensor $M_2(2)$. The full gradient is computed by summing term shown where the red tensor is replaced by $M_2(2)000$, the term where the red tensor is replaced by $000\mathbb{I}$, and their transposes. Here, $\mathds{1}$ is an identity matrix of the appropriate dimension and $\mathbb{I}$ is the $2\times 2$ identity matrix.
• Figure 3: The ground state energy of the bosons on a $d=258\times258\times258$ cubic lattice, computed with the optimization of local tensors with $\chi=4$ using the L-BFGS algorithm (10 modes, with 10 iterations) zhuAlgorithm778LBFGSB1997. We set the potential energies to unity and assigned a weight of $0.1$ to the nearest-neighbor interaction. on the nearest neighbor interaction. The dashed line represents the trivial bound on the ground state energy, $d/2$. The optimization converged to within $\approx 0.8\%$ of the exact ground state energy $E^{\mathrm{SD}}_0\approx 8.3226\times10^{6}$.
• Figure 4: Utilizing the majorization relation in Eq. \ref{['sympstielf']}, we show here that Sopt can retrieve the first three symplectic eigenvalues of a 12 mode symplectic spectrum, where we numerically solve the optimization in Eq. \ref{['sympsteilfunitri']}, for different $k$ values.
• Figure 5: Sopt accurately retrieves ground-state energies on dipole coupled lattices, amended with diagonal position momentum coupling, with coupling strength given by $c$. All lattice calculations are shown at $\rho=1.9$, with $N=5\times 5\times 5$ in (a), $N=10\times 10$ in (b), $N=2\times 7\times 7$ in (c) and $N=3\times 5 \times 5$ in (d).

Theorems & Definitions (2)

• Lemma 1
• Lemma 2