Time Dependent Variational Principle for Tree Tensor Networks

We present a generalization of the Time Dependent Variational Principle (TDVP) to any finite sized loop-free tensor network. The major advantage of TDVP is that it can be employed as long as a representation of the Hamiltonian in the same tensor network structure that encodes the state is available. Often, such a representation can be found also for long-range terms in the Hamiltonian. As an application we use TDVP for the Fork Tensor Product States tensor network for multi-orbital Anderson impurity models. We demonstrate that TDVP allows to account for off-diagonal hybridizations in the bath which are relevant when spin-orbit coupling effects are important, or when distortions of the crystal lattice are present.


Introduction
The development of the Density Matrix Renormalization Group (DMRG) [1,2] was an immensely important milestone in our understanding of one-dimensional quantum systems.The subsequent realizations that DMRG produces Matrix Product States [3] (MPS) and that it can be formulated as a variational method [4], ultimately led to the development of numerous approaches using not only MPS but also general Tensor Networks to handle quantum systems.Notable examples are the Projected Entangled Pair States (PEPS) [5,6], the Multi-scale Entanglement Renormalization Ansatz (MERA) [7] and so-called Tree-Tensor Networks (TTN) [8][9][10][11][12][13][14][15][16][17] including also the recently developed Fork Tensor Product States (FTPS) method [18,19].Among the most important properties of tensor networks is whether their graph is loop-free, i.e., whether there exists only a single path from one tensor to any other.While PEPS and MERA are not loop-free, the TTNs and MPS are.Cutting any edge of a loop-free network, results in two separated segments and therefore gives a notion of left and right with respect to this edge.This in turn allows a controlled truncation scheme based on the Schmidtdecomposition of quantum states in the spirit of DMRG.One of the major reasons behind the success of tensor networks are the celebrated area laws of entanglement [20] stating that the entanglement of ground states of gapped Hamiltonians with short-range couplings is proportional to the surface area connecting the two regions.MPS in 1-d and PEPS in 2-d efficiently encode quantum states obeying these area laws and are hence efficient parametrizations.In addition, MPS-based time evolution for one dimensional systems is an important method to calculate dynamical properties [21][22][23][24].Approaches to perform the real-time evolution include, among others, the Time-dependent Density Matrix Renormalization Group (tDMRG) [25,26], the closely related Time Evolving Block Decimation (TEBD) [27,28] as well as the Time Dependent Variational Principle (TDVP) [29][30][31].An in depth comparison of several time evolution algorithms performed in Ref. [32] came to the conclusion that while all approaches have strengths and weaknesses, TDVP is among the most reliable methods to perform the time evolution.While time evolution approaches for MPS are well established, much less has been done for general tensor networks.So far, mostly TEBD (and variations) have been used, for example for the MERA network [33], for PEPS [34][35][36] and for TTNs [9,17,18,37,38].The advantage of TEBD is its relative simplicity, since it effectively boils down to a repeated application of short range operators obtained from a Suzuki-Trotter decomposition [39] of the full timeevolution operator.However, one of the major disadvantages of TEBD is that it can become difficult to implement for more complicated Hamiltonians, especially when long-range couplings are present.One approach to treat such couplings is an MPO-based approach introduced by Zaletel et al. [44] in which a MPO approximation of the time-evolution operator is constructed.Alternatively, TDVP circumvents this problem by only demanding a Hamiltonian represented in the same tensor network structure as the state which is often easy to find.Additionally, TDVP in its single-site variant exactly respects conserved quantities of the Hamiltonian like energy or magnetization [30].Although some works applied TDVP to more general tensor networks [45,46], it is not obvious how these algorithms work in detail and how it can be generalized.A notable exception is Ref. [47] which introduces TDVP for binary TTNs.Parallel to these developments in the tensor network physics community, very similar approaches to TDVP have been developed in quantum chemistry under the name of Multi-layer Multi-Configurational Time-Dependent Hartree approach [40][41][42][43].These methods effectively generate tensor networks by repeatedly grouping degrees of freedom together and transforming them with (time dependent) basis transformations into new degrees of freedom.A more practical motivation for the formulation of TDVP for TTNs are Dynamical Mean-Field Theory (DMFT) calculations using the FTPS tensor network.So far, this approach has been used for so-called diagonal hybridizations only.On the other hand, real materials often exhibit off-diagonal hybridizations, which can for example come from spin-orbit coupling, or from distortions of the crystal lattice.For off-diagonal hybridizations, the TEBD approach used so far [18,19] is difficult to generalize and we hence choose to use TDVP in these situations.Although part of the motivation for this work comes from the FTPS tensor network, in this paper we formulate TDVP for general loop-free and finite-size tensor networks.After establishing the relevant concepts of TTNs in Sec. 2, we generalize TDVP to these networks in Sec. 3. Finally in Sec. 4 we show how this approach can be used for the FTPS tensor network and that it can be applied to off-diagonal hybridizations.6).Cutting the link between sites 4 and 5, as indicated by the dashed line, results in two disconnected tensor network segments and defines a notion of left and right at each link.In this example, sites 1, 2, 3, 6 and 7 are the leaves of the TTN.

