Variational approach to open quantum systems with long-range competing interactions

Competition between short- and long-range interactions underpins many emergent phenomena in nature. Despite rapid progress in their experimental control, computational methods capable of accurately simulating open quantum many-body systems with complex long-ranged interactions at scale remain scarce. Here, we address this limitation by introducing an efficient and scalable approach to dissipative quantum lattices in one and two dimensions, combining matrix product operators and time-dependent variational Monte Carlo. We showcase the versatility, effectiveness, and unique methodological advantages of our algorithm by simulating the non-equilibrium dynamics and steady states of spin-$\frac{1}{2}$ lattices with competing algebraically-decaying interactions for as many as $N=200$ sites, revealing the emergence of spatially-modulated magnetic order far from equilibrium. This approach offers promising prospects for advancing our understanding of the complex non-equilibrium properties of a diverse variety of experimentally-realizable quantum systems with long-ranged interactions, including Rydberg atoms, ultracold dipolar molecules, and trapped ions.

💡 Research Summary

The manuscript introduces a novel computational framework for simulating open quantum many‑body systems with competing long‑range interactions. By marrying matrix‑product operators (MPO) with time‑dependent variational Monte Carlo (t‑VMC), the authors construct a variational manifold for the density matrix and project the exact Lindblad dynamics onto this manifold using the Dirac–Frenkel principle. The resulting equations of motion for the MPO parameters involve a metric tensor and a set of “variational forces,” both of which are expressed as expectation values over configurations sampled from the probability distribution defined by the current MPO. Monte Carlo sampling thus replaces the exponentially costly exact evaluation of inner products, while the MPO structure enables efficient contraction of the local Lindblad estimator even for highly non‑local terms.

Key technical innovations include: (i) a direct tensor‑contraction scheme for the local Lindblad estimator that avoids the auxiliary sampling required in neural‑network VMC; (ii) a decomposition of the Lindblad superoperator into a sum of n‑local components, allowing the algorithm to handle arbitrary algebraic decay (1/r^α) of interaction strengths; (iii) the use of a second‑order Heun integrator with adaptive step sizing, which together with a signal‑to‑noise regularization stabilizes the time evolution; and (iv) the complete elimination of Trotter errors because the method never discretizes the propagator e^{ℒt}. Consequently, the approach yields highly accurate steady‑state observables even at long times.

The authors benchmark the method against time‑evolving matrix‑product states (t‑MPS) with a second‑order Suzuki‑Trotter decomposition. For a 1D anisotropic Heisenberg chain with N = 200 sites and bond dimension χ = 20, the t‑VMC+MPO results for bulk magnetization, nearest‑ and next‑nearest‑neighbor spin‑spin correlations match the t‑MPS data across the entire evolution. In a separate test on a 4 × 4 2D square lattice, the same level of agreement is observed. Moreover, when evaluating relaxation to the steady state of a dissipative transverse‑field Ising chain, the stochastic t‑VMC+MPO dynamics (using a simple forward‑Euler step) achieves a smaller steady‑state error than t‑MPS even when the latter employs a four‑fold smaller time step, highlighting the advantage of avoiding Trotter accumulation.

The most striking application concerns systems with competing long‑range Ising interactions. The authors study a model where two power‑law couplings coexist: an antiferromagnetic dipolar term (J₁ < 0, α₁ = 3) and a ferromagnetic van‑der‑Waals term (J₂ > 0, α₂ = 6). Simulations of both 1D chains (N = 200) and 2D lattices (4 × 4) reveal that, after a transient, the non‑equilibrium steady state develops spatially modulated magnetic order—essentially a stripe‑like pattern of alternating magnetization. This demonstrates that frustration induced by competing algebraic interactions can survive in the presence of drive and dissipation, a phenomenon previously explored mainly in closed systems.

Scalability is a central claim: the computational cost scales with the MPO bond dimension χ and the number of Monte Carlo samples, not with the system size, enabling simulations of hundreds of spins in one dimension and modest two‑dimensional lattices. The method also works for periodic boundary conditions, which are often problematic for conventional tensor‑network time‑evolution algorithms.

In summary, the paper delivers a powerful, versatile, and scalable tool for tackling open quantum lattice models with arbitrary long‑range couplings. By combining the compactness of MPO representations with the stochastic efficiency of variational Monte Carlo, it overcomes the limitations of both Trotter‑based tensor‑network methods and neural‑network VMC approaches. The demonstrated ability to capture emergent spatial order in driven‑dissipative settings positions this technique as a promising avenue for future theoretical studies of Rydberg‑atom arrays, dipolar molecule platforms, trapped‑ion crystals, and other experimental systems where long‑range interactions and environmental coupling coexist.