"Scaling of contraction costs for entanglement renormalization algorithms including tensor Trotterization and variational Monte Carlo"
T. Barthel and Q. Miao arXiv:2407.21006,
pdf
The multi-scale entanglement renormalization ansatz (MERA) is a hierarchical class of tensor network states motivated by the real-space renormalization group. It is used to simulate strongly correlated quantum many-body systems. For prominent MERA structures in one and two spatial dimensions, we determine the optimal scaling of contraction costs as well as corresponding contraction sequences and algorithmic phase diagrams. This is motivated by recent efforts to employ MERA in hybrid quantum-classical algorithms, where the MERA tensors are Trotterized, i.e., chosen as circuits of quantum gates, and observables as well as energy gradients are evaluated by sampling causal-cone states. We investigate whether tensor Trotterization and/or variational Monte Carlo (VMC) sampling can lead to quantum-inspired classical MERA algorithms that perform better than the traditional optimization of full MERA based on the exact evaluation of energy gradients. Algorithmic phase diagrams indicate the best MERA method depending on the scaling of the energy accuracy and the number of Trotter steps with the bond dimension. The results suggest substantial gains due to VMC for two-dimensional systems.

"Equivalence of cost concentration and gradient vanishing for quantum circuits: An elementary proof in the Riemannian formulation"
Q. Miao and T. Barthel arXiv:2402.07883,
pdf, Quantum Sci. Technol. 9, 045039 (2024)
The optimization of quantum circuits can be hampered by a decay of average gradient amplitudes with increasing system size. When the decay is exponential, this is called the barren plateau problem. Considering explicit circuit parametrizations (in terms of rotation angles), it has been shown in Arrasmith et al (2022 Quantum Sci. Technol. 7 045015) that barren plateaus are equivalent to an exponential decay of the variance of cost-function differences. We show that the issue is particularly simple in the (parametrization-free) Riemannian formulation of such optimization problems and obtain a tighter bound for the cost-function variance. An elementary derivation shows that the single-gate variance of the cost function is strictly equal to half the variance of the Riemannian single-gate gradient, where we sample variable gates according to the uniform Haar measure. The total variances of the cost function and its gradient are then both bounded from above by the sum of single-gate variances and, conversely, bound single-gate variances from above. So, decays of gradients and cost-function variations go hand in hand, and barren plateau problems cannot be resolved by avoiding gradient-based in favor of gradient-free optimization methods.

"Driven-dissipative Bose-Einstein condensation and the upper critical dimension"
Y. Zhang and T. Barthel arXiv:2311.13561,
pdf, Phys. Rev. A 109, L021301 (2024)
Driving and dissipation can stabilize Bose-Einstein condensates. Using Keldysh field theory, we analyze this phenomenon for Markovian systems that can comprise on-site two-particle driving, on-site single-particle and two-particle loss, as well as edge-correlated pumping. Above the upper critical dimension, mean-field theory shows that pumping and two-particle driving induce condensation right at the boundary between the stable and unstable regions of the non-interacting theory. With nonzero two-particle driving, the condensate is gapped. This picture is consistent with the recent observation that, without symmetry constraints beyond invariance under single-particle basis transformations, all gapped quadratic bosonic Liouvillians belong to the same phase. For systems below the upper critical dimension, the edge-correlated pumping penalizes high-momentum fluctuations, rendering the theory renormalizable. We perform the one-loop renormalization group analysis, finding a condensation transition inside the unstable region of the non-interacting theory. Interestingly, its critical behavior is determined by a Wilson-Fisher-like fixed point with universal correlation-length exponent nu=0.6 in three dimensions.

"Criteria for Davies irreducibility of Markovian quantum dynamics"
Y. Zhang and T. Barthel arXiv:2310.17641,
pdf, J. Phys. A: Math. Theor. 57, 115301 (2024)
The dynamics of Markovian open quantum systems are described by Lindblad master equations, generating a quantum dynamical semigroup. An important concept for such systems is (Davies) irreducibility, i.e., the question whether there exist non-trivial invariant subspaces. Steady states of irreducible systems are unique and faithful, i.e., they have full rank. In the 1970s, Frigerio showed that a system is irreducible if the Lindblad operators span a self-adjoint set with trivial commutant. We discuss a more general and powerful algebraic criterion, showing that a system is irreducible if and only if the multiplicative algebra generated by the Lindblad operators $L_a$ and the operator $K=iH+\sum_a L^+_aL_a$, involving the Hamiltonian $H$, is the entire operator space. Examples for two-level systems, show that a change of Hamiltonian terms as well as the addition or removal of dissipators can render a reducible system irreducible and vice versa. Examples for many-body systems show that a large class of spin chains can be rendered irreducible by dissipators on just one or two sites. Additionally, we discuss the decisive differences between (Davies) reducibility and Evans reducibility for quantum channels and dynamical semigroups which has lead to some confusion in the recent physics literature, especially, in the context of boundary-driven systems. We give a criterion for quantum reducibility in terms of associated classical Markov processes and, lastly, discuss the relation of the main result to the stabilization of pure states and argue that systems with local Lindblad operators cannot stabilize pure Fermi-sea states.

"Machine learning with tree tensor networks, CP rank constraints, and tensor dropout"
H. Chen and T. Barthel arXiv:2305.19440,
pdf, IEEE Trans. Pattern Anal. Mach. Intell. 1 (2024)
Tensor networks approximate order-N tensors with a reduced number of degrees of freedom that is only polynomial in N and arranged as a network of partially contracted smaller tensors. As suggested in [arXiv:2205.15296] in the context of quantum many-body physics, computation costs can be further substantially reduced by imposing constraints on the canonical polyadic (CP) rank of the tensors in such networks. Here we demonstrate how tree tensor networks (TTN) with CP rank constraints and tensor dropout can be used in machine learning. The approach is found to outperform other tensor-network based methods in Fashion-MNIST image classification. A low-rank TTN classifier with branching ratio 4 reaches test set accuracy of 90.3 percent with low computation costs. Consisting of mostly linear elements, tensor network classifiers avoid the vanishing gradient problem of deep neural networks. The CP rank constraints have additional advantages: The number of parameters can be decreased and tuned more freely to control overfitting, improve generalization properties, and reduce computation costs. They allow us to employ trees with large branching ratios which substantially improves the representation power.

"Isometric tensor network optimization for extensive Hamiltonians is free of barren plateaus"
Q. Miao and T. Barthel arXiv:2304.14320,
pdf, Phys. Rev. A 109, L050402 (2024)
We explain why and numerically confirm that there are no barren plateaus in the energy optimization of isometric tensor network states (TNS) for extensive Hamiltonians with finite-range interactions. Specifically, we consider matrix product states, tree tensor network states, and the multiscale entanglement renormalization ansatz. The variance of the energy gradient, evaluated by taking the Haar average over the TNS tensors, has a leading system-size independent term and decreases according to a power law in the bond dimension. For a hierarchical TNS with branching ratio b, the variance of the gradient with respect to a tensor in layer t scales as (br)^t, where r is the second largest eigenvalue of the Haar-average doubled layer-transition channel and decreases algebraically with increasing bond dimension. The observed scaling properties of the gradient variance bear implications for efficient initialization procedures.