Tree Tensor Networks Basics
In this section, we discuss concepts of TTNs relevant for the formulation of TDVP.All these properties are generalizations of the corresponding concepts for MPS.Although these have been discussed previously in several publications (see for example Refs.[14,17,48]), here, we present them in a format that will suit us for the subsequent formulation of the TDVP algorithm.

TTNs
x Any state |ψ of a quantum system consisting of N sites with local basis states |s i on site i can be expanded in the corresponding product basis: The coefficient c Tensor networks represent this rank-N tensor as a product over tensors of much smaller rank: Each tensor has a set of auxiliary indices Q i = {q k : q k is attached to node i } such that each auxiliary index is part of exactly two tensors.We call r i = |Q i | the number Figure 3: With an SVD, we can orthogonalize a tensor towards one of its neighbors with which it shares an index.Top: Tensor T is reshaped into a matrix and the U -matrix of its SVD is used as new tensor.S • V † is absorbed into the neighboring tensor.Middle: Graphical representation of U † •U = 1.Bottom: A tensor that is normalized towards one of its neighbors is depicted as a triangle pointing in the direction of this neighbor.The picture shows tensor 4 of the TTN in Fig. 1 orthogonalized towards tensor 2. Let us call the index connecting tensor 4 with tensor 2, q 2 .In this case, we denote tensor 4 as of indices of the tensor on site i.Additionally, we attached to each tensor a physical index as for example in the FTPS tensor network.While for general TTNs not all tensors have a physical index, the following results can be straightforwardly generalized by just removing the physical index from the notation.Alternatively, every tensor without a physical index could be interpreted as having a dummy index with just a single entry corresponding to a single state, say |0 , onto which the Hamiltonian acts as an identity H |0 = |0 .Note that if all sites have a physical index, the number of links is L = N − 1.In the following, we will often omit sums over auxiliary indices q 1 •••q L and assume Einstein convention for the summations.
An example for a TTN with N = 7 sites and L = 6 auxiliary indices (links) is shown in Fig. 1.The property distinguishing a TTN from a general tensor network is that the graph of a TTN is loop-free, i.e., to move from one site to any other there is only one unique path along the links.This also implies that by cutting any link, the tensor networks splits into two disconnected segments.Therefore, at each link there is a notion of left and right which is a first hint towards the capability of TTNs to access the Schmidt decomposition and with it also the reduced density matrix as demonstrated below.We also define the leaves of the TTN as all tensors with just a single link index.For convenience, we assume site N to be a leave of the TTN.Since TTNs are loop-free, one can also define a measure of distance d ij between two sites i and j given by the number of links one has to traverse to move from site i to site j.

