PaperlayerView exact revision on arXiv

Source-linked AI summary

Quantum many-body dynamics in two dimensions with artificial neural networks

Markus Schmitt, Markus Heyl

arXiv:1912.08828v3cond-mat.str-elcond-mat.dis-nnquant-ph

TL;DR

Accurate real-time simulation of two-dimensional quantum dynamics remains challenging. This paper uses neural-network wave functions to evolve the transverse-field Ising model, reporting long time scales and collapse-and-revival oscillations comparable to or exceeding tensor-network capabilities.

  • Problem

    Accurately simulating post-quench dynamics can be difficult when small neural networks produce fluctuations across random initializations.

  • Method

    The approach encodes quantum many-body wave functions with a convolutional neural network and controls evolution error using adaptive time steps.

  • Results

    Ensemble averages over initial conditions showed very good agreement with iPEPS reference data for a quench to the critical point.

  • Takeaways & Limitations

    The reported critical-point dynamics are reliable when results are averaged over ensembles of initial conditions.

Abstract

from arXiv · show

The efficient numerical simulation of nonequilibrium real-time evolution in isolated quantum matter constitutes a key challenge for current computational methods. This holds in particular in the regime of two spatial dimensions, whose experimental exploration is currently pursued with strong efforts in quantum simulators. In this work we present a versatile and efficient machine learning inspired approach based on a recently introduced artificial neural network encoding of quantum many-body wave functions. We identify and resolve some key challenges for the simulation of time evolution, which previously imposed significant limitations on the accurate description of large systems and long-time dynamics. As a concrete example, we study the dynamics of the paradigmatic two-dimensional transverse field Ising model, as recently also realized experimentally in systems of Rydberg atoms. Calculating the nonequilibrium real-time evolution across a broad range of parameters, we, for instance, observe collapse and revival oscillations of ferromagnetic order and demonstrate that the reached time scales are comparable to or exceed the capabilities of state-of-the-art tensor network methods.

Supplemental material to Quantum many-body dynamics in two dimensions with artificial neural networks

The supplied supplemental-material excerpt contains bibliographic metadata rather than substantive scientific findings, identifying the authors and the arXiv version date.

  • The listed authors are Markus Schmitt and Markus Heyl.
  • The supplemental material corresponds to arXiv:1912.08828v3.

I. CONVOLUTIONAL NEURAL NETWORK WAVE FUNCTION

The CNN wave-function ansatz maps basis configurations to complex coefficients through alternating affine, convolutional, and nonlinear layers. Translation and point-symmetry structure, filter locality, and specialized activations are used to encode correlations and improve optimization.

  • Architecture: The CNN maps a basis configuration s to its wave-function coefficient ψη(s) through alternating nonlinear and affine maps arranged in layers.The linear affine component resembles a convolution of the preceding layer.
  • Architecture: Channels contain NT vertices associated with lattice-translation orbits, while filters connect activations across channels.The filters involve NF ≤ NT activations from each channel of the previous layer.
  • Parameterization: Complex filters and biases provide variational parameters capable of representing the complex-valued wave function.The parameter set is η ≡ (F, b), comprising filters F and biases b.
  • Symmetries: The ansatz can enforce translation and point symmetries of the wave function through the CNN construction.The text specifically targets invariance under translations and point symmetries π ∈ P.
  • Implementation: For multilayer CNNs, filters connect square dF × dF patches, with dF × L required to exceed the system’s linear extent for maximal-distance correlations.Single-layer CNNs are fully connected, whereas multilayer networks use local square patches.
  • Activations and optimization: An even polynomial activation incorporates the model’s Z2 symmetry by setting first-layer biases to zero, while later-layer derivatives help mitigate vanishing gradients.Lower-order polynomials substantially reduce network expressivity, motivating the chosen first-layer activation.

II. EXPLOITING LOCALITY AND CAUSALITY WITH THE CNN · III. TIME-DEPENDENT VARIATIONAL PRINCIPLE

