We introduce the concept of directed loops in stochastic series expansion and path-integral quantum Monte Carlo methods. Using the detailed balance rules for directed loops, we show that it is possible to smoothly connect generally applicable simulation schemes (in which it is necessary to include backtracking processes in the loop construction) to more restricted loop algorithms that can be constructed only for a limited range of Hamiltonians (where backtracking can be avoided). The "algorithmic discontinuities" between general and special points (or regions) in parameter space can hence be eliminated. As a specific example, we consider the anisotropic S=1/2 Heisenberg antiferromagnet in an external magnetic field. We show that directed-loop simulations are very efficient for the full range of magnetic fields (zero to the saturation point) and anisotropies. In particular, for weak fields and anisotropies, the autocorrelations are significantly reduced relative to those of previous approaches. The back-tracking probability vanishes continuously as the isotropic Heisenberg point is approached. For the XY model, we show that back tracking can be avoided for all fields extending up to the saturation field. The method is hence particularly efficient in this case. We use directed-loop simulations to study the magnetization process in the two-dimensional Heisenberg model at very low temperatures. For LxL lattices with L up to 64, we utilize the step structure in the magnetization curve to extract gaps between different spin sectors. Finite-size scaling of the gaps gives an accurate estimate of the transverse susceptibility in the thermodynamic limit: chi( perpendicular )=0.0659+/-0.0002.

A cluster update (the "operator-loop") is developed within the framework of a numerically exact quantum Monte Carlo method based on the power series expansion of exp(−βH) (stochastic series expansion). The method is generally applicable to a wide class of lattice Hamiltonians for which the expansion is positive definite. For some important models the operator-loop algorithm is more efficient than loop updates previously developed for "worldline" simulations. The method is here tested on a two-dimensional anisotropic Heisenberg antiferromagnet in a magnetic field.The path integral formulation of quantum statistical mechanics is a useful starting point for numerical studies of interacting many-body systems in cases where positive-definiteness can be assured. Monte Carlo algorithms based on the "Trotter decomposition" 1,2 in discrete imaginary time, commonly referred to as "worldline" methods, have been used extensively for studies of quantum spins and bosons, as well as fermions in one dimension (in higher dimensions the fermion path integral is not positive definite). 3 Recently, two important technical developments have lead to significantly more efficient simulation algorithms. A generalization 4 of cluster updates used in classical Monte Carlo simulations 5 can reduce the autocorrelation times of some simulations by orders of magnitude, 6,7 thereby enabling studies of models in parameter regimes where standard local updating schemes do not efficiently explore the configuration space. Algorithms have also been constructed that work directly in the imaginary time continuum, 8-11 thus producing results free of systematic errors without the extrapolations to zero discretization which are required in order to obtain numerically exact results using the Trotter decomposition.There are, however, still unresolved issues for these improved algorithms. For some important models the loop schemes do not take into account all interactions in the system, and hence an a posteriori acceptance probability has to be assigned after the loop-clusters have been constructed. 7,12 This can seriously affect the efficiency of simulations. Some loop algorithms also break down due to "freezing", 4,13 when the probability is high for a single cluster to encompass the whole system. It is also often a highly non-trivial task to construct an algorithm for a new Hamiltonian -it would clearly be desirable to have a simple recipe for an arbitrary model.In this Communication, a general loop-type updating scheme is constructed within the "stochastic series expansion" (SSE) 8,14 framework. This approach to quantum simulations is based on sampling the diagonal matrix elements of the power series expansion of exp (−βH) [where H is the Hamiltonian and β the inverse temperature] and is related to a less general method proposed by Handscomb. 15 The SSE scheme is as general in applicability as the worldline method, and like the continuous time variant, it is numerically exact (there is also a strong relationship between the two methods 11 ). SSE algor...

The ground state parameters of the two-dimensional S = 1/2 antiferromagnetic Heisenberg model are calculated using the Stochastic Series Expansion quantum Monte Carlo method for L × L lattices with L up to 16. The finite-size results for the energy E, the sublattice magnetization M , the long-wavelength susceptibility χ ⊥ (q = 2π/L), and the spin stiffness ρs, are extrapolated to the thermodynamic limit using fits to polynomials in 1/L, constrained by scaling forms previously obtained from renormalization group calculations for the nonlinear σ model and chiral perturbation theory. The results are fully consistent with the predicted leading finite-size corrections and are of sufficient accuracy for extracting also subleading terms. The subleading energy correction (∼ 1/L 4 ) agrees with chiral perturbation theory to within a statistical error of a few percent, thus providing the first numerical confirmation of the finite-size scaling forms to this order. The extrapolated ground state energy per spin, E = −0.669437(5), is the most accurate estimate reported to date. The most accurate Green's function Monte Carlo (GFMC) result is slightly higher than this value, most likely due to a small systematic error originating from "population control" bias in GFMC. The other extrapolated parameters are M = 0.3070(3), ρs = 0.175(2), χ ⊥ = 0.0625(9), and the spinwave velocity c = 1.673(7). The statistical errors are comparable with those of the best previous estimates, obtained by fitting loop algorithm quantum Monte Carlo data to finite-temperature scaling forms. Both M and ρs obtained from the finite-T data are, however, a few error bars higher than the present estimates. It is argued that the T = 0 extrapolations performed here are less sensitive to effects of neglected higher-order corrections and therefore should be more reliable.

Using ground-state projector quantum Monte Carlo simulations in the valence-bond basis, it is demonstrated that nonfrustrating four-spin interactions can destroy the Néel order of the two-dimensional S=1/2 Heisenberg antiferromagnet and drive it into a valence-bond solid (VBS) phase. Results for spin and dimer correlations are consistent with a single continuous transition, and all data exhibit finite-size scaling with a single set of exponents, z=1, nu=0.78+/-0.03, and eta=0.26+/-0.03. The unusually large eta and an emergent U(1) symmetry, detected using VBS order parameter histograms, provide strong evidence for a deconfined quantum critical point.