"Absence of barren plateaus and scaling of gradients in the energy optimization of isometric tensor network states"
T. Barthel and Q. Miao arXiv:2304.00161,
pdf
Vanishing gradients can pose substantial obstacles for high-dimensional optimization problems. Here we consider energy minimization problems for quantum many-body systems with extensive Hamiltonians and finite-range interactions, which can be studied on classical computers or in the form of variational quantum eigensolvers on quantum computers. Barren plateaus correspond to scenarios where the average amplitude of the energy gradient decreases exponentially with increasing system size. This occurs, for example, for quantum neural networks and for brickwall quantum circuits when the depth increases polynomially in the system size. Here we prove that the variational optimization problems for matrix product states, tree tensor networks, and the multiscale entanglement renormalization ansatz are free of barren plateaus. The derived scaling properties for the gradient variance provide an analytical guarantee for the trainability of randomly initialized tensor network states (TNS) and motivate certain initialization schemes. In a suitable representation, unitary tensors that parametrize the TNS are sampled according to the uniform Haar measure. We employ a Riemannian formulation of the gradient based optimizations which simplifies the analytical evaluation.

"Convergence and quantum advantage of Trotterized MERA for strongly-correlated systems"
Q. Miao and T. Barthel arXiv:2303.08910,
pdf
Strongly-correlated quantum many-body systems are difficult to study and simulate classically. We recently proposed a variational quantum eigensolver (VQE) based on the multiscale entanglement renormalization ansatz (MERA) with tensors constrained to certain Trotter circuits. Here, we extend the theoretical analysis, testing different initialization and convergence schemes, determining the scaling of computation costs for various critical spin models, and establishing a quantum advantage. For the Trotter circuits being composed of single-qubit and two-qubit rotations, it is experimentally advantageous to have small rotation angles. We find that the average angle amplitude can be reduced substantially with negligible effect on the energy accuracy. Benchmark simulations show that choosing TMERA tensors as brick-wall circuits or parallel random-pair circuits yields very similar energy accuracies.

"Tensor network states with low-rank tensors"
H. Chen and T. Barthel arXiv:2205.15296,
pdf
Tensor networks are used to efficiently approximate states of strongly-correlated quantum many-body systems. More generally, tensor network approximations may allow to reduce the costs for operating on an order-N tensor from exponential to polynomial in N, and this has become a popular approach for machine learning. We introduce the idea of imposing constraints on the canonical polyadic rank of the tensors that compose the tensor network. With this modification, the time and space complexities for the network optimization can be substantially reduced while maintaining high accuracy. We detail this idea for tree tensor network states (TTNS) and projected entangled-pair states. Simulations of spin models on Cayley trees with low-rank TTNS exemplify the effect of rank constraints on the expressive power. We find that choosing the tensor rank r to be on the order of the bond dimension m, is sufficient to obtain high-accuracy groundstate approximations and to substantially outperform standard TTNS computations. Thus low-rank tensor networks are a promising route for the simulation of quantum matter and machine learning on large data sets.

"Criticality and phase classification for quadratic open quantum many-body systems"
Y. Zhang and T. Barthel arXiv:2204.05346,
pdf, Phys. Rev. Lett. 129, 120401 (2022)
We study the steady states of translation-invariant open quantum many-body systems governed by Lindblad master equations, where the Hamiltonian is quadratic in the ladder operators, and the Lindblad operators are either linear or quadratic and Hermitian. These systems are called quasifree and quadratic, respectively. We find that steady states of one-dimensional systems with finite-range interactions necessarily have exponentially decaying Green's functions. For the quasifree case without quadratic Lindblad operators, we show that fermionic systems with finite-range interactions are noncritical for any number of spatial dimensions and provide bounds on the correlation lengths. Quasifree bosonic systems can be critical in D>1 dimensions. Last, we address the question of phase transitions in quadratic systems and find that, without symmetry constraints beyond invariance under single-particle basis and particle-hole transformations, all gapped Liouvillians belong to the same phase.

"Solving quasi-free and quadratic Lindblad master equations for open fermionic and bosonic systems"
T. Barthel and Y. Zhang arXiv:2112.08344,
pdf, J. Stat. Mech. 113101 (2022)
The dynamics of Markovian open quantum systems are described by Lindblad master equations. For fermionic and bosonic systems that are quasi-free, i.e., with Hamiltonians that are quadratic in the ladder operators and Lindblad operators that are linear in the ladder operators, we derive the equation of motion for the covariance matrix. This determines the evolution of Gaussian initial states and the steady states, which are also Gaussian. Using ladder super-operators (a.k.a. third quantization), we show how the Liouvillian can be transformed to a many-body Jordan normal form which also reveals the full many-body spectrum. Extending previous work by Prosen and Seligman, we treat fermionic and bosonic systems on equal footing with Majorana operators, shorten and complete some derivations, also address the odd-parity sector for fermions, give a criterion for the existence of bosonic steady states, cover non-diagonalizable Liouvillians also for bosons, and include quadratic systems. In extension of the quasi-free open systems, quadratic open systems comprise additional Hermitian Lindblad operators that are quadratic in the ladder operators. While Gaussian states may then evolve into non-Gaussian states, the Liouvillian can still be transformed to a useful block-triangular form, and the equations of motion for k-point Green's functions form a closed hierarchy. Based on this formalism, results on criticality and dissipative phase transitions in such models are discussed in a companion paper.

"Quantum-classical eigensolver using multiscale entanglement renormalization"
Q. Miao and T. Barthel arXiv:2108.13401,
pdf, Phys. Rev. Research 5, 033141 (2023)
We propose a variational quantum eigensolver (VQE) for the simulation of strongly correlated quantum matter based on a multiscale entanglement renormalization ansatz (MERA) and gradient-based optimization. This MERA quantum eigensolver can have substantially lower computation costs than corresponding classical algorithms. Due to its narrow causal cone, the algorithm can be implemented on noisy intermediate-scale quantum (NISQ) devices and still describe large systems. It is particularly attractive for ion-trap devices with ion-shuttling capabilities. The number of required qubits is system-size independent and increases only to a logarithmic scaling when using quantum amplitude estimation to speed up gradient evaluations. Translation invariance can be used to make computation costs square-logarithmic in the system size and describe the thermodynamic limit. We demonstrate the approach numerically for a MERA with Trotterized disentanglers and isometries. With a few Trotter steps, one recovers the accuracy of the full MERA.