The CNN architecture exploits locality by mediating long-distance correlations through deep layers, allowing more parameters to contribute when correlations remain local. The section also introduces the TDVP derivation from the Fubini–Study distance and examines noise in the resulting equation.

  • II. EXPLOITING LOCALITY AND CAUSALITY WITH THE CNN: The CNN mediates long-distance correlations through its deep layers.
  • II. EXPLOITING LOCALITY AND CAUSALITY WITH THE CNN: The correlation-growth analysis illustrates CNN and RBM parameters during a one-dimensional transverse-field Ising quench to criticality.
  • II. EXPLOITING LOCALITY AND CAUSALITY WITH THE CNN: Compared with a sparsely connected architecture, the CNN engages a larger fraction of variational parameters despite local correlations.
  • II. EXPLOITING LOCALITY AND CAUSALITY WITH THE CNN: This parameter utilization can make the CNN more efficient when correlations remain local.
  • III. TIME-DEPENDENT VARIATIONAL PRINCIPLE: The TDVP equation is reviewed through a derivation beginning from the Fubini–Study distance.
  • III. TIME-DEPENDENT VARIATIONAL PRINCIPLE: The TDVP discussion also examines the role of noise and includes supporting data.

A. TDVP from Fubini-Study distance

The TDVP is formulated by minimizing the short-time Fubini–Study distance between the variationally updated state and its unitary evolution, yielding a linear equation for the parameter velocities. The resulting Hamiltonian dynamics conserves energy under accurate force estimates and sufficiently small time steps, including hermitian regularized inversions.

  • A. TDVP from Fubini-Study distance: The TDVP minimizes the short-time Fubini–Study distance between the updated variational state and the state generated by unitary evolution.The distance is expanded consistently to second order in the short-time step, while first-order terms suffice in an intermediate expression because second-order contributions cancel.
  • A. TDVP from Fubini-Study distance: Stationarity of the expanded distance with respect to the parameter velocities yields the linear TDVP equation.The associated residual vector is R_k = S_k,k′ ˙η_k′ − F_k, and the error expression is obtained from a second-order expansion in the time step.
  • A. TDVP from Fubini-Study distance: The TDVP equation defines Hamiltonian dynamics on the variational manifold and therefore conserves the energy expectation value.The derivation uses the hermiticity of the S matrix.
  • A. TDVP from Fubini-Study distance: Accurate estimates of F_k and a sufficiently small time step τ preserve energy for any hermitian S^-1, including hermitian regularizations of the inversion.Regularization does not affect energy conservation provided hermiticity is maintained.

B. The role of noise for the time-dependent variational principle

The section analyzes Monte Carlo noise in the TDVP equation using a transformed eigenbasis and analytical signal-to-noise expressions under a Gaussian approximation. Simulations show relative Monte Carlo fluctuations remain consistent across orders of magnitude.

  • Noise in the S-matrix: Transforming the random variables into the eigenbasis of S provides a basis for analyzing Monte Carlo noise in the TDVP equation.The transformed variables have vanishing mean and covariance, with variance ⟨|Q_k|^2⟩ = σ_k^2.
  • Noise in the S-matrix: The analysis presents analytical signal-to-noise ratios for key TDVP quantities based on a Gaussian approximation, justified later in the section.The Gaussian assumption is expressed through ⟨|Q_k|^4⟩ = 3σ_k^4.
  • Noise in the S-matrix: Simulations show that relative Monte Carlo fluctuations are the same over all orders of magnitude.This observation is reported as consistent with the analytical noise analysis and illustrated in Fig. 2.

2. Noise in the F-vector … B. Network initialization