Tensor Gauge and Orthogonality Center
The representation of a quantum state as a tensor network is highly non-unique.This gauge degree of freedom can be used to obtain useful representations of the same quantum state as a TTN with certain properties, which can speed up calculations dramatically.As shown in Fig. 2, at each link one can insert an identity 1 = G • G −1 for any invertible matrix G.By absorbing G into one tensor and G −1 into the other, a different representation of the same state is reached.In this part, we make use of this gauge degree of freedom to define an orthogonality center of the TTN.A tensor T s i Q i can be orthogonalized towards one of its neighbors with which it shares link q k as follows: with rows (s i , Q i \ q k ) and column (q k ).
• Perform an SVD (a QR decomposition is faster): as the new local tensor on site i and absorb S • V † into the corresponding neighbor by multiplying S • V † onto it (formally also relabel α → q k ).
The SVD as well as the QR decomposition guarantees that the new site tensor has the property (see Fig. 3) For tensors orthogonalized towards their neighbor along link q k we introduce the notation (see Fig. 3 bottom).As already mentioned, in TTNs there is a unique path between any two tensors.Therefore, by orthogonalizing a tensor towards one of its neighbors, we also orthogonalize it towards all other tensors, which can be reached via this neighbor.For example, to orthogonalize tensors 1, 2 and 6 in Fig. 1 towards tensor 3, we orthogonalize all of them towards tensor 4 with the procedure described above.Next, let us introduce orthogonality centers.Site i is an orthogonality center with tensor C s i Q i if all tensors of all other sites are orthogonalized towards site i.To obtain such an orthogonality center, we can use the following algorithm: 1. Find the maximum distance d max between site i and any other site in the TTN.
2. Initialize d = d max and perform the following steps until d = 0 • Orthogonalize all sites j that are at distance d from site i towards site i, i.e., towards the single neighbor on the path from j to i.
• Reduce d by one d → d − 1.  1 is placed on site 4, the center tensor C has four link indices q 1 • • • q 4 .Each of these links corresponds to one of four mutually orthogonal set of states |q 4 .This orthogonality is a direct result of the orthogonality property of the U -matrices of the SVD used on all sites except site 4 (see also Fig. 3).The inset is a reminder of the TTN used in this section.
For example, to orthogonalize the TTN shown in Fig. 1 towards site 4, we first orthogonalize sites 3 and 7 towards site 5 and then sites 1, 2, 6 and 5 towards site 4.
The wave function of a TTN with orthogonality center C s i Q i can be written as: Here, S i q k is the segment of the tensor network that is obtained by cutting index q k and which does not contain site i.The states |q k form an orthogonal basis and are defined in Fig. 4 for the TTN of Fig. 1 with orthogonality center on site i = 4. Orthogonality centers allow to easily calculate local observables acting on the orthogonality center.For example, the expectation value of the operator Â = s i ,s i A s i ,s i |s i s i | acting non-trivially only on site i, reduces to: where the bar denotes complex conjugation.Orthogonality centers hence reduce the costly contraction over the whole tensor network, to a simple contraction over the center tensor C only.
Figure 5: Definition of the vertical subspace for site 4 in Eq. 9 with end point site 7.The space of all tensors defining the kernel of the map from the tangent space to the physical hilbert space, is spanned by a matrix X q l q l (yellow squares) for each link of the TTN.The physical index s 4 is labeled separately to distinguish it from the bond indices.The inset is a reminder of the TTN used in this section.

Truncation of TTNs
A TTN with an orthogonality center allows to calculate any Schmidt decomposition of the quantum state with respect to the two parts of the system defined by cutting any of the links of the orthogonality center.To do so, we reshape the center tensor C into a matrix with physical index s i and one of the links q k , combined into the row index and all other indices into the column indices, i.e., C s i q 1 ,q 2 ,••• ,qr i = C (s i ,q k ),(Q i \q k ) .The Schmidt decomposition then follows from an SVD of this matrix: Note that if one is interested solely in the truncation of a given orthogonality center and keeping it at the same site, an even more efficient approach would be to perform an SVD on the matrix C (q k ),(s i Q i \q k ) .Again, the orthogonality of the states |q k and the orthogonality of the U and V † matrices guarantee that the left and right vectors also form an orthogonal basis and hence Eq. 5 is a true Schmidt decomposition.Note that this Schmidt decomposition separates all sites in segment S i q k as well as site i to the rest of the lattice.From there it is straightforward to calculate the reduced density matrix for one of these two subsystems and approximate states by keeping only the largest eigenvalues in the spirit of DMRG.

TDVP equations for Tree Tensor Networks
In this section, we generalize the derivation of the tangent space projector presented for MPS in Ref. [29] to general TTNs.While the overall approach is very similar to the derivation for MPS, the lack of a clear start and end point in the tensor network geometry will make the derivation and the subsequent integration of the equations quite different from standard MPS.TDVP amounts to the solution of the Schrödinger equation in the space spanned by the tensor network without ever leaving this manifold (at least in its single-site variant).In TDVP, one solves a modified Schrödinger equation by projecting its right-hand side onto the so-called tangent space: In the following, we want to find a representation of the tangent space projection operator , which not only depends on the current state |ψ[T ] but importantly also on the structure of the TTN.