"On the closedness and geometry of tensor network state sets"
T. Barthel, J. Lu, and G. Friesecke arXiv:2108.00031,
pdf, Lett. Math. Phys. 112, 72 (2022), Special Issue "Mathematical Physics and Numerical Simulation of Many-Particle Systems"Lett. Math. Phys. 2022
Tensor network states (TNS) are a powerful approach for the study of strongly correlated quantum matter. The curse of dimensionality is addressed by parametrizing the many-body state in terms of a network of partially contracted tensors. These tensors form a substantially reduced set of effective degrees of freedom. In practical algorithms, functionals like energy expectation values or overlaps are optimized over certain sets of TNS. Concerning algorithmic stability, it is important whether the considered sets are closed because, otherwise, the algorithms may approach a boundary point that is outside the TNS set and tensor elements diverge. We discuss the closedness and geometries of TNS sets, and we propose regularizations for optimization problems on non-closed TNS sets. We show that sets of matrix product states (MPS) with open boundary conditions, tree tensor network states, and the multiscale entanglement renormalization ansatz are always closed, whereas sets of translation-invariant MPS with periodic boundary conditions (PBC), heterogeneous MPS with PBC, and projected entangled pair states are generally not closed. The latter is done using explicit examples like the W state, states that we call two-domain states, and fine-grained versions thereof.

"Superoperator structures and no-go theorems for dissipative quantum phase transitions"
T. Barthel and Y. Zhang arXiv:2012.05505,
pdf, Phys. Rev. A 105, 052224 (2022)
In the thermodynamic limit, the steady states of open quantum many-body systems can undergo nonequilibrium phase transitions due to a competition between coherent and driven-dissipative dynamics. Here, we consider Markovian systems and elucidate structures of the Liouville superoperator that generates the time evolution. In many cases of interest, an operator-basis transformation can bring the Liouvillian into a block-triangular form, making it possible to assess its spectrum. The spectral gap sets the asymptotic decay rate. The superoperator structure can be used to bound gaps from below, showing that, in a large class of systems, dissipative phase transitions are actually impossible and that the convergence to steady states follows an exponential temporal decay. Furthermore, when the blocks on the diagonal are Hermitian, the Liouvillian spectra obey Weyl ordering relations. The results apply, for example, to Davies generators and quadratic systems and are also demonstrated for various spin models.

"Eigenstate entanglement scaling for critical interacting spin chains"
Q. Miao and T. Barthel arXiv:2010.07265,
pdf, Quantum 6, 642 (2022)
With increasing subsystem size and energy, bipartite entanglement entropies of energy eigenstates cross over from the groundstate scaling to a volume law. In previous work, we pointed out that, when strong or weak eigenstate thermalization (ETH) applies, the entanglement of all or, respectively, almost all eigenstates follow universal scaling functions which are determined by the subsystem entropy of thermal states. This was demonstrated by field-theoretical arguments and by analysis of large systems of non-interacting fermions and bosons. Here, we further substantiate such scaling properties for integrable and non-integrable interacting spin-1/2 chains at criticality using exact diagonalization. In particular, we analyze XXZ and transverse-field Ising models with and without next-nearest-neighbor interactions. We first confirm that the crossover for subsystem entropies in thermal ensembles can be described by a universal scaling function following from conformal field theory. Then, we analyze the validity of ETH for entanglement in these models. Even for the relatively small system sizes that can be simulated, the distributions of eigenstate entanglement entropies are sharply peaked around the subsystem entropies of the corresponding thermal ensembles.

"Low-energy physics of isotropic spin-1 chains in the critical and Haldane phases"
M. Binder and T. Barthel arXiv:2005.03643,
pdf, Phys. Rev. B 102, 014447 (2020)
Using a matrix product state algorithm with infinite boundary conditions, we compute high-resolution dynamic spin and quadrupolar structure factors in the thermodynamic limit to explore the low-energy excitations of isotropic bilinear-biquadratic spin-1 chains. Haldane mapped the spin-1 Heisenberg antiferromagnet to a continuum field theory, the nonlinear sigma model (NLSM). We find that the NLSM fails to capture the influence of the biquadratic term and provides only an unsatisfactory description of the Haldane phase physics. But several features in the Haldane phase can be explained by noninteracting multimagnon states. The physics at the Uimin-Lai-Sutherland point is characterized by multisoliton continua. Moving into the extended critical phase, we find that these excitation continua contract, which we explain using a field-theoretic description. New excitations emerge at higher energies and, in the vicinity of the purely biquadratic point, they show simple cosine dispersions. Using block fidelities, we identify them as elementary one-particle excitations and relate them to the integrable Temperley-Lieb chain.

"Scaling functions for eigenstate entanglement crossovers in harmonic lattices"
T. Barthel and Q. Miao arXiv:1912.10045,
pdf, Phys. Rev. A 104, 022414 (2021)
For quantum matter, eigenstate entanglement entropies obey an area law or log-area law at low energies and small subsystem sizes and cross over to volume laws for high energies and large subsystems. The crossover can be described by scaling functions. We demonstrate this for the harmonic lattice model which describes quantized lattice vibrations and is a regularization for free scalar field theories, describing, e.g., spin-0 bosonic particles. In one dimension, the groundstate entanglement obeys a log-area law. For dimensions d>1, it displays area laws, even at criticality, because excitation energies vanish only at a single point in momentum space. We sample excited states. The distribution of their entanglement entropies are sharply peaked around subsystem entropies of corresponding thermodynamic ensembles in accordance with the eigenstate thermalization hypothesis. In this way, we determine entanglement scaling functions which correspond to the quantum critical regime of the model. We show how infrared singularities of the system can be handled and how to access the thermodynamic limit using a perturbative trick for the covariance matrix. Eigenstates for quasi-free bosonic systems are not Gaussian. We resolve this problem by considering appropriate squeezed states instead. For these, entanglement entropies can be evaluated efficiently.

"Eigenstate entanglement: Crossover from the ground state to volume laws"
Q. Miao and T. Barthel arXiv:1905.07760,
pdf, Phys. Rev. Lett. 127, 040603 (2021)
For the typical quantum many-body systems that obey the eigenstate thermalization hypothesis (ETH) we argue that the entanglement entropy of (almost) all energy eigenstates is described by a single crossover function. The ETH implies that the crossover functions can be deduced from subsystem entropies of thermal ensembles and have universal properties. These functions capture the full crossover from the ground-state entanglement regime at low energies and small subsystem size (area or log-area law) to the extensive volume-law regime at high energies or large subsystem size. For critical one-dimensional systems, a universal scaling function follows from conformal field theory and can be adapted for nonlinear dispersions. We use it to also deduce the crossover scaling function for Fermi liquids in d>1 dimensions. The analytical results are complemented by numerics for large noninteracting systems of fermions in1