The method quantifies Monte Carlo noise in the F-vector using signal-to-noise estimates that remain accurate across a broad range, while motivating the Gaussian approximation heuristically. Numerical initialization adapts the computational basis, network depth, and ground-state search to obtain stable starting states.

  • 2. Noise in the F-vector: The signal-to-noise ratio depends on the index k through the ratio of σ_k^2 and |ρ_k|^2.
  • 2. Noise in the F-vector: 10% is the maximum relative deviation between the Gaussian-based and explicit Monte Carlo variance estimates, despite SNR(ρ_k) spanning four orders of magnitude.Equation (24) therefore provides a decent quantitative approximation under the Gaussian assumption.
  • 3. Joint distribution of Qk and ¯Eloc: The Gaussian behavior of Q_k and ¯E_loc is heuristically attributed to summing many essentially random terms, with the local energy comprising N identically distributed contributions.The argument invokes a central-limit-theorem mechanism rather than a strict proof.
  • A. Computational basis: The quantization axis follows the initial polarization: the z-basis is used for x-polarized states, and the x-basis for z-polarized states.This makes small-random-weight networks close to the desired initial states; the z-polarized state in the z-basis is computationally pathological because the S-matrix vanishes.
  • B. Network initialization: Single-layer networks use random weights drawn uniformly from [−w, w] with w = 10^-3.
  • B. Network initialization: Deep networks initialize F^(l)_{c,c′,k} uniformly on [−w^(l), w^(l)] with w^(l) = (N_F(α_{l−1} + α_l))^-1/2 to control gradient behavior.This choice also relies on odd activation functions in deep layers.
  • B. Network initialization: The initialized network undergoes Stochastic Reconfiguration ground-state search using a Hamiltonian tailored to the desired initial state.The search terminates when the energy variance density falls below 10^-7.

C. Adaptive time step · D. Monte Carlo sampling

The method combines a second-order Heun integrator with adaptive step-size control based on estimated integration errors and an S-matrix-induced norm. Monte Carlo sampling uses single-spin flips with periodic global flips to cover both Z2-related modes.

  • C. Adaptive time step: The simulations use a second-order consistent integration scheme, reducing Monte Carlo sampling requirements relative to higher-order adaptive methods such as Dormand–Prince.Errors are estimated by varying integration step sizes.
  • C. Adaptive time step: The Heun method is chosen as the second-order consistent integrator for evolving the variational equations.The method advances an ODE solution from y(t) using a step size τ.
  • C. Adaptive time step: Varying the step size τ enables estimation of the integration error, using the scheme’s O(τ^3) error scaling.The procedure compares one step of size τ with two steps of size τ/2.
  • C. Adaptive time step: The step size is adjusted to meet a desired tolerance ϵ, using the difference between solutions obtained with alternative step sizes.This difference provides the integration error δ used for adaptation.
  • C. Adaptive time step: Errors are measured with the norm induced by the S-matrix, weighting updates according to their significance for the physical state and normalizing by the parameter count P.The S-matrix serves as the metric tensor of the variational manifold.
  • C. Adaptive time step: Adaptive steps can vary over one order of magnitude, and selecting the maximal permissible value at each time substantially reduces computational cost.The time-step evolution is illustrated for a quench to h = hc/10 with α = (8; 10).
  • D. Monte Carlo sampling: Monte Carlo sampling uses simple Markov Chain Monte Carlo with single-spin flip updates to sample |ψ(s)|2.Every 200th proposed update is instead a global spin flip.
  • D. Monte Carlo sampling: Periodic global spin flips prevent sampling from being confined to only one of the two Z2-related modes of the distribution.A global flip is proposed every 200 updates.

E. Regularization scheme

The update vectors are regularized using the signal-to-noise ratio (SNR) of Monte Carlo estimates, with a tunable SNR cutoff and soft-cutoff treatment alongside an adaptive integrator. The SNR can be estimated from Monte Carlo samples, while Gaussian and binning-based alternatives address computational cost and sample dependence.

  • Regularization scheme: Update vectors are regularized based on the SNR of the Monte Carlo estimate of ρ_k.This regularization introduces the SNR cutoff λ_SNR as an additional simulation hyperparameter.
  • Regularization scheme: Soft cutoffs are preferred over hard cutoffs when combined with the adaptive integrator.The stated motivation is to avoid hard cutoff behavior during the regularized update.
  • SNR estimation: The SNR is estimated from a Monte Carlo sample of size N_MC.The supplied passage introduces this estimate but does not include the defining expression.
  • SNR estimation: When the Gaussian approximation is accurate, Eq. (24) provides a cheaper SNR estimate, although its higher-cost alternative was never computationally relevant in these simulations.The comparison concerns the cost of the two SNR-estimation expressions.
  • SNR estimation: Binning analysis can estimate fluctuations when independence of individual samples is doubtful, but this situation did not occur in the presented simulations.The authors therefore did not need the binning alternative for their reported simulations.