Tangent Space Projector
Any element of the tangent space |Θ[B] is parametrized by a set of tensors Importantly, for each summand we use the representation of the state |ψ[T ] in which site i is the orthogonality center, i.e., all tensors T s j Q j are orthogonalized towards site i such that The gauge degree of freedom in the TTN, reflects itself in the tangent space that not all linearly independent choices of B s i Q i result in different tangent vectors.Ref. [29] solves this problem by first defining the so-called vertical subspace, i.e., all tensors B s i Q i that give the zero-state |Θ[B] = 0 and hence define the kernel of the map from the tensors to the physical Hilbert space.Then, imposing a gauge prescription, they fix this kernel to a single element which guarantees that the resulting parametrization is unique.In order to arrive at a result that resembles the MPS algorithm, we first need to define a fixed end point of the TTN with the restriction that it should be a leave.Note however that any site of the tensor network can be used as end point.Without loss of generality, we choose site N as end point.The vertical subspace, i.e., all tensors B s i Q i for which |Θ[B] = 0 can then be parametrized by matrices X q k q k such that: (T N [q l ] ) s i Q i is the unique tensor of the state |ψ[T ] with site i orthogonalized towards the neighbor on the other end of the link q l .This definition of the vertical subspace is depicted in Fig. 5 for the tensor B s 4 Q 4 .The factor sgn(q l → N ) is 1 if link q l points towards the end point and −1 otherwise.This construction guarantees that for any choice of B s i Q i in the vertical subspace, |Θ[B] = 0, because the single term with positive sign (q l → N ) is exactly canceled by one negative term of its neighbor (since there (q l N )).Note that this definition of the vertical subspace reduces in the case of MPS to the definition used in Ref. [29] if the right-most site of the MPS is chosen as the end point.To uniquely specify the kernel, we impose the following matrix-valued (with indices q k and q k ) gauge fixing condition for the B-tensors of the tangent space: Again, the bar denotes complex conjugation.Above, q k is the single index pointing towards the end point N .These are N − 1 matrix-valued constraints, for the X-matrices living on L = N − 1 indices.This implies that no ambiguity is left in the definition of the kernel, if we choose B-tensors according to Eq. 10.
It also guarantees that the overlap between two tangent vectors reduces to a contraction over local tensors only: under the constraints given by Eq. 10.With Eq. 8 and using a Lagrange multipliers λ [i] q k q k to account for the constraints, the minimization can be reformulated as: with The solution to this minimization problem can be found by setting the derivative with respect to Bs i Q i as well as λ [i] q k q k to zero.Using some algebra we find the minimum for all sites i = N : Figure 6: Definition of the states used in the projection operator onto the link q k defined in the second line of Eq. 15 for link q 1 connecting sites i = 4 and j = 5.
This allows us to obtain a representation of the tangent space projector |Θ[B] = P T |ψ[T ] |Ξ as: Where <i,j>q k denotes a sum over all nearest neighbors i and j with the corresponding index q k connecting these two sites.The graphical representation of the states in the second line of the tangent space projector for the bond connecting sites i = 4 and j = 5 is shown in Fig. 6.Formally, this result resembles the projection operator obtained for MPS [29].The first term with positive sign corresponds to the forward time propagation of the site tensor.The second term on the other hand is the evolution backwards in time of the bonds between two site tensors and is a direct consequence of the gauge fixing of the tangent vectors used in Eq. 10.