"The matrix product approximation for the dynamic cavity method"
T. Barthel arXiv:1904.03312,
pdf, J. Stat. Mech. 013217 (2020)
Stochastic dynamics of classical degrees of freedom, defined on vertices of locally tree-like graphs, can be studied in the framework of the dynamic cavity method which is exact for tree graphs. Such models correspond for example to spin-glass systems, Boolean networks, neural networks, and other technical, biological, and social networks. The central objects in the cavity method are edge messages - conditional probabilities of two vertex variable trajectories. In this paper, we discuss a rather pedagogical derivation for the dynamic cavity method, give a detailed account of the novel matrix product edge message (MPEM) algorithm for the solution of the dynamic cavity equation as introduced in Phys. Rev. E 97, 010104(R) (2018), and present optimizations and extensions. Matrix product approximations of the edge messages are constructed recursively in an iteration over time. Computation costs and precision can be tuned by controlling the matrix dimensions of the MPEM in truncations. Without truncations, the dynamics is exact. Data for Glauber-Ising dynamics shows a linear growth of computation costs in time. In contrast to Monte Carlo simulations, the approach has a much better error scaling. Hence, it gives for example access to low probability events and decaying observables like temporal correlations. We discuss optimized truncation schemes and an extension that allows to capture models which have a continuum time limit.

"Optimized Lie-Trotter-Suzuki decompositions for two and three non-commuting operators"
T. Barthel and Y. Zhang arXiv:1901.04974,
pdf, Ann. Phys. 418, 168165 (2020)
Lie-Trotter-Suzuki decompositions are an efficient way to approximate operator exponentials exp(tH) when H is a sum of n (non-commuting) terms which, individually, can be exponentiated easily. They are employed in time-evolution algorithms for tensor network states, digital quantum simulation protocols, path integral methods like quantum Monte Carlo, and splitting methods for symplectic integrators in classical Hamiltonian systems. We provide optimized decompositions up to sixth order. The leading error term is expanded in nested commutators (Hall bases) and we minimize the 1-norm of the coefficients. For two terms, several of the optima we find are close to those in [McLachlan, SlAM J. Sci. Comput. 16, 151 (1995)]. For three terms, our results substantially improve over unoptimized decompositions by Forest, Ruth, Yoshida, and Suzuki. We explain why these decompositions are sufficient to efficiently simulate one- and two-dimensional systems with finite-range interactions. This follows by solving a partitioning problem for the interaction graph.

"Infinite boundary conditions for response functions and limit cycles in iDMRG, demonstrated for bilinear-biquadratic spin-1 chains"
M. Binder and T. Barthel arXiv:1804.09163,
pdf, Phys. Rev. B 98, 235114 (2018)
Response functions <A(x,t)B(y,0)> for one-dimensional strongly correlated quantum many-body systems can be computed with matrix product state (MPS) techniques. Especially, when one is interested in spectral functions or dynamic structure factors of translation-invariant systems, the response for some range |x-y|<L is needed. We demonstrate how the number of required time-evolution runs can be reduced substantially: (a) If finite-system simulations are employed, the number of time-evolution runs can be reduced from L to 2sqrt(L). (b) To go beyond, one can employ infinite MPS (iMPS) such that two evolution runs suffice. To this purpose, iMPS that are heterogeneous only around the causal cone of the perturbation are evolved in time, i.e., the simulation is done with so-called infinite boundary conditions. Computing overlaps of these states, spatially shifted relative to each other, yields the response functions for all distances |x-y|. As a specific application, we compute the dynamic structure factor for ground states of bilinear-biquadratic spin-1 chains with very high resolution and explain the underlying low-energy physics. To determine the initial uniform iMPS for such simulations, infinite-system density-matrix renormalization group (iDMRG) can be employed. We discuss that, depending on the system and chosen bond dimension, iDMRG with a cell size n may converge to a non-trivial limit cycle of length m. This then corresponds to an iMPS with an enlarged unit cell of size mn.

"Fundamental limitations for measurements in quantum many-body systems"
T. Barthel and J. Lu arXiv:1802.04378,
pdf, Phys. Rev. Lett. 121, 080406 (2018)
Dynamical measurement schemes are an important tool for the investigation of quantum many-body systems, especially in the age of quantum simulation. Here, we address the question whether generic measurements can be implemented efficiently if we have access to a certain set of experimentally realizable measurements and can extend it through time evolution. For the latter, two scenarios are considered (a) evolution according to unitary circuits and (b) evolution due to Hamiltonians that we can control in a time-dependent fashion. We find that the time needed to realize a certain measurement to a predefined accuracy scales in general exponentially with the system size -- posing a fundamental limitation. The argument is based, on the construction of epsilon-packings for manifolds of observables with identical spectra and a comparison of their cardinalities to those of epsilon-coverings for quantum circuits and unitary time-evolution operators. The former is related to the study of Grassmann manifolds.

"Typical one-dimensional quantum systems at finite temperatures can be simulated efficiently on classical computers"
T. Barthel arXiv:1708.09349,
pdf
It is by now well-known that ground states of gapped one-dimensional (1d) quantum-many body systems with short-range interactions can be studied efficiently using classical computers and matrix product state techniques. A corresponding result for finite temperatures was missing. For 1d systems that can be described by a local 1+1d field theory, it is shown here that the cost for the classical simulation at finite temperatures grows in fact only polynomially with the inverse temperature and is system-size independent -- even for quantum critical systems. In particular, we show that the thermofield double state (TDS), a purification of the equilibrium density operator, can be obtained efficiently in matrix-product form. The argument is based on the scaling behavior of Rényi entanglement entropies in the TDS. At finite temperatures, they obey the area law. For quantum critical systems, the entanglement is found to grow only logarithmically with inverse temperature, S~log(beta). The field-theoretical results are confirmed by quasi-exact numerical simulations of quantum magnets and interacting bosons.

"Symmetric minimally entangled typical thermal states for canonical and grand-canonical ensembles"
M. Binder and T. Barthel arXiv:1701.03872,
pdf, Phys. Rev. B 95, 195148 (2017)
Based on the density matrix renormalization group (DMRG), strongly correlated quantum many-body systems at finite temperatures can be simulated by sampling over a certain class of pure matrix product states (MPS) called minimally entangled typical thermal states (METTS). When a system features symmetries, these can be utilized to substantially reduce MPS computation costs. It is conceptually straightforward to simulate canonical ensembles using symmetric METTS. In practice, it is important to alternate between different symmetric collapse bases to decrease autocorrelations in the Markov chain of METTS. To this purpose, we introduce symmetric Fourier and Haar-random block bases that are efficiently mixing. We also show how grand-canonical ensembles can be simulated efficiently with symmetric METTS. We demonstrate these approaches for spin-1/2 XXZ chains and discuss how the choice of the collapse bases influences autocorrelations as well as the distribution of measurement values and, hence, convergence speeds.