F. Ensemble averaging of initial conditions · G. Summary of simulation parameters · V. CONVERGENCE CHECKS

Small random-initialization ensembles can produce fluctuating critical-quench dynamics, but ensemble averaging agrees well with iPEPS reference data. The simulations used varied hyperparameters without requiring precise fine-tuning, alongside distinct regularization choices for different figures.

  • F. Ensemble averaging of initial conditions: Small network sizes produced noticeable fluctuations in post-quench dynamics across random initializations.This occurred despite the high accuracy of the initial ground-state search for quenches to the critical point.
  • F. Ensemble averaging of initial conditions: Ensemble averages of initial conditions showed very good agreement with iPEPS reference data.The agreement was demonstrated exemplarily in Fig. 4.
  • G. Summary of simulation parameters: The hyperparameters used for the Fig. 2 simulations are summarized in Table II.The table collects the various settings used for those simulations.
  • G. Summary of simulation parameters: Fig. 1 used no SNR-based regularization.Instead, Tikhonov regularization was combined with generalized cross-validation to determine an adaptive regularization parameter.
  • G. Summary of simulation parameters: Variation among the hyperparameters in Table II does not imply that precise fine-tuning is required.The authors state that the variations reflect experimentation with different settings.
  • G. Summary of simulation parameters: Different hyperparameter settings nonetheless led to consistent results.This consistency is presented as evidence that the observed variation reflects experimentation rather than a requirement for exact tuning.

A. Network size

The study tests simulation accuracy across neural-network sizes using fully connected single-layer CNNs. Increasing the size from α = 1 to α = 3 substantially reduces TDVP error and yields converged order-parameter results by α = 3.

  • A. Network size: The order-parameter results for α = 3 and α = 5 fully coincide, while Fisher information density shows small deviations.This comparison provides a systematic accuracy check across network sizes.
  • A. Network size: Fully connected single-layer CNNs with sizes α = 1, α = 3, and α = 5 are compared using observables and integrated TDVP error.The tested quantities include the order parameter, Fisher information density fQ(t), and integrated TDVP error R2(t).
  • A. Network size: Going from α = 1 to α = 3 substantially reduces the integrated TDVP error, with corresponding improvements in both observables.The smallest network produces inaccurate results and a much larger Fisher information density.

B. Finite size effect for quench to critical point

At the critical point, finite system size affects local-observable dynamics on relatively short timescales, while agreement with infinite-system iPEPS dynamics extends systematically as system size increases.

  • B. Finite size effect for quench to critical point: Finite-size effects affect local-observable dynamics on timescales exceeded by iPEPS when the system is smaller than N = 8 × 8.This is shown for transverse magnetization after quenching from a paramagnetic product state to h = hc across N = 6 × 6, N = 7 × 7, N = 8 × 8, and N = 10 × 10.
  • B. Finite size effect for quench to critical point: The time interval agreeing with infinite-system iPEPS dynamics extends systematically with increasing system size.The comparison concerns transverse magnetization after the quench to the critical point.

VI. COMPUTATIONAL COMPLEXITY AND PARALLEL COMPUTE PERFORMANCE

The time-evolution algorithm has overall complexity O(NMC × max(N^2, P) × P), with network evaluations dominating an example GPU-parallelized workload. Its computational structure supports hybrid distributed- and shared-memory parallelism, achieving perfect speedup with up to 256 MPI processes.

  • Computational complexity: When NMC > P, the overall complexity is O(NMC × max(N^2, P) × P), with most compute time spent on network evaluations.This distribution was observed using the GPU-parallelized algorithm.
  • Parallel implementation: Hybrid parallelization exploits independent Monte Carlo chains across compute nodes and shared-memory operations across CPU cores or GPUs.Batched network evaluation further improves arithmetic intensity, while global communication is required only once for summing matrix elements.
  • Parallel performance: Perfect speedup was achieved with up to 256 MPI processes when using N = 8×10^4 samples.The implementation used MPI for distributed memory, OpenMP/MKL for CPUs, and CUDA for GPU accelerators.
Loading 1912.08828v3…