Single-Site TDVP
With the representation of the projection operator in Eq. 15, we can go back to the projected time dependent Schrödinger equation (Eq.6) and integrate each term one by one using Trotter breakups [39].First, let us discuss a first-order update, which can later easily be modified to perform a second order integration.Since each term in the projection operator keeps all but one tensor fixed, the integration can be performed locally.Therefore, we define effective Hamiltonians for the sites i and for the links q k : and solve equations of the form: where A is either a site-tensor or a link tensor and k ) (positive sign).In matrix form, the exponential of these effective Hamiltonians can be efficiently calculated using Krylov exponentiation.
A full TDVP step is then given by a series of N − 1 local updates of a site tensor and the corresponding link tensor connecting the site to the end point as shown below.The single local update on site i and link q k is • Orthogonalize the TTN such that site i is the orthogonality center.
• Calculate the one-site effective Hamiltonian H eff = H (s i Q i ),(s i Q i ) (Eq. 16a) and forward time evolve (negative sign) according to Eq. 17 with A = C s i Q i .If site i is the chosen end point, stop here; otherwise continue.
• Reshape the time evolved tensor into a matrix . As usual, take the U -tensor as new tensor on site i.
• Calculate the link effective Hamiltonian k ) (Eq. 16b) for link q k .To do so, use the time evolved tensor obtained in the previous step for site i.Then evolve tensor from the previous step backwards in time (positive sign) according to Eq. 17.Finally, absorb the C-tensor onto the neighbor of site i along q k by multiplying it onto its site tensor.
A full TDVP time step can then achieved by the following sweeping procedure: 1. Choose a start and an end point; initialize site i as the chosen start point.

Perform the following steps until i is the chosen end point:
• Find the link q k ∈ Q i that connects site i to the end point.
• If any tensor attached to the other links Q i \ q k has not been updated, choose one of these links and choose one of the leaves attached to the corresponding segment of the TTN as new site i.
• Otherwise, perform a local update on site i as described above and choose the neighbor of site i along link q k as new site i.
3. Perform one last local update for the endpoint i = N as described above.
A depiction of the sweeping order for the TTN in Fig. 1 is shown in Fig. 7.The procedure described above defines a first-order time step.A second-order method can easily be obtained by performing the first order time step with ∆t 2 and then simply perform the exact same steps in reverse order corresponding to repeated second order Trotter breakups e τ (A+B+C) = e Importantly, this means that during the local update, the link update has to be performed before the site update (see also caption of Fig. 7).Note that for a given TTN, there can be several versions of this algorithm depending on the sequence of chosen indices in step 2. Very often though, the TTN structure itself defines some natural order when to time evolve which sites, as we will see in the next section for the FTPS tensor network.