"Matrix product purifications for canonical ensembles and quantum number distributions"
T. Barthel arXiv:1607.01696,
pdf, Phys. Rev. B 94, 115157 (2016)
Matrix product purifications (MPPs) are a very efficient tool for the simulation of strongly correlated quantum many-body systems at finite temperatures. When a system features symmetries, these can be used to reduce computation costs substantially. It is straightforward to compute an MPP of a grand-canonical ensemble, also when symmetries are exploited. This paper provides and demonstrates methods for the efficient computation of MPPs of canonical ensembles under utilization of symmetries. Furthermore, we present a scheme for the evaluation of global quantum number distributions using matrix product density operators (MPDOs). We provide exact matrix product representations for canonical infinite-temperature states, and discuss how they can be constructed alternatively by applying matrix product operators to vacuum-type states or by using entangler Hamiltonians. A demonstration of the techniques for Heisenberg spin-1/2 chains explains why the difference in the energy densities of canonical and grand-canonical ensembles decays as 1/L.

"Finite-temperature effects on interacting bosonic one-dimensional systems in disordered lattices"
L. Gori, T. Barthel, A. Kumar, E. Lucioni, L. Tanzi, M. Inguscio, G. Modugno, T. Giamarchi, C. D'Errico, and G. Roux arXiv:1512.04238,
pdf, Phys. Rev. A 93, 033650 (2016)
We analyze the finite-temperature effects on the phase diagram describing the insulating properties of interacting 1D bosons in a quasi-periodic lattice. We examine thermal effects by comparing experimental results to exact diagonalization for small-sized systems and to density-matrix renormalization group (DMRG) computations. At weak interactions, we find short thermal correlation lengths, indicating a substantial impact of temperature on the system coherence. Conversely, at strong interactions, the obtained thermal correlation lengths are significantly larger than the localization length, and the quantum nature of the T=0 Bose glass phase is preserved up to a crossover temperature that depends on the disorder strength. Furthermore, in the absence of disorder, we show how quasi-exact finite-T DMRG computations, compared to experimental results, can be employed to estimate the temperature, which is not directly accessible in the experiment.

"A matrix product algorithm for stochastic dynamics on networks, applied to non-equilibrium Glauber dynamics"
T. Barthel, C. De Bacco, and S. Franz arXiv:1508.03295,
pdf, Phys. Rev. E 97, 010104(R) (2018)
We introduce and apply a novel efficient method for the precise simulation of stochastic dynamical processes on locally tree-like graphs. Networks with cycles are treated in the framework of the cavity method. Such models correspond, for example, to spin-glass systems, Boolean networks, neural networks, or other technological, biological, and social networks. Building upon ideas from quantum many-body theory, the new approach is based on a matrix product approximation of the so-called edge messages -- conditional probabilities of vertex variable trajectories. Computation costs and accuracy can be tuned by controlling the matrix dimensions of the matrix product edge messages (MPEM) in truncations. In contrast to Monte Carlo simulations, the algorithm has a better error scaling and works for both, single instances as well as the thermodynamic limit. We employ it to examine prototypical non-equilibrium Glauber dynamics in the kinetic Ising model. Because of the absence of cancellation effects, observables with small expectation values can be evaluated accurately, allowing for the study of decay processes and temporal correlations.

"Phase diagram of an extended quantum dimer model on the hexagonal lattice"
T. M. Schlittler, T. Barthel, G. Misguich, J. Vidal, and R. Mosseri arXiv:1507.04643,
pdf, Phys. Rev. Lett. 115, 217202 (2015)
We introduce a quantum dimer model on the hexagonal lattice that, in addition to the standard three-dimer kinetic and potential terms, includes a competing potential part counting dimer-free hexagons. The zero-temperature phase diagram is studied by means of quantum Monte Carlo simulations, supplemented by variational arguments. It reveals some new crystalline phases and a cascade of transitions with rapidly changing flux (tilt in the height language). We analyze perturbatively the vicinity of the Rokhsar-Kivelson point, showing that this model has the microscopic ingredients needed for the devil's staircase scenario [Fradkin et al. Phys. Rev. B 69, 224415 (2004)], and is therefore expected to produce fractal variations of the ground-state flux.

"Phase diagram of the hexagonal lattice quantum dimer model: Order parameters, ground-state energy, and gaps"
T. M. Schlittler, R. Mosseri, and T. Barthel arXiv:1501.02242,
pdf, Phys. Rev. B 96, 195142 (2017)
The phase diagram of the quantum dimer model on the hexagonal (honeycomb) lattice is computed numerically, extending on earlier work by Moessner, Sondhi, and Chandra. The different ground state phases are studied in detail using several local and global observables. In addition, we analyze imaginary-time correlation functions to determine ground state energies as well as gaps to the first excited states. This leads in particular to a confirmation that the intermediary so-called plaquette phase is gapped -- a point which was previously advocated with general arguments and some data for an order parameter, but required a more direct proof. On the technical side, we describe an efficient world-line quantum Monte Carlo algorithm with improved cluster updates that increase acceptance probabilities by taking account of potential terms of the Hamiltonian during the cluster construction. The Monte Carlo simulations are supplemented with variational computations.

"Minimally entangled typical thermal states versus matrix product purifications for the simulation of equilibrium states and time evolution"
M. Binder and T. Barthel arXiv:1411.3033,
pdf, Phys. Rev. B 92, 125119 (2015)
For the simulation of equilibrium states and finite-temperature response functions of strongly-correlated quantum many-body systems, we compare the efficiencies of two different approaches in the framework of the density matrix renormalization group (DMRG). The first is based on matrix product purifications. The second, more recent one, is based on so-called minimally entangled typical thermal states (METTS). For the latter, we highlight the interplay of statistical and DMRG truncation errors, discuss the use of self-averaging effects, and describe schemes for the computation of response functions. For critical as well as gapped phases of the spin-1/2 XXZ chain and the one-dimensional Bose-Hubbard model, we assess computation costs and accuracies of the two methods at different temperatures. For almost all considered cases, we find that, for the same computation cost, purifications yield more accurate results than METTS -- often by orders of magnitude. The METTS algorithm becomes more efficient only for temperatures well below the system's energy gap. The exponential growth of the computation cost in the evaluation of response functions limits the attainable timescales in both methods and we find that in this regard, METTS do not outperform purifications.

