{'fraction_non_alphanumeric': 0.03474903474903475, 'fraction_numerical': 0.01045981045981046, 'mean_word_length': 4.833742833742834, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Variational neural network ansatz for steady states in open quantum systems

We present a general variational approach to determine the steady state of open quantum lattice systems via a neural network approach. The steady-state density matrix of the lattice system is constructed via a purified neural network ansatz in an extended Hilbert space with ancillary degrees of freedom. The variational minimization of cost functions associated to the master equation can be performed using a Markov chain Monte Carlo sampling. As a first application and proof-of-principle, we apply the method to the dissipative quantum transverse Ising model.

We present a general variational approach to determine the steady state of open quantum lattice systems via a neural network approach. The steady-state density matrix of the lattice system is constructed via a purified neural network ansatz in an extended Hilbert space with ancillary degrees of freedom. The variational minimization of cost functions associated to the master equation can be performed using a Markov chain Monte Carlo sampling. As a first application and proof-of-principle, we apply the method to the dissipative quantum transverse Ising model.
In spite of the tremendous experimental progress in the isolation of quantum systems, a finite coupling to the environment [1] is unavoidable and certainly plays a crucial role in the practical implementation of quantum information and quantum simulation protocols [2]. Moreover, through an active control of the environment via the so-called reservoir engineering, an open quantum manybody system can be prepared in non-trivial phases [3][4][5] with also possible quantum applications [6,7]. The theoretical description of open quantum manybody systems is in general out-of-the equilibrium and much less developed than for equilibrium systems. A mixed state with a finite entropy can be described by a density matrix, whose evolution is described by a master equation.
Recently, a few theoretical methods have been developed to solve the master equation of open quantum manybody systems, including analytical approaches based on the Keldysh formalism [8,9], numerical algorithms based on matrix product operator and tensor-network techniques [10][11][12][13][14], cluster mean-field methods [15,16], corner-space renormalization [17][18][19], Gutzwiller mean-field [20], full configuration-interaction Monte Carlo [21]. The research in the field is very active, since the different methods are optimal for different specific regimes. For example, the corner-space renormalization method is best suited for systems with moderate entropy, while matrix product operator techniques to systems with short-range quantum correlations.
In the last decade, the field of artificial neural networks has enjoyed a dramatic expansion and success thanks to remarkable applications in the recognition of complex patterns such as visual images or human speech (for a recent review see, e.g., [22]). The optimization (supervised learning) of the network is obtained by tuning the weights quantifying the connections between neural units via a variational minimization of a properly defined cost function. The wavefunction of a manybody system is in general a complex quantity, which is hard to be recognized. Recent works have proposed to exploit artificial neural networks to construct trial wavefunctions, where the connection weights in the network play the role of variational parameters [23,24]. Neural network approaches have already been succesffuly applied to a wide number (see e.g. [25][26][27][28][29]) of close Hamiltonian systems. However, they have not yet been generalized to the important quantum manybody problem of open systems.
In this Letter, we present a theoretical approach based on a variational neural network ansatz in order to determine the steady state of the master equation of open quantum lattice systems. We construct the ansatz for the mixed density matrix starting from a Restricted Boltzman Machine ansatz for a pure many-body wavefunction in an extended Hilbert space. We determine the optimal variational parameters by minimizing a cost function which involves the Liouvillian superoperator associated to the master equation for the density matrix. As a first application, we have considered the dissipative tranverse field quantum Ising model. We present a proofof-principle demonstration by benchmarking the neural network calculations of the steady state against numerically exact simulations performed by quantum trajectories in the full Hilbert space [30]. Our minimization of the cost function is performed by Markov chain Monte Carlo sampling of the gradient and is thus scalable to a large number of lattice sites. Perspectives of the present approach are discussed in the conclusions.
The general task that we wish to solve is the determination of the steady state of an open quantum system described by the Lindblad master equation [1] for the system reduced density matrixρ, which reads (setting = 1): where L is the so-called Liouvillian superoperator depending on the system Hamiltonian operatorĤ. The coupling to the environment is represented by interaction channels with the reservoir characterized by dissipation rates γ j and jump operatorsL j acting on the system. The steady-state density matrix is defined as arXiv:1902.10104v1 [quant-ph] 26 Feb 2019 ρ ss = lim t→∞ρ (t) and can be shown to be the eigenmatrix corresponding to the 0−eigenvalue of L. Because of the exponential growth of the Hilbert space with the number of lattice sites, describing the full density matrix requires exponentially many complex numbers, which in practice can be done exactly only for a small number of sites. If one wants to attack the problem within a variational framework, the density matrix can be represented by an ansatzρ v depending on a set of variational parameters v. If {|σ = |σ 1 , σ 2 , ..., σ N } denotes a basis of states for the system Hilbert space, the density matrix can be expressed in the form In order to construct our neural network ansatz for the density matrix, we consider an extended Hilbert space H = H S ⊗ H A where H S,A represents respectively the system and ancillary Hilbert spaces. Such extended space is spanned by the basis set {|σ, a } where a = (a 1 , a 2 , ..., a Na ) labels the ancillary degrees of freedom. We start by considering a pure state in the extended Hilbert space, represented by the wavefunction ψ v (σ, a). In this framework the reduced density matrix of the system S is obtained by tracing out the ancillary degrees freedom [31], namely The next step is to represent ψ v (σ, a) via a neural network ansatz. This purified procedure automatically ensures thatρ v is Hermitian and positive semi-definite, as required for a density matrix. In a recent paper, Torlai and Melko [31] proposed to describe purified wavefunctions as Both the amplitude P v A (σ, a) and phase-related function P v θ (σ, a) of the purified wavefunction are given by the Boltmann-like expression where the associated dimensionless energy reads The rectangular matrix W ν weighs the connections between the system variables (visible layer) to the auxiliary variables (hidden layer), while the weight matrix U ν quantifies the connection between the system variables and the ancillary ones (ancillary layer). Such neural network ansatz is represented by a tri-partite Restricted Boltzmann Machine depicted in fig. 1. In other words, there are two independent artificial neural networks, one for the amplitude (ν = A) and one ν , Wν , Uν ). One network is used for the representation of the amplitude (ν = A) of the purified wavefunction in Eq. (3), while another independent network with the same topology is used to represent the phase (ν = θ).
for the phase (ν = θ). By substituting those formulas into Eq. (3) and carrying out the sum over the ancillary degrees of freedom one obtains a closed formula for the entries of the density matrix: where the expression of Γ +/− and Π can be found in Ref. [31]. It is worth pointing out that this scheme is not specific to this network topology, but relies only on the general fact that if two visible layers are connected by a shallow ancillary layer, the ancilla can be traced out analitically and an efficient neural-network description of the density matrix can be obtained.
Having defined a variational ansatzρ(v), we have now to define a variational principle to determine the optimal parameters. In particular, we have to recast the search for the steady state into a minimization problem with a positive-definite cost function C(v) which has a global minimum when the master equation Lρ = 0 is satisfied. Moreover, in order to be able to deal with large Hilbert spaces, we need a quantity which can be sampled and computed efficiently. These requirements are met by the following cost function expressed in terms of the 2-norm of the time derivative of the density-matrix: as (i) C(v ss ) = 0 ⇐⇒ρ(v ss ) =ρ ss and (ii) C(v) ≥ 0. It is useful to rewrite Eq. (7) as a sum over the whole space of bounded operators on the Hilbert space: where p v (σ,σ) = |ρ v (σ,σ)| 2 /Z corresponds to a probability distribution as Z = σ,σ |ρ v (σ,σ)| 2 . The local contribution reads: In order to find the global minimum of the cost function (7) by means of gradient-based iterative schemes we need to compute its gradient where we have defined the log-derivatives of the den- , which can be efficiently computed for the considered neural network.
The computational complexity of evaluating ∇C(v) exactly grows exponentially with the size of the system. This cost can be considerably reduced if one only uses an estimate of ∇C(v) obtained by sampling the values (σ,σ) according to the probability p v (σ,σ). Because the normalisation factor Z is not fixed, we cannot sample the distribution directly and have to resort to a Markov Chain Monte Carlo [32] with Metropolis update rules [33]. At every sampling step, we propose to update the configuration (σ,σ) → (σ ,σ ) by switching a random number of spins and accept the new configuration with probability min(exp[p v (σ,σ)/p v (σ ,σ )], 1).
Finally, in order to find the global minimum of the cost function we employ a standard Stochastic Gradient Descent algorithm [34]. In order to improve the performance of the Stochastic Gradient Descent (i.e. to reduce the number of iterations needed to converge to the global minima of the cost function) we update the variational parameters according to to the metric of the space of density matrices exploiting the Stochastic Reconfiguration Approach [35]. During the optimization procedure we sample the physical observables of interest through another Markov chain as Θ = Tr ρΘ beingσ α j the Pauli matrices (α ∈ {x, y, z}) acting on the j-th site. The first term represents the nearestneighbor spin-spin interaction depending only on the zcomponents, being V the coupling strength. The second term accounts for a local and uniform magnetic field along the transverse direction x. We consider local dissipative spin-flip processes described by the site-dependent jump operatorL , which fully determine the master equation in Eq. (1).
Numerical results for steady-state observables of the dissipative quantum transverse Ising model on a 1D periodic chain are reported in Fig. 2. In particular, we report the spatial components of the averaged magnetization as a function of the magnetic field g (in units of the dissipation rate γ) for V /γ = 2. For N = 16 lattice sites the predictions of the neural-network variational method (circles) are compared to the results obtained with a brute-force exact integration of the master equation in the whole Hilbert space, showing a remarkable agreement for g 1γ and g 2.5γ. Interestingly, for 1 g/γ 2.5, we note that the gradient-descent procedure requires more iterations. This region corresponds to the range of g/γ where the smallest nonzero eigenvalue of L † L decreases significantly [36]. Finally, as an example of convergence, the top panel of Fig. 3 depicts a typical evolution of the cost function in the iterative minimization procedure for a fixed set of parameters (g/γ = 1), showing a good convergence towards the global minimum. In the bottom panel of Fig. 3, the convergence of the x-component of the averaged magnetization is also reported.
In conclusion, we have presented a general variational approach for the steady-state density matrix of open quantum manybody systems based on an artificial neural network scheme. Our method is scalable since the cost function associated to the Liouvillian of the master equation can be calculated via Monte Carlo sampling. We have demonstrated a proof-of-principle of the theoretical scheme by a successful benchmarking to brute-force finite-size simulations in the full Hilbert space for arrays of spins described by the dissipative quantum transverse Ising model. We would like to point out that the present approach does not depend on the specific network topology. Indeed, the variational procedure presented in this Letter is general and can be applied to many other neural-networks or physically-inspired variational ansätze. There are many future developments at the horizon, including the study of dynamical properties, the use of deep neural networks and/or alternative cost functions, comparison with other existing techniques as well as the study of disordered systems without translational invariance. The neural network approach has the potential to pave the way to the theoretical study of a wide spectrum of open quantum manybody systems. We thank G. Carleo and V. Savona for discussions. We acknowledge support from ERC (via Consolidator Grant CORPHO No. 616233). This work was granted access to the HPC resources of TGCC under the allocation 2018-A0050510601 attributed by GENCI (Grand Equipement National de Calcul Intensif).
Note: after presenting this work at the conference "At the Crossroad of Physics and Machine Learning" (KITP, February 11, 2019), we were orally informed by participants of related independent theoretical works that have been carried on in parallel. Their results have been recently submitted to the arXiv [37][38][39].