Two-Site TDVP
It is also straightforward to generalize the single-site TDVP approach presented above to a two-site TDVP integration scheme which allows to dynamically adapt the necessary bond dimensions.To do so, we need to define the two-site effective Hamiltonian H 2-site for two sites i and j connected by the index q l : With this, only small modifications to the algorithm presented above are necessary.The single update for sites i and j sharing link q l becomes: • Orthogonalize the TTN such that site i is the orthogonality center.
• Calculate the two-site effective Hamiltonian H 2-site according to Eq. 18. and forward time evolve (negative sign) with • Reshape the time evolved tensor into a matrix A (s i Q i \q l ),(s j Q j \q l ) and perform an SVD . In this step one can also truncate the smallest Schmidt values.As usual, keep the U -tensor to update site i and C s j Q j as site tensor on site j, shifting the orthogonality center to site j.If site j is the chosen end point, stop here; otherwise continue.
• Calculate the one-site effective Hamiltonian H eff = H (s j Q j ),(s j Q j ) (Eq. 16a) for site j.
To do so, use the time-evolved tensor obtained in the previous step for site i.Then evolve tensor A = C s j Q j backwards in time (positive sign in Eq. 17).
A full two-site TDVP step can then be performed by: 1. Choose a start and end point.Initialize site i as the chosen start point.
2. Perform the following steps until i is the chosen end point: • Find the link q k and the corresponding neighbor j that connects site i to the end point.
• If any tensor attached to the other links Q i \ q k has not been updated, choose one of these links and choose one of the leaves attached to the corresponding segment of the TTN as new site i.
• Otherwise, perform a local update on site i and j as described above and go to site j, i.e., i → j.
A depiction of the necessary sweeping order of this two-site scheme is shown in Fig. 8 4 TDVP for FTPS An FTPS is a special TTN designed to efficiently encode states of multi-orbital Anderson Impurity Models (AIMs).An AIM consists of an interacting impurity coupled to a bath of free fermions with Hamiltonian c † mσk (c mσk ) creates (annihilates) an electron in chain m with spin σ on site k, where k = 0 denotes the impurity site (see Fig. 9 (b)).n mσk are the corresponding particle number operators.H int is the interaction Hamiltonian that only couples impurity degrees of freedom and for which we choose the Kanamori Interaction [18,49] without the spin-flip and pairhopping terms parametrized by two interaction strengths U and J where σ is the opposite spin direction of σ.In the following, we will use a combined index l = (mσ) to denote the orbital and spin-degrees of freedom.For a single orbital, an FTPS reduces to a MPS, while for multiple orbitals it has tensors with three link indices as depicted in Fig. 9 for a two-orbital model.It consists of a single MPS-like chain for the bath tensors of each orbital/spin and impurity tensors connecting the different chains.An FTPS for a N orb -orbital AIM has a total of N C = 2N orb chains.For simplicity we assume that each chain has the same number of bath sites N b .
According to the algorithm presented in the previous section, we first need to choose a start and end point.We choose to start at the outermost bath site of the first chain (site 4 in Fig. 9 (a)) and the outermost bath site of the last chain as end point (site 16 in Fig. 9 (a)).To actually perform the time evolution, we choose to employ a hybrid TDVP scheme using 2-site TDVP for the bath tensors as well as for the bath-impurity link, and 1-site TDVP for the impurity tensors itself and the corresponding impurity-impurity links.We choose to use 1-site TDVP for the impurity links, since 2-site TDVP becomes computationally expensive, since one would have to deal with tensors with four link indices.This leads to the following algorithm for a single time step: 1.For l = 1 : N C − 1 perform the following steps: -Perform a two-site step on sites i = (l, k) and j = (l, k − 1) (see Fig. 9 for the definition of the site-labeling).
• Perform a one-site step on the impurity tensor i = (l, 0); q k connects site (l + 1, 0) Impurity site Bath site l i n c r e a s e k increase Figure 9: Graphical representation of a FTPS tensor network for a two orbital model.For each orbital, we use two chains, one for each spin-species.(a) one way to label the sites is just to numerate them in ascending order.(b) a different way to label sites is to specify the chain (orbital m and spin σ) as well as an index (bath index k).This way to label sites resembles the labels used for the operators of the Hamiltonian in Eq. 19.
2. l = N C , for k = 0 : N b − 1 perform the following steps: • Perform a two-site step on sites i = (l, k) and j = (l, k + 1).
For the actual calculations, we apply the second order version of this algorithm by using only the half time step and reapplying each step in reverse order.Again, this also means that the order in the local updates changes.Note that the backwards propagation during the two-site involving an impurity site cancels with the subsequent forwards time evolution of the one-site step on the same impurity tensor.Therefore, these two steps can be omitted.As a first demonstration of this algorithm, let us compare the TDVP time evolution to the TEBD-like approach used in Refs.[18,19].Therefore, we look at the greater Greens function of the impurity defined by: |ψ 0 is the ground state of Hamiltonian H with ground state energy E 0 .For a degenerate two orbital model, Fig. 10 shows that the TDVP time evolution indeed produces the correct result.In a recent publication, the authors have shown, that for diagonal hybridizations, TDVP has larger errors than TEBD for the bath geometry chosen here [50].This means that for such systems, the TEBD approach is most likely preferable over TDVP.For more involved baths on the other hand, TEBD can become difficult to formulate as discussed next.
One of the major advantages of TDVP is that it allows to perform the time evolution for arbitrary couplings in the Hamiltonian between the sites, as long as an MPO with the same tensor network structure as the state can be found.Eq. 19 is in fact not the most general AIM, since the bath only couples diagonally to its impurity.Often, one is also interested in so-called off-diagonal hybridizations which can be encoded as hoppings from impurity l to a different bath l .Therefore, we can account for off-diagonal hybridizations by replacing the for a two orbital model between the TEBD time-evolution approach used in Ref. [18] and TDVP.The calculation was performed for a spin-and orbital degenerate model using 9 bath sites per orbital and spin with parameters k = −0.8,−0.6, • • • , 0.6, 0.8 and V [k] l = 0.1 ∀ k, l.Therefore, only diagonal entries of the Green's function are non-zero and for the diagonals there is only one independent function, i.e., G > l ,l (t) = δ l l G > (t).Interaction parameters were U = 1, J = 0.1.The impurity on-site energy was chosen to obtain particle hole symmetry, i.e., mσ0 = − 3U −5J