"Bound states and entanglement in the excited states of quantum spin chains"
J. Mölter, T. Barthel, U. Schollwöck, and V. Alba arXiv:1407.0066,
pdf, J. Stat. Mech. P10029 (2014), Special Issue "Quantum Entanglement in Condensed Matter Physics"J. Stat. Mech. 2014
We investigate the entanglement properties of the excited states of the spin-1/2 Heisenberg (XXX) chain with isotropic antiferromagnetic interactions, by exploiting the Bethe ansatz solution of the model. We consider eigenstates obtained from both real and complex solutions ('strings') of the Bethe equations. Physically, the former are states of interacting magnons, whereas the latter contain bound states of groups of particles. We first focus on the situation with few particles in the chain. Using exact results and semiclassical arguments, we derive an upper bound Smax for the entanglement entropy. This exhibits an intermediate behaviour between logarithmic and extensive, and it is saturated for highly-entangled states. As a function of the eigenstate energy, the entanglement entropy is organized in bands. Their number depends on the number of blocks of contiguous Bethe-Takahashi quantum numbers. In the presence of bound states a significant reduction in the entanglement entropy occurs, reflecting that a group of bound particles behaves effectively as a single particle. Interestingly, the associated entanglement spectrum shows edge-related levels. At a finite particle density, the semiclassical bound Smax becomes inaccurate. For highly-entangled states the subsystem entropy is proprtional to the chord length, signalling the crossover to extensive entanglement. Finally, we consider eigenstates containing a single pair of bound particles. No significant entanglement reduction occurs, in contrast with the few-particle case.

"Multispinon continua at zero and finite temperature in a near-ideal Heisenberg chain"
B. Lake, D. A. Tennant, J.-S. Caux, T. Barthel, U. Schollwöck, S. E. Nagler, and C. D. Frost arXiv:1307.4071,
pdf, Phys. Rev. Lett. 111, 137205 (2013)
The space- and time-dependent response of many-body quantum systems is the most informative aspect of their emergent behaviour. The dynamical structure factor, experimentally measurable using neutron scattering, can map this response in wavevector and energy with great detail, allowing theories to be quantitatively tested to high accuracy. Here, we present a comparison between neutron scattering measurements on the one-dimensional spin-1/2 Heisenberg antiferromagnet KCuF3, and recent state-of-the-art theoretical methods based on integrability and density matrix renormalization group simulations. The unprecedented quantitative agreement shows that precise descriptions of strongly correlated states at all distance, time and temperature scales are now possible, and highlights the need to apply these novel techniques to other problems in low-dimensional magnetism.

"Domain-wall melting in ultracold-boson systems with hole and spin-flip defects"
J. C. Halimeh, A. Wöllert, I. P. McCulloch, U. Schollwöck, and T. Barthel arXiv:1307.0513,
pdf, Phys. Rev. A 89, 063603 (2014)
Quantum magnetism is a fundamental phenomenon of nature. As of late, it has garnered a lot of interest because experiments with ultracold atomic gases in optical lattices could be used as a simulator for phenomena of magnetic systems. A paradigmatic example is the time evolution of a domain-wall state of a spin-1/2 Heisenberg chain, the so-called domain-wall melting. The model can be implemented by having two species of bosonic atoms with unity filling and strong on-site repulsion U in an optical lattice. In this paper, we study the domain-wall melting in such a setup on the basis of the time-dependent density matrix renormalization group (tDMRG). We are particularly interested in the effects of defects that originate from an imperfect preparation of the initial state. Typical defects are holes (empty sites) and flipped spins. We show that the dominating effects of holes on observables like the spatially resolved magnetization can be taken account of by a linear combination of spatially shifted observables from the clean case. For sufficiently large U, further effects due to holes become negligible. In contrast, the effects of spin flips are more severe as their dynamics occur on the same time scale as that of the domain-wall melting itself. It is hence advisable to avoid preparation schemes that are based on spin-flips.

"Algebraic versus exponential decoherence in dissipative many-particle systems"
Z. Cai and T. Barthel arXiv:1304.6890,
pdf, Phys. Rev. Lett. 111, 150403 (2013), covered by Phys.org "Quantum particles find safety in numbers"phys.org, 10/16/2013
The interplay between dissipation and internal interactions in quantum many-body systems gives rise to a wealth of novel phenomena. Here we investigate spin-1/2 chains with uniform local couplings to a Markovian environment using the time-dependent density matrix renormalization group (tDMRG). For the open XXZ model, we discover that the decoherence time diverges in the thermodynamic limit. The coherence decay is then algebraic instead of exponential. This is due to a vanishing gap in the spectrum of the corresponding Liouville superoperator and can be explained on the basis of a perturbative treatment. In contrast, decoherence in the open transverse-field Ising model is found to be always exponential. In this case, the internal interactions can both facilitate and impede the environment-induced decoherence.

"Precise evaluation of thermal response functions by optimized density matrix renormalization group schemes"
T. Barthel arXiv:1301.2246,
pdf, New J. Phys. 15, 073010 (2013)
This paper provides a study and discussion of earlier as well as novel more efficient schemes for the precise evaluation of finite-temperature response functions of strongly correlated quantum systems in the framework of the time-dependent density matrix renormalization group (tDMRG). The computational costs and bond dimensions as functions of time and temperature are examined at the example of the spin-1/2 XXZ Heisenberg chain in the critical XY phase and the gapped Néel phase. The matrix product state purifications occurring in the algorithms are in one-to-one relation with corresponding matrix product operators. This notational simplification elucidates implications of quasi-locality on the computational costs. Based on the observation that there is considerable freedom in designing efficient tDMRG schemes for the calculation of dynamical correlators at finite temperatures, a new class of optimizable schemes, as recently suggested in arXiv:1212.3570, is explained and analyzed numerically. A specific novel near-optimal scheme that requires no additional optimization reaches maximum times that are typically increased by a factor of two, when compared against earlier approaches. These increased reachable times make many more physical applications accessible. For each of the described tDMRG schemes, one can devise a corresponding transfer matrix renormalization group (TMRG) variant.

"Scaling of the thermal spectral function for quantum critical bosons in one dimension"
T. Barthel, U. Schollwöck, and S. Sachdev arXiv:1212.3570,
pdf
We present an improved scheme for the precise evaluation of finite-temperature response functions of strongly correlated systems in the framework of the time-dependent density matrix renormalization group. The maximum times that we can reach at finite temperatures T are typically increased by a factor of two, when compared against the earlier approaches. This novel scheme, complemented with linear prediction, allows us now to evaluate dynamic correlators for interacting bosons in one dimension. We demonstrate that the considered spectral function in the quantum critical regime with dynamic critical exponent z=2 is captured by the universal scaling form S(k,omega)=(1/T)*Phi(k/sqrt(T),omega/T) and calculate the scaling function precisely.

"Quasi-locality and efficient simulation of Markovian quantum dynamics"
T. Barthel and M. Kliesch arXiv:1111.4210,
pdf, Phys. Rev. Lett. 108, 230504 (2012)
We consider open many-body systems governed by a time-dependent quantum master equation with short-range interactions. With a generalized Lieb-Robinson bound, we show that the evolution in this very generic framework is quasilocal; i.e., the evolution of observables can be approximated by implementing the dynamics only in a vicinity of the observables' support. The precision increases exponentially with the diameter of the considered subsystem. Hence, time evolution can be simulated on classical computers with a cost that is independent of the system size. Providing error bounds for Trotter decompositions, we conclude that the simulation on a quantum computer is additionally efficient in time. For experiments and simulations in the Schrödinger picture, our result can be used to rigorously bound finite-size effects.

"Solving condensed matter ground state problems by semidefinite relaxations"
T. Barthel and R. Hübener arXiv:1106.4966,
pdf, Phys. Rev. Lett. 108, 200404 (2012)
We present a new generic approach to the condensed-matter ground-state problem which is complementary to variational techniques and works directly in the thermodynamic limit. Relaxing the ground-state problem, we obtain semidefinite programs (SDP). These can be solved efficiently, yielding strict lower bounds to the ground-state energy and approximations to the few-particle Green's functions. As the method is applicable for all particle statistics, it represents in particular a novel route for the study of strongly correlated fermionic and frustrated spin systems in D>1 spatial dimensions. It is demonstrated for the XXZ model and the Hubbard model of spinless fermions. The results are compared against exact solutions, quantum Monte Carlo, and Anderson bounds, showing the competitiveness of the SDP method.

"Dissipative quantum Church-Turing theorem"
M. Kliesch, T. Barthel, C. Gogolin, M. Kastoryano, and J. Eisert arXiv:1105.3986,
pdf, Phys. Rev. Lett. 107, 120501 (2011), see also D. Browne "Viewpoint: Quantum simulation hits the open road"Physics 4, 72 (2011)
We show that the time evolution of an open quantum system, described by a possibly time dependent Liouvillian, can be simulated by a unitary quantum circuit of a size scaling polynomially in the simulation time and the size of the system. An immediate consequence is that dissipative quantum computing is no more powerful than the unitary circuit model. Our result can be seen as a dissipative Church-Turing theorem, since it implies that under natural assumptions, such as weak coupling to an environment, the dynamics of an open quantum system can be simulated efficiently on a quantum computer. Formally, we introduce a Trotter decomposition for Liouvillian dynamics and give explicit error bounds. This constitutes a practical tool for numerical simulations, e.g., using matrix-product operators. We also demonstrate that most quantum states cannot be prepared efficiently.

"Real-space renormalization yields finitely correlated states"
T. Barthel, M. Kliesch, and J. Eisert arXiv:1003.2319,
pdf, Phys. Rev. Lett. 105, 010502 (2010)
Real-space renormalization approaches for quantum lattice systems generate certain hierarchical classes of states that are subsumed by the multi-scale entanglement renormalization ansatz (MERA). It is shown that, with the exception of one dimension, MERA states can be efficiently mapped to finitely correlated states, also known as projected entangled pair states (PEPS), with a bond dimension independent of the system size. Hence, MERA states form an efficiently contractible class of PEPS and obey an area law for the entanglement entropy. It is shown further that there exist other efficiently contractible schemes violating the area law.

"Contraction of fermionic operator circuits and the simulation of strongly correlated fermions"
T. Barthel, C. Pineda, and J. Eisert arXiv:0907.3689,
pdf, Phys. Rev. A 80, 042333 (2009), also in Virtual Journal of Nanoscale Science and Technology20, Issue 20 (2009), and in Virtual Journal of Quantum Information9, Issue 11 (2009)
A fermionic operator circuit is a product of fermionic operators of usually different and partially overlapping support. Further elements of fermionic operator circuits (FOCs) are partial traces and partial projections. The presented framework allows for the introduction of fermionic versions of known qudit operator circuits (QUOC), important for the simulation of strongly correlated d-dimensional systems: the multiscale entanglement renormalization ansaetze (MERA), tree tensor networks (TTN), projected entangled pair states (PEPS), or their infinite-size versions (iPEPS etc.). After the definition of a FOC, we present a method to contract it with the same computation and memory requirements as a corresponding QUOC, for which all fermionic operators are replaced by qudit operators of identical dimension. A given scheme for contracting the QUOC relates to an analogous scheme for the corresponding fermionic circuit, where additional marginal computational costs arise only from reordering of modes for operators occurring in intermediate stages of the contraction. Our result hence generalizes efficient schemes for the simulation of d-dimensional spin systems, as MERA, TTN, or PEPS to the fermionic case.

"Unitary circuits for strongly correlated fermions"
C. Pineda, T. Barthel, and J. Eisert arXiv:0905.0669,
pdf, Phys. Rev. A 81, 050303(R) (2010)
We introduce a scheme for efficiently describing pure states of strongly correlated fermions in higher dimensions using unitary circuits. A local way of computing local expectation values is presented. We formulate a dynamical reordering scheme, corresponding to time-adaptive Jordan-Wigner transformation that avoids non-local string operators and only keeps suitably ordered the causal cone. Primitives of such a reordering scheme are highlighted. Fermionic unitary circuits can be contracted with the same complexity as in the spin case. The scheme gives rise to a variational description of fermionic models that does not suffer from a sign problem. We present a numerical example on a 9x9 fermionic lattice model to show the functioning of the approach.

"Spectral functions in one-dimensional quantum systems at finite temperature using the density matrix renormalization group"
T. Barthel, U. Schollwöck, and S. R. White arXiv:0901.2342,
pdf, Phys. Rev. B 79, 245101 (2009)
We present for the first time time-dependent density-matrix renormalization-group simulations (t-DMRG) at finite temperatures. It is demonstrated how a combination of finite-temperature t-DMRG and time-series prediction allows for an easy and very accurate calculation of spectral functions in one-dimensional quantum systems, irrespective of their statistics, for arbitrary temperatures. This is illustrated with spin structure factors of XX and XXX spin-1/2 chains. For the XX model we can compare against an exact solution and for the XXX model (Heisenberg antiferromagnet) against a Bethe Ansatz solution and quantum Monte Carlo data.