2
. The time step for TEBD was ∆t = 0.01 and for TDVP ∆t = 0.1, since TDVP generally allows to use larger time steps [50].Truncated weight (sum of all truncated Schmidt values) for TEBD was 10 −12 without restricting the bond dimension, and for TDVP it was 10 −9 for all links except the impurity-impurity links which were not truncated, but restricted to a maximal dimension of 50.hybridization terms in Eq. 19 with: It turns out that for each k, the matrix V ll can be chosen as a lower-triangular matrix.This means, that for a spin-symmetric, two-orbital model there are three free parameters for each value of k (instead of two for the diagonal hybridization).
As a second demonstration of the TDVP approach for FTPS we calculate the 2 × 2 matrix of the greater Green's function of such a spin symmetric two-orbital model.Since the TEBD approach we compared with in Fig. 10 is difficult to generalize to such off-diagonal hybridizations, we perform the calculation in the non-interacting case U = J = 0 and note that for tensor network based approaches this is a highly non-trivial situation.This is because the bipartitions defined by the links of the FTPS structure have non-trivial entanglement also for non-interacting systems, and the off-diagonal hoppings V

[k]
ll for l = l introduce entanglement between the orbitals, i.e., non-trivial links between the impurities.The results of such a comparison can be seen in Fig. 11.Having access to the exact solution, we also plot the difference between the exact and numerical Green's functions in the bottom panels.Again we find very good agreement between TDVP and the reference calculations.Finally let us demonstrate that the results indeed converge with respect to the control parameters.The left plot of Fig. 12 shows the scaling of the error as a function of ∆t and we indeed observe the expected ∼ ∆t 2 behavior at larger values of ∆t.The deviation of this behavior at smaller ∆t can be understood from the additional errors due to the truncations of the tensor network in the ground state as well as during time evolution.Additionally, we frequently observed that when TDVP is used in the star-geometry representation of the bath (with long-range couplings V mσ ), a good balance between truncation and ∆t is necessary.Surprisingly we found that it is often advantageous to use rather large time steps compared to what one would use in TEBD calculations.In the right plot of Fig. 12 we show the convergence of the error as a function of dimension of the impurity-impurity links which is usually the bottle-neck of FTPS calculations as these links need to transport the entanglement between the different orbitals l.Also here, we observe convergence with the control parameter, showing that TDVP indeed can be efficiently used to account for off-diagonal hybridizations.

Conclusion
We presented a generalization of the Time Dependent Variational Principle (TDVP) to general loop-free tensor networks (TTNs).The major advantage of TDVP over the commonly used TEBD approach is that the latter is often difficult to implement if long-range couplings are present in the Hamiltonian.TDVP on the other hand allows to perform the time evolution (either in imaginary-or real-time) for any Hamiltonian for which a representation in the same TTN structure can be found, which is often possible for long-range couplings.Using a similar derivation as in Ref. [29], we were able to find the projection operator onto the tangent space for any TTN -the central object in TDVP.Integrating the terms in the tangent space projector one after the other, equivalent to a Suzuki Trotter breakup, we were able to formulate TDVP in its single-site as well as two-site variant.We then applied TDVP to the FTPS tensor network which is a TTN especially suited for multi-orbital Anderson impurity In each segment, the top panel shows the Green's function itself, while the bottom panel shows the absolute value of the difference (t) = |G > exact (t) − G > TDVP (t)|, i.e., the numerical error.We obtained the exact solution from diagonalization of the hopping matrix at U = J = 0.The calculation was performed for a spindegenerate model using 9 bath sites per orbital and spin.We allowed off-diagonal hopping terms only between the orbital degrees of freedom, i.e., V [k] mm = 0.05 ∀ k for m = 2, m = 1.These off-diagonal terms correspond to a hopping processes from the impurity of orbital 2 to the bath of orbital 1.The TDVP time step was chosen ∆t = 0.1 and on-site energies were mσ0 = 0. Note that the off-diagonal hybridizations break the orbital degeneracy, albeit for the parameters chosen only slightly and the differences between the two orbitals are barely visible.Truncated weight (sum of all truncated Schmidt values) during DMRG and the time evolution was 10 −9 , except during the time evolution of the impurity impurity links where no truncation was performed.During DMRG as well as time evolution, the impurity-impurity links were restricted to 140.