"Magnetism, coherent many-particle dynamics, and relaxation with ultracold bosons in optical superlattices"
T. Barthel, C. Kasztelan, I. P. McCulloch, and U. Schollwöck arXiv:0809.5141,
pdf, Phys. Rev. A 79, 053627 (2009)
We study how well magnetic models can be implemented with ultracold bosonic atoms of two different hyperfine states in an optical superlattice. The system is captured by a two-species Bose-Hubbard model, but realizes in a certain parameter regime actually the physics of a spin-1/2 Heisenberg magnet, describing the second-order hopping processes. Tuning of the superlattice allows for controlling the effect of fast first-order processes versus the slower second-order ones. Using the density-matrix renormalization-group method, we provide the evolution of typical experimentally available observables. The validity of the description via the Heisenberg model, depending on the parameters of the Hubbard model, is studied numerically and analytically. The analysis is also motivated by recent experiments [S. Foelling et al., Nature (London) 448, 1029 (2007); S. Trotzky et al., Science 319, 295 (2008)] where coherent two-particle dynamics with ultracold bosonic atoms in isolated double wells were realized. We provide theoretical background for the next step, the observation of coherent many-particle dynamics after coupling the double wells. Contrary to the case of isolated double wells, relaxation of local observables can be observed. The tunability between the Bose-Hubbard model and the Heisenberg model in this setup could be used to study experimentally the differences in equilibration processes for nonintegrable and Bethe ansatz integrable models. We show that the relaxation in the Heisenberg model is connected to a phase averaging effect, which is in contrast to the typical scattering driven thermalization in nonintegrable models. We discuss the preparation of magnetic ground states by adiabatic tuning of the superlattice parameters.

"Quasiperiodic Bose-Hubbard model and localization in one-dimensional cold atomic gases"
G. Roux, T. Barthel, I. P. McCulloch, C. Kollath, U. Schollwöck, and T. Giamarchi arXiv:0802.3774,
pdf, Phys. Rev. A 78, 023628 (2008)
We compute the phase diagram of the one-dimensional Bose-Hubbard model with a quasiperiodic potential by means of the density-matrix renormalization group technique. This model describes the physics of cold atoms loaded in an optical lattice in the presence of a superlattice potential whose wavelength is incommensurate with the main lattice wavelength. After discussing the conditions under which the model can be realized experimentally, the study of the density vs the chemical potential curves for a nontrapped system unveils the existence of gapped phases at incommensurate densities interpreted as incommensurate charge-density-wave phases. Furthermore, a localization transition is known to occur above a critical value of the potential depth V_2 in the case of free and hard-core bosons. We extend these results to soft-core bosons for which the phase diagrams at fixed densities display new features compared with the phase diagrams known for random box distribution disorder. In particular, a direct transition from the superfluid phase to the Mott-insulating phase is found at finite V_2. Evidence for reentrances of the superfluid phase upon increasing interactions is presented. We finally comment on different ways to probe the emergent quantum phases and most importantly, the existence of a critical value for the localization transition. The latter feature can be investigated by looking at the expansion of the cloud after releasing the trap.

"Dephasing and the steady state in quantum many-particle systems"
T. Barthel and U. Schollwöck arXiv:0711.4896,
pdf, Phys. Rev. Lett. 100, 100601 (2008), also in Virtual Journal of Quantum Information8, Issue 3 (2008)
We discuss relaxation in bosonic and fermionic many-particle systems. For integrable systems, time evolution can cause a dephasing effect, leading for finite subsystems to steady states. We explicitly derive those steady subsystem states and devise sufficient prerequisites for the dephasing to occur. We also find simple scenarios, in which dephasing is ineffective and discuss the dependence on dimensionality and criticality. It follows further that, after a quench of system parameters, entanglement entropy will become extensive. This provides a way of creating strong entanglement in a controlled fashion.

"Entanglement entropy in collective models"
J. Vidal, S. Dusuel, and T. Barthel arXiv:cond-mat/0610833,
pdf, J. Stat. Mech. P01015 (2006)
We discuss the behaviour of the entanglement entropy of the ground state in various collective systems. Results for general quadratic two-mode boson models are given, yielding the relation between quantum phase transitions of the system (signalled by a divergence of the entanglement entropy) and the excitation energies. Such systems naturally arise when expanding collective spin Hamiltonians at leading order via the Holstein-Primakoff mapping. In a second step, we analyse several such models (the Dicke model, the two-level Bardeen-Cooper-Schrieffer model, the Lieb-Mattis model and the Lipkin-Meshkov-Glick model) and investigate the properties of the entanglement entropy over the whole parameter range. We show that when the system contains gapless excitations the entanglement entropy of the ground state diverges with increasing system size. We derive and classify the scaling behaviours that can be met.

"Entanglement entropy beyond the free case"
T. Barthel, S. Dusuel, and J. Vidal arXiv:cond-mat/0606436,
pdf, Phys. Rev. Lett. 97, 220402 (2006), also in Virtual Journal of Nanoscale Science and Technology14, Issue 24 (2006)
We present a perturbative method to compute the ground state entanglement entropy for interacting systems. We apply it to a collective model of mutually interacting spins in a magnetic field. At the quantum critical point, the entanglement entropy scales logarithmically with the subsystem size, the system size, and the anisotropy parameter. We determine the corresponding scaling prefactors and evaluate the leading finite-size correction to the entropy. Our analytical predictions are in perfect agreement with numerical results.

"Entanglement scaling in critical two-dimensional fermionic and bosonic systems"
T. Barthel, M.-C. Chung, and U. Schollwöck arXiv:cond-mat/0602077,
pdf, Phys. Rev. A 74, 022329 (2006)
We relate the reduced density matrices of quadratic fermionic and bosonic models to their Green's function matrices in a unified way and calculate the scaling of the entanglement entropy of finite systems in an infinite universe exactly. For critical fermionic two-dimensional (2D) systems at $T=0$, two regimes of scaling are identified: generically, we find a logarithmic correction to the area law with a prefactor dependence on the chemical potential that confirms earlier predictions based on the Widom conjecture. If, however, the Fermi surface of the critical system is zero-dimensional, then we find an area law with a sublogarithmic correction. For a critical bosonic 2D array of coupled oscillators at $T=0$, our results show that the entanglement entropy follows the area law without corrections.

"Entanglement and boundary critical phenomena"
H.-Q. Zhou, T. Barthel, J. O. Fjærestad, and U. Schollwöck arXiv:cond-mat/0511732,
pdf, Phys. Rev. A 74, 050305(R) (2006), also in Virtual Journal of Nanoscale Science and Technology14, Issue 23 (2006)
We investigate boundary critical phenomena from a quantum information perspective. Bipartite entanglement in the ground state of one-dimensional quantum systems is quantified using the Renyi entropy S_alpha, which includes the von Neumann entropy (alpha=1) and the single-copy entanglement (alpha=infinity) as special cases. We identify the contribution from the boundary entropy to the Renyi entropy, and show that there is an entanglement loss along boundary renormalization group (RG) flows. This property, which is intimately related to the Affleck-Ludwig g-theorem, can be regarded as a consequence of majorization relations between the spectra of the reduced density matrix along the boundary RG flows. We also point out that the bulk contribution to the single-copy entanglement is half of that to the von Neumann entropy, whereas the boundary contribution is the same.