Figure 1 :
Figure 1: Example of a TTN with 7 tensors with different numbers of link-indices on each site.Each dot represents a tensor and each line an index, where a connected line implies summation over this index.The open lines are the physical indices s 1 • • • s N (N = 7), while the connected lines are the link indices q 1 • • • q L (L = 6).Cutting the link between sites 4 and 5, as indicated by the dashed line, results in two disconnected tensor network segments and defines a notion of left and right at each link.In this example, sites 1, 2, 3, 6 and 7 are the leaves of the TTN.

Figure 2 :
Figure 2: Gauge degree of freedom in tensor networks.At each link, one can insert an identity 1 = G • G −1 without changing the physical state |ψ .By absorbing G into one tensor and G −1 into the other, we obtain a different representation of the same state |ψ .

Figure 4 :
Figure 4: If the orthogonality center of the TTN depicted in Fig. 1 is placed on site 4, the center tensor C has four link indices q 1 • • • q 4 .Each of these links corresponds to one of four mutually orthogonal set of states |q Similar to MPS, we can now reformulate the projection problem of an arbitrary state |Ξ onto the tangent space |Θ[B] = P T |ψ[T ] |Ξ as a minimization problem: min B || |Θ[B] − |Ξ || 2 ,

τ 2 C e τ 2 B e τ A e τ 2 B e τ 2 C
used on Eq. 15.

Figure 7 :
Figure7: Update sequence to perform a first-order single-site TDVP time step from time t to t + ∆t for the TTN shown in Fig.1.Start point is site 1 and end point site 7. Yellow denotes tensor that are updated in the current step.Red and blue tensors indicate whether this tensor is taken at time t + ∆t (red) or t (blue).Triangles indicate the orthogonalization of each tensor.Updates on site-tensors are in forward direction (negative sign), while updates on bond-tensors are backwards time evolutions (positive sign in Eq 17) For a second order update, first perform all steps (1) → (13) with time step dt 2 in the order shown and then reapply them in the reverse order (13) → (1), again with time step dt 2 .

Figure 10 :
Figure10: Comparison of the impurity greater Greens function G > (t) for a two orbital model between the TEBD time-evolution approach used in Ref.[18] and TDVP.The calculation was performed for a spin-and orbital degenerate model using 9 bath sites per orbital and spin with parameters k = −0.8,−0.6, • • • , 0.6, 0.8 and V

Figure 11 :
Figure 11: Comparison of impurity greater Greens function G > ll (t) for a two-orbital model between TDVP and the exact solution.In each segment, the top panel shows the Green's function itself, while the bottom panel shows the absolute value of the difference (t) = |G > exact (t) − G > TDVP (t)|, i.e., the numerical error.We obtained the exact solution from diagonalization of the hopping matrix at U = J = 0.The calculation was performed for a spindegenerate model using 9 bath sites per orbital and spin.We allowed off-diagonal hopping terms only between the orbital degrees of freedom, i.e., V (m σ ) = δ σσ V [k]mm and therefore G (mσ)(m σ ) = δ σσ G mm .The parameters were k = −0.8,−0.6, • • • , 0.6, 0.8 for all orbitals, diagonal hybridizations V [k] mm = 0.1 ∀ k, l and off-diagonal hybridizations V 1. Start point is site 1 and end point site 7. Yellow denotes tensor that are updated in the current step.Red and blue tensors indicate whether this tensor is taken at time t + ∆t (red) or t (blue).Triangles indicate the orthogonalization of each tensor.Updates on site-tensors are in forward direction (negative sign), while updates on bond-tensors are backwards time evolutions (positive sign in Eq 17) For a second order update, first perform all steps (1) → (13) with time step dt 2 in the order shown and then reapply them in the reverse order (13) → (1), again with time step dt 2 .Update sequence to perform a first-order two-site TDVP time step from time t to t + ∆t for the TTN shown in Fig.1.Start point is site 1 and end point site 7. Yellow denotes tensors that are updated in the current step.Red and blue tensors indicate whether this tensor is taken at time t + ∆t (red) or t (blue).Triangles indicate the orthogonalization of each tensor.Updates on two-site tensors are in forward direction (negative sign), while updates on a single site are backwards time evolutions (positive sign in Eq. 17) For a secondorder update, first perform all steps (1) → (11) with time step dt 2 in the order shown and then reapply them in the reverse order (11) → (1), again with time step dt 2 .