Hierarchical generalization of dual unitarity

Quantum dynamics with local interactions in lattice models display rich physics, but is notoriously hard to study. Dual-unitary circuits allow for exact answers to interesting physical questions in clean or disordered one- and higher-dimensional quantum systems. However, this family of models shows some non-universal features, like vanishing correlations inside the light-cone and instantaneous thermalization of local observables. In this work we propose a generalization of dual-unitary circuits where the exactly calculable spatial-temporal correlation functions display richer behavior, and have non-trivial thermalization of local observables. This is achieved by generalizing the single-gate condition to a hierarchy of multi-gate conditions, where the first level recovers dual-unitary models, and the second level exhibits these new interesting features. We also extend the discussion and provide exact solutions to correlators with few-site observables and discuss higher-orders, including the ones after a quantum quench. In addition, we provide exhaustive parametrizations for qubit cases, and propose a new family of models for local dimensions larger than two, which also provides a new family of dual-unitary models.


Introduction
One of the pivotal problems in quantum manybody physics is understanding the dynamics of extended systems with local interactions.Although the local interactions are simple, they generate complex dynamics, which is in general too hard to describe.The dynamics can be characterized with different probes such as local correlation functions, entanglement spreading, and other quantum information quantities.The un-derstanding of this type of dynamics is currently at the center of attention in many fields spanning from nonequilibrium statistical mechanics, quantum information, condensed matter to highenergy physics and quantum gravity.
However, the complexity of quantum dynamics, both analytically and numerically, presents a significant hurdle.For example, the bond dimension of Matrix Product States (MPS) typically increases exponentially due to the linear growth in entanglement entropy [1,2].This necessitates the use of solvable models to unravel many-body behavior.The most well-known examples are noninteracting (Gaussian) systems such as free fermions or bosons, Clifford circuits, and Betheansatz interacting integrable models [3].Unfortunately, all of these models are not chaotic, in contrast to generic examples.If one is prepared to average over an ensemble of systems, random unitary circuits [4] provide examples of solvable chaotic dynamics.However, the averaged results miss a lot of relevant physics and are less relevant for translationally invariant systems.
A recently discovered solvable family of models, known as dual-unitary quantum circuits [5], has distinctly different solvable structure which does not require averaging, and moreover contains both integrable and chaotic examples.The basic property, which enables solvability is the unitarity of the local gates in the space direction.In this paper, we generalize this to a condition on a few gates and unravel new families of solvable models.
Despite the success of the dual-unitary circuits (and its extension) in describing physical properties of non-equibrium dynamics, they exhibit some non-universal features.Firstly, in dual-unitary circuits the non-vanishing correlation functions exist only on light-cones edges.Secondly, the thermalization of local observables is instantaneous.The fundamental reason behind these non-universal features is that their single gate's conditions are too restrictive.In fact, this family of models doesn't even include many gates that are solvable by other methods, e.g. the Identity, Controlled-Not and Controlled-Z.This raises the question: Can the dual-unitary condition be relaxed to allow for richer physics while still maintaining the solvability of the spatial-temporary correlation function?
In this paper, we answer these questions in the affirmative by relaxing the dual-unitary condition and extending it to a hierarchy of conditions that contains more and more local gates.The dualunitary condition on one gate forms the the first level of the Hierarchy denoted by L 1 , whereas the circuit with two-gate condition at the second level of the Hierarchy L 2 allows for the exact calculation of two-point correlation functions, which are richer than for dual unitaries, i.e. non-vanishing at the same site and different times.Moreover, when quenched from the solvable initial states, the L 2 circuit exhibits non-trivial thermalization of local observables.We go beyond L 2 and show that higher levels of hierarchical circuits limit the maximum speed of information spreading.We provide complete parameterization for the L 2 and L 3 circuits in the qubit case.For the larger local Hilbert space dimensions, we propose a new method to construct a class of circuits using the Clifford group that is analytically trackable.This method can also be used to construct new families of dual-unitary circuit.
The paper is structured as follows: In Sec. 2, we introduce the notation.Subsec.2.1 reviews the dual-unitary circuits and introduces our parametrization method using Clifford groups.The hierarchical generalization is outlined in subsec.2.3.After that, we dive into details of the L 2 and L 3 circuits in subsecs.2.4 and 2.5, including their parametrization.In Sec. 3, we discuss the physical applications of the different levels of Hierarchical circuits.Subsec.3.1 considers the two-point correlation functions for L 1 , L 2 and L 3 circuits.In subsec.3.2, we extend our discussion to the correlation functions of multisite observables, and three-point correlation functions.Subsec.3.3 discusses the evolution of an initial state from a quantum quench.We generalize the solvable initial states [6] for the L 2 circuit and explore the relationship between quench dynamics and quantum thermalization.In Sec. 4, we summarize the main results of the paper and discuss future directions.

Hierarchical generalization of dual unitarity
In this paper we consider a chain comprised of L cells, with each cell containing 2 sites at integers and odd-half integer sites.At each site there is a Hilbert space with a local dimension D. Consequently, the corresponding total Hilbert space is H = (C D ) 2L .The local basis is denoted by |j⟩ with j = 0, 1, • • • , D − 1.The chain's dynamics is governed by a brickwall Floquet circuit where T 2L is a periodic translation operator on 2L sites, and u a local gate.Here, for simplicity, we assume translational invariance of the circuit and introduce periodic boundary conditions.However, our result can be easily generalized to non-uniform cases and open boundary conditions.Above, we graphically represented local unitary gates with dimension D 2 × D 2 by a box with in-coming and outgoing legs, satisfying unitarity conditions ( Our results can be more succinctly expressed in the folded picture, where an operator over (C D ) 2L is vectorized to a vector in (C D ) 4L by the linear map on the basis The time evolution in Schrodinger picture can also be vectorized to Graphically, u † is folded back behind u, thereby forming a joint operator w ≡ u ⊗ u * .It is also convenient to denote the vectorized identity operator in the folded picture as an empty bullet With these notations, the unitarity condition (3) is graphically expressed as = and = .

Dual Unitarity
Understanding the dynamics of extended locally interacting systems is at the core of quantum many-body physics.However, this problem is usually analytically intractable and numerically exponentially hard.To make progress, we need some additional structure.One possibility is to demand the so-called dual-unitarity condition [5] mentioned in the introduction, which enables various exact calculations even for chaotic dynamics.Dual unitarity demands that the gate u is unitary even if we exchange the roles of space and time.This switching corresponds to changing A gate is dual-unitary [5] if both u and ũ are unitary, so in addition to (3) we also require which in the folded graphical language yields The family of models defined in this way encompasses free, interacting integrable and chaotic models [5].The parametrization of dual-unitary gates for D = 2 has been fully determined [5].Despite a lack of the complete parametrization of dual-unitary gates for D ≥ 3, several families have been proposed [29,30,31,32,33,34].We proceed to add another family to the list, resulting in a novel extensive family of dual-unitary gates in higher dimensions.

New Parametrization of Qudit Gates
In the following we provide a non-exhuasive parameterization of qudit gates, which will prove useful for constructing examples at D > 2, both for dual-unitary gates and their generalizations.This parametrization first appeared in Ref. [45] in the study of Operator-Schmidt decomposition, and is a special case of a novel, more general framework we report in App.A; however, since we do not need this more general framework in the following analysis, we restrict ourselves to the special case in the main text.
Consider the following family of two-qudit unitary gates: where v 1 , v 2 , v 3 , v 4 are single site unitary gates, and u 0 is defined as Here {θ p,q } 0≤p,q≤D−1 is a collection of U(1) phases, and {|ψ p,q ⟩} 0≤p,q≤D−1 is an orthonormal basis for the 2-qudit Hilbert space (C D )2 defined as where σ, τ are the D × D dimensional generators of the Clifford group satisfying the relations with ω = e 2πi/D a D-th root of unity.The matrices σ, τ generate the full matrix algebra M D (C), and we can always choose σ to be diagonal and τ to be real.Explicitly, they are defined as where |D⟩ ≡ |0⟩.
(16) Notice that the single site unitary gates v 1 , v 2 , v 3 , v 4 do not appear in the above expression.We simplify Eq. ( 16) further with the following relations satisfied by τ p,q τ p,q τ r,s = ω qr τ p+r,q+s , They follow from Eqs. (13) and (14) by straightforward computation.Simplifying Eq. (16) using Eq.(17), and comparing the coefficients of both sides using the fact that {τ p,q } 0≤p,q≤D−1 forms a basis of the matrix algebra M D (C), we obtain 0≤p,q≤D−1 θ * p,q θ p+k,q+l = 0, for (k, l) ̸ = (0, 0).(18) In this way, the original dual unitarity condition, which involves 2D 4 equations and 2D 4 − 1 real unknowns simplifies to a set of D 2 − 1 equations with D 2 − 1 real unknowns (notice that we can set θ 0,0 = 1 without loss of generality).A simple yet nontrivial ansatz for θ p,q is 1 where µ ∈ Z, and λ, ν ∈ Z if D is odd while λ, ν ∈ Z/2 if D is even (which guarantees that θ p,q is periodic both in p and q with period D).This ansatz also results in the perfect tensors in odd dimensions which are found in [32].Inserting the ansatz (19) into Eq.(18), we see that dual-unitarity requires that k = l = 0 is the only solution to the following system of equations For example, when D = 3, λ = µ = 1, ν = −1 satisfies this condition.In later sections we will use the ansatz Eq. (11) and Eq. ( 19) to find examples of hierarchical generalizations of dual-unitary gates.
In this subsection we recapped the basics of dual-unitarity and introduced a novel family of dual-unitary models for D > 2. A particular subset of solutions from this family appeared before in [34].This leaves us in a good position to introduce the generalization in the next subsection.

Hierarchical Generalization
As mentioned in the Introduction, dual unitarity imposes conditions on only a single gate, which restricts the possible physical behaviours.It also excludes certain fundamental and well-known gates, such as the Identity and the Controlled-Not gates, which are solvable yet not dual-unitary.To unveil more intricate quantum dynamics and include these Clifford gates in a more general notion of solvability, we define a hierarchy of conditions.This gives us new families of models.
Since only one green box plays a part in the dual-unitary condition (9), we will call dualunitarity also the first level of the Hierarchy and denote it as L 1 .In the subsequent subsections 2.4 and 2.5, we extend the concept of dual unitarity (L 1 ) to conditions involving two and three gates, resulting in the second level L 2 and the third level L 3 of the Hierarchy.

Second level of the Hierarchy
In this subsection, we introduce the gates from L 2 , which are more general than dual-unitary gates.L 2 contains CNOT and identity, as well as a large family of non-trivial gates whose dynamics cannot be solved by any previous techniques and reveals richer physics.The gates from this family fulfill a condition involving two gates, which is weaker than a dual-unitary condition.Here we focus on the case when this condition holds both from the left and the right, but we comment on the non-symmetric case in Sec.3.1.3.The two conditions are: (21) Algebraically, we can express the condition as This, together with the fact that the identity is in In the following we focus on the gates that are in L 2 but not in L 1 , the set we denote as L 2 = L 2 − L 1 .Note that a necessary condition for a gate u to be in L 2 is that ũ is not invertible, since otherwise one can contract both sides of Eq. ( 21) with ũ−1 and recover the dual unitary condition.
Similarly to the L 1 case [5], we can figure out the complete parametrization of L 2 for qubits.As explained in App.C.1, we used analytical analysis with the numerical help of Mathematica.When D = 2, an exhaustive parametrization of 2-qubit gates is Here σ j are Pauli matrices and v 1 , v 2 , v 3 , v 4 are all single site gates from SU(2), and 0 ≤ J j < π/2.
We may simplify the gate structure by setting v 3 = v 4 = I D without any loss of generality3 .The trivial example from L 2 is a tensor product of two single-site operators.Apart from that, the L 2 condition fixes J z = π 4 , J x = J y = 04 and v 1 , v 2 to be elements of the set where U (r, θ, ϕ) is defined as e ir(cos θ σz+sin θ cos ϕ σx+sin θ sin ϕ σy) , representing a SU(2) on the Bloch sphere.Geometrically, this specific combination of r, θ, ϕ represents a rotation that maps σ z to the x − y plane.
The dimension of L 2 can be counted as follows.Out of 12 parameters determining the 4 local SU(2) gates, e.g. the Euler angles, two are redundant because the rotation around the z−axis commutes with the Ising interaction resulting from J z = π 4 , J x = J y = 0. Further, Eq. ( 24) provides 2 constraints.After considering the global phase, the total independent parameters to characterize a qubit L 2 circuit is 12 − 2 − 2 + 1.Therefore, we have defined a new 9-dimensional family of solvable models which are not part of 12-dimensional set of L 1 gates [33].
Note that the control not gate (CNOT) can be decomposed into the form of Eq. ( 23) (with v 4 different than identity) as: with J x = 0, J y = 0, J z = π 4 and an addition Hadamard gate.To check that this satisfy (24), we include v 4 back in the previous layer of the gates therefore obtaining combined(c) In the case of bigger local dimensions D, we do not yet possess a complete parametrization for the L 2 case.Nevertheless, we can discern two distinct and rich families.The first family is associated with generalized Controlled-NOT gate in higher dimension surrounded by 4 single site with C τ = i |i⟩ ⟨i| ⊗ τ i .Following a similar argument as below Eq. ( 23), we set v 3 = v 4 = I D .In this case, v 1 and v 2 must satisfy These two equations share a symmetry of exchanging columns and rows between themselves.
The second family is derived using the Clifford group method from subsection 2.1.Utilizing the proposed ansatz from Eqs. ( 10) and ( 11), we set v 3 = v 4 = I D and simplify Eq. ( 22) to: Here b is a shorthand for 0≤p b ,q b ≤D−1 .The above equation should hold for ∀(s, t) ̸ = (0, 0) and (k, l) ̸ = (0, 0).If all terms in the first sum vanish separately, we obtain the L 1 .From these nonlinear equations, we can derive a family of L 2 , which is just one of the many possible solutions.The family is defined for D = 4k + 2 as Another nontrivial example is given by In both of the two examples, the so far unspecified v 1 and v 2 belong to a non-trivial subset of SU(D).Finding all possible v 1 and v 2 is, in general, hard.One can try guessing good candidates and check if the conditions in Eq. ( 28) are satisfied.There is a way of simplifying the conditions using the structure of the Clifford ansatz, which help in deducing v 1 and v 2 , see App.C.2 for the details.Different choices lead to both ergodic and non-ergodic dynamics.
Before concluding this subsection, we would like to highlight that the two equalities in Eq. ( 21) and unitarity additionally imply . i.e., a Hierarchical condition along with unitarity implies its 180 • rotation version.The proof is shown in App.D. Eq. (31) will play an important role in computing the spatiotemporary correlation functions.

Third level of the Hierarchy
Following the principles from subsection 2.4, we define the third level hierarchical condition for L 3 as An immediate observation reveals that a gate characterized as L 2 is also the L 3 .Nonetheless, we are again interested in the special subset of L 3 which does not belong to L 2 , designated as We again use the complete parameterization of 2-qubit gates (23) and wlog set v 3 = v 4 = I 2 .The condition defining L 3 is satisfied either for all diagonal gates or for the case where J x = J y = 0 and any J z with v i satisfying v i = cos ϕ i σ x + sin ϕ i σ y , i ∈ {1, 2}.
To get some examples for D > 2, we use our Clifford gate parametrization method from Secion 2.1.The algebraic equation of θ p,q can be found in the App.B. To obtain some examples, we take the single-site operators v i to be the identity.Some classes of the solutions obtained in this way are shown below.
In principle, nothing stops us from going beyond the L 3 , by demanding even more general condition with even more gates.We expect the examples to be constructed in a similar way.

Spatio-temporary correlator functions
In this subsection we focus on the spatio-temporal correlation functions, which are the most common objects to characterize the dynamics.In particular, they provide information about the thermalization and ergodicity of the system.
In most cases, the exact non-perturbative calculation of the spatio-temporary correlators is only available in free models and to some extent in interacting integrable ones [46,47,48].Important progress has been made in understanding the correlations also in chaotic models, in particular dual-unitary (L 1 ) circuits which we extend here.
Due to the trivial propagation of an identity operator, we are only interested in the correlation function between two traceless Hilbert-Schmidt normalized operators a i , b j 5 .Working in the Heisenberg picture, the spatio-temporal correlation function of normalized local operators can be expressed as Employing the time unitarity enables us to simplify the circuit from the bottom and top, yielding

Dual unitarity
For completeness, here we briefly summarize the result for dual-unitary circuit from [5].Intuitively, the time unitarity and space unitarity both determine a light cone outside which the correlation function vanishes.Therefore, correlators can solely manifest at the intersection of these two cones, forming a 1-dimensional straight line that precisely bisects the temporal and spatial directions with other correlators vanishing.That other correlations vanish can be seen by repeatedly applying (9) to expression (36), which results in ⟨ | ⟩ = 0 since the operators a and b are traceless.In Eq. (37) each time step is just a quantum channel over D × D Hilbert space.Therefore, the correlators for the L 1 circuits can always be calculated efficiently and propagates only along two directions with maximal speed.

L 2 circuits
Moving beyond the L 1 , we are interested in which new features appear in the L 2 circuits.Said differently, we are interested in what happens if the weaker condition (21) defining L 2 is satisfied but dual unitarity condition (9) is not.
We apply Eq. ( 21) to Eq. ( 36) and further simplify the circuit to Lastly, Eq. ( 31) is utilised to address the corner of the path, leading to This correlator vanishes because the discontinuous path will be simplified to Tra i Trb j and both are traceless according to our assumption.
Therefore, the existence of nonvanishing correlators is limited to three possible directions, either at the light cone or at velocity zero.
These expressions can be written using four single qudit channels: To simplify the analysis, we assume j to be an integer, and the other case follows analogously.Thus ) This C ij (t) behaves differently than that in the dual-unitary case [5], as the circuits from L 2 allow for an additional non-vanishing direction along the time axis.
Let us mention here the connection with triunitaries circuits proposed in [39], where the correlation function also exclusively manifest in the same three directions.In fact, we can group the legs of two 2-qubit gates into a 3-qubit gates as However, in the tri-unitary case, the condition is which is a much stronger condition than Eq. ( 21).Interestingly, in the context of qubits (D = 2), it is impossible to observe all in principle allowed physical behaviors.The correlations along the light rays vanish since channels ϵ L and ϵ R correspond to the total depolarizing channel.Nevertheless, when D > 2, there are examples manifesting all of the properties discussed above, i.e. nonvanishing correlations in all three directions at all times.In other words, both the correlations in Eq. ( 40) are nontrivial.A such example is given in Eq. ( 29) which is also shown in Fig. 2(a).In this figure, the operator has support on two nearest neighbor sites, to eliminate the odd/even effect (for details see subsec.3.2).

L 3 and higher levels
In the case where the gate is classified as L 3 , the correlation function is reduced to Intriguingly, this correlator does not vanish within the entire light cone, and no closed expression for it can be derived with a scaling polynomial in system size.Nevertheless, the hierarchical conditions still imply that some of the correlations are zero.As a general rule, for a kth level of Hierarchical circuit L k = L k − L k−1 , the correlations are nonzero along both the middle part as shown in Eq. ( 45) and the maximal velocity light rays.If we focus on the middle part and consider it as the inner light cone, its inner light cone velocity will be suppressed to for k ≥ 2, see Fig. 1.Notice that a Hierarchical condition along with unitarity always imply its 180 • rotation, as explained in Eq. ( 21), Eq.Finally, we also illustrate the possibility of choosing different types of condition (or no condition at all) from the left and right directions, resulting in an asymmetric inner light cone.This possibility is shown and discussed in Fig. 1.

Bigger operators and higher orders
In contrast to previous research, which concentrated primarily on correlators supported on a single site, exploring operators with multi-site support sheds light on more intricate underlying physical phenomena.Specifically, for correlators supported on multiple sites, their behavior resembles that described in Eq. (42), where correlation functions manifest exclusively along three directions.Remarkably, in the case of qubits, these correlation functions persist over time unlike the single-site supported ones.
Here we present examples of the correlators on nearest neighbor sites.The derivation is essentially the same as in the previous section with some special attention to the simplifications around the operators.The location of a multisite operator is indexed by its leftmost site.For simplicity, we assume j is an even number and represent the nearest neighbor two-site operator by a black square = .
The correlation function in Eq. ( 46) along the light cone is exactly the same as Eq.(42).By implementing the quantum channel defined as we can express Eq.(47) in a compact analytical Let us have a look at the correlation function's intriguing temporal decay.Generally speaking, the long-term behavior of the correlation function will be dominated by the largest eigenvalue λ of the quantum channel, evolving as ∼ λ t6 .If |λ| = 1, the correlation function will persist without decay.This behavior is referred to as nonergodic.A simple example of L 2 circuit in higher dimension, Eq. ( 29) with D = 6, v 1 = v 2 = I D , falls into this class, as illustrated in Fig. 2(a).Conversely, if |λ| < 1, the correlation function  29).(b) The correlation function for the qubit case supported on 3 sites.The gate is an ergodic element of L 2 , with parameters given in the App.E. The asymmetry between the left and right sides results from even/odd effects.In both figures, a i = b j = h, with h a random normalized traceless Hermitian operator.j is fixed at 20.The location of a multi-site operator is defined as the location of its left end.
will be ergodic and exhibit an exponential decay, a point we will come back later in Sec.3.3 in the context of quantum quenches.An ergodic example can also be constructed from Eq. ( 29) by choosing D = 6, v i = σ x ⊗ κ i , i ∈ {1, 2} for almost any κ i ∈ SU(3).
The scenario with three sites operators follows the preceding discussions without any additional difficulty.The correlators along the time axis can be expressed with the quantum channel Its largest non-trivial eigenvalue is typically smaller than 1.Fig. 2(b) showcases this ergodic qubit circuit.We further examine the correlators at i = j and i = j + t in Fig. 3.
Moving beyond the scope of 2−point correlation functions, 3−point correlation functions provide more information of the non-equilibrium dy-  namics.They are defined as follows: If i and j are on the same side of k, the 3−point correlation functions become trivial for L 2 circuits, i.e., either vanishes or reduces to 2−point correlations.Therefore, without loss of generality, we can assume i < k < j such that (51) Unlike the L 1 , there are non-trivial correlation even when both a and b are strictly inside the light cone.
Nevertheless, despite the fact that the L 2 condition greatly simplifies the circuit complexity, it does not fully resolve all computational challenges.Viewing from Eq. ( 51) we obtain that the maximum number of qudits (legs) we need to store when contracting the graph diagonally is l+ 1  2 with l labeled in Eq. ( 51), thus the computational complexity scales as e O(|(t 1 −k+i)−(t 2 −j+k)|) .

Quantum Quench in L 2 circuits
In this subsection, we examine the correlation functions for L 2 circuits following a quantum quench, i.e. originating from an initial density matrix ρ L (0).Here ρ L (0) can either be a pure state or a mixed state with a local purification.Therefore we can write it with a local purifying Matrix Product State (MPS) as [43,6], where Here γ i is the purification index, which we sum over in ρ L (0).Without additional specifications, the gates in this subsection are assumed to be from L 2 .Graphically, the vectorized density matrix can be represented as (53) A physical density matrix is normalized in the thermodynamic limit where E(0) is the space transfer matrix This implies that E(0) has a unique nondegenerate left and right fixed point whose eigenvalue is one with ⟨△ |□⟩ = 1.

1−point correlation functions
The most basic information about the dynamic from a quantum quench is contained in the 1−point correlation function, which specifies the relaxation and thermalization of local observables.It is defined as Here, E(t) = E 1 (t) and E O 1 (t) are appropriate space transfer matrices In this subsection we always assume periodic boundary conditions with large enough L such that the transfer matrix can be replaced by its fixed point.
If the initial state satisfies the following with being the right eigenvector of E(0), 7 we can find an eigenstate of the transfer matrix with eigenvalue 1 Since lim L→∞ Tr (ρ L (t)) = 1 due to the normalization, the largest eigenvalue of the transfer matrix must be 1 and non-degenerate.Therefore, the 1−point correlation function can be analytically calculated with this eigenvector.
It is worth to note that if the initial state satisfies the condition = , Eq. ( 59) is automatically satisfied.This implies that we have identified a larger solvable class than both the pure solvable initial states for dual-unitary evolution [6] and in general mixed initial states for open 3-way unital evolution [43].
The correlation function for 2−site observables after a quench is very similar to Eq. (47).The only difference is the substitution of at the 7 Note that by contracting the left top leg with an empty bullet in both sides of Eq. ( 59), we can directly show that Eq. ( 59) implies is the right eigenvector of E(0).

base for
.
where ⟨△| and |□⟩ are the left and right fixed points from Eq. ( 56).Following the same argument as in Sec.3.2, the long time behavior is dictated by the largest eigenvalue of the quantum channel, Q.The eigenspectrum of Q can be completely deduced for qubits.Using the general parametrization of gates from L 2 in accordance with Eqs. ( 23) and ( 24), whereby , only two non-zero eigenvalues, {1, λ}, exist.Eigenvalue 1 corresponds to the trivial eigenvector, i.e. identity and λ is given by (62) Here r 1 , r 2 are real numbers from the interval ].The dynamics is non-ergodic with λ = 1 if ϕ 1 = ϕ 2 + 0, π and either r 1 = π − r 2 or cos 2r 1 = cos 2r 2 = 0. On the other hand, the circuit is nonergodic also with λ = −1 if cos(ϕ 1 − ϕ 2 ) = 0 and the last term in Eq. ( 62) equals 1, for example, by r 1 = π 2 , r 2 = π 4 .This gives all possible nonergodic circuits for two-qubit gates from L 2 , apart from a tensor product of single-site gates.All other examples are ergodic and show exponential decay.
The most straightforward solution of initial states from Eq. ( 59) is the pure state with a bond dimension 1.One example is the qubit Bell state which we use in Fig. 4, where we show both thermalizing and non-thermalizing behaviors.In contrast to the situation here, for L 1 (dual-unitarity) the expectation values of local observables trivially vanish at all times, even for clearly nonergodic circuits such as circuits made out of SWAPs.63).The dynamics is governed by a L 2 circuit (symbols).We show its linear fit by solid lines.The blue and green data points show an exponential decay of the expectation value of observable, indicative of thermalization dynamics.In contrast, the red points remains constant at long times implying nonthermalization and non-ergodicity.The parameters for these three circuits can be found in App.E.

2−point correlation functions
Finally, we are moving to the 2−point correlation functions after a quench, which are non trivial even for the L 1 [6].Diagrammatically they are represented as By leveraging time unitarity, we can deduce the emergence of two backward propagating light cones originating from the operators, interconnected by the initial state's transfer matrix E(0) The number of repetitions of the transfer matrices E(0) in the central region depends on the separation, given by j − i − (t 1 + t 2 ) + 1 2 .The suitable L 2 2-point correlator solvability condition for this correlation function is and similar expression from left Here |△⟩ and |□⟩ are the left and right unique fixed point from Eq. ( 56).Note that this condition is stronger than the solvability condition for the 1-point correlators.Following the same argument in [43], it can be shown that any MPDO in the local purification form fulfilling solvability condition can be cast in right canonical form, which we use.Thus is replaced by (vectorized identity).
Supposing that the length L is large enough, we can replace the outer parts of Eq. (65) by Eq. (56).After this simplification, we obtain The leftmost and rightmost corners can be further reduced by the first equations from (66) and (67).Then with the help of the second equations from (66), (67) and the contraction properties of the L 2 , we arrive at Applying the first equation from (66), we finally obtain a considerably simplified expression Interestingly, if the separation j − i is shorter than t 1 + t 2 , the correlation function will vanish for traceless operators, which is the same as for L 1 [6].
Solvability condition from Eq. (66) can also be formulated in an algebraic form.If we define an operator K γ as the solvability condition for an initial state in the right canonical form can be written as The same discussion holds from the left side, except that we need to keep |△⟩ instead of | ⟩.
For example, where the initial state is pure, that is, without summation over γ, the first equation implies that K is proportional to a unitary operator.This leads to the conclusion that U is a L 1 , as per the second equation.Consequently, a solvable initial state for L 2 can only be a mixed state.As an example, if we choose U = CNOT and the bond dimension for MPS to be 1, a nontrivial solution is

Conclusions and perspectives
In this paper, we have generalized dual-unitary circuits to larger solvable classes referred to as levels of the Hierarchy.These families of circuits with comparatively relaxed conditions exhibit richer physical properties than the dualunitary ones.Most intriguingly, their correlation functions are non-zero at the initial site and later times.Furthermore, these circuits show nontrivial thermalization of local operators, in contrast with the standard dual-unitary circuit.
In our work we provide the complete parametrization of L 1 , L 2 , L 3 for qubit circuits.For local dimensions bigger than two, we propose a systematic approach using the Clifford group to construct novel examples, which include both ergodic and nonergodic instances.
It is worth mentioning that even though we focus on the translational invariant dynamics, the results and methods presented do not depend crucially on this property and can be extended to time and space inhomogeneous settings.This study, therefore, provides a strong foundation for further explorations in the realm of solvable models and their real-world applications.They can be directly implemented on current NISQ devices, and used as benchmarks in the realms where ordinary classical simulations are no longer possible.
We mentioned that one interesting example is the CNOT gate and its extensions.They are a particular cases of Floquet East model, which shows interesting physics, e.g.hydrophobicity and localization, both in classical [49] and quantum realm [50].Interestingly, for a subset of the solvable L 2 connected to CNOT, exact solutions can also be obtained using the so called zipper equations [51] by devising an exact eigenvectors of the transfer matrix [52].
Of course, the ideas presented here can be extended to conditions on bunch of gates in different configurations from the ones presented here.One direct way is generalizing the dual-unitarity round-a-face [33].It is feasible to define the L 2 condition in this context as We anticipate these circuits to exhibit similar features as the L 2 circuits discussed here.For examples, their correlation functions might exclusively exist along three directions.However, parametrization of these types of circuits and their subsequent physical applications represent an intriguing path for future research.One can also be more bold, and demand for instance the following conditions: It can be shown, using a simplification procedure similar to that in Sec.3.1, that this condition leads to the exact solvability of two point spatial-temporal correlation functions in the region −t/3 < x < 0. This is interesting since it gives us a 2D space-time region where the spatialtemporal correlation functions are solvable and non-vanishing, while all dual-unitary circuits and their L 2 generalization have vanishing spatialtemporal correlations except on a few lines.Unfortunately, to date, we have not yet been able to find a nontrivial unitary gate satisfying this condition, due to the computational difficulty of solving this equation.We leave this possibility as an open direction for the future.

A A more general parametrization of dual-unitarity and hierarchy gates based on finite groups
In this appendix we give a further generalization of the parametrization of dual-unitarity and hierarchy gates presented in Sec.2.1, based on projective representations of finite groups.Let G be a finite group and let Γ be an irreducible, unitary finite dimensional (with dimension D) projective representation of G, i.e.Γ † a = Γ a −1 , ∀a ∈ G and where We use Eq.(10) again for the parametrization of 2-qudit unitary gate u, but now u 0 is defined as where {θ a } a∈G is a collection of U(1) phases, and {|ψ a ⟩} a∈G is an orthonormal basis for the 2-qudit Hilbert space (C D ) 2 defined as We proceed to investigate conditions on the parameters {θ a } a∈G for u to be a dual-unitary gate.After a space-time reshuffling of indices defined in Eq. (7), , where Then the unitarity condition (8) on ũ is equivalent to a,b∈G We simplify Eq. (81) further with Eq. (75) and, using the fact that {Γ a } a∈G forms a basis of the matrix algebra M D (C), we obtain where ω is a k-th root of unity.The projective representation Γ is defined as The parametrization presented in Sec.2.1 corresponds to the special case n = 1.

B The parametrization of the Hierarchical gates in higher dimensions
In this appendix, we provide further details about deriving Eq. (28).In particular, inserting Eqs.(10) and (11) into the first equation of (21) with With the help of Eq. ( 17), this can be simplified to a,b,c,d We can relabel the dummy variables p a = p b + k, q a = q b + l, p c = p d + s, q c = q d + t so that the independent Clifford group matrix is separated as where we have used the equality τ pc,qc = τ p d ,q d τ s,t ω −q d s .Note that Eq. ( 88) is automatically satisfied for (k, l) = (0, 0).Therefore, for ∀(k, l) ̸ = (0, 0), the left hand side of Eq. (88) must vanish, namely, either 0, which gives the first line of Eq. ( 28) in the main text.The second line of Eq. ( 28) containing v 2 simplifies in a similar way using the second equation of (21).
The derivation of the L 3 condition follows exactly same procedure as above but may involve lots of tedious calculations.Here we just skip these details and give the result directly.With the parametrization method using Clifford group, Eq. ( 32) is equivalent to C Examples of L 2 Hierarchical circuits

C.1 Exhaustive parametrization for qubits
In the qubit case, the full parametrization of a two-qubit gate Eq. ( 23) is used with v 3 = v 4 = I.The Hierarchical condition Eq. ( 21) or Eq. ( 22) can be transformed into a series of algebraic equations of J x , J y , J z , r i , θ i , ϕ i with v i = exp{ir i (cos θ i σ z +sin θ i cos ϕ i σ x +sin θ i sin ϕ i σ y )}, i ∈ {1, 2}.However, these equations are very complicated and virtually impossible to be analytically solved.
Our strategy is to combine the analytical analysis with the numerical help of Mathematica.Those algebraic equations can be represented as a cost function f ≥ 0, designed to be minimized.The equations are solved for f = 0.The minimization is done by the Mathematica function NMinimize and we apply different constraints to it.For example, we require that all three Js are nonvanishing, i.e. not equal to nπ/2, n ∈ Z, and also avoid the dual unitary case.Extensive numerical experiments very strongly suggest that it is impossible to achieve a minimum value close to 0 under these constraints, leading to the conclusion that a vanishing J is necessary to satisfy the condition of L 2 .Consequently, we proceed by setting J x = 0 in the subsequent analysis.This process of simplification and analysis continues until the resultant algebraic equations become tractable for manual examination.

C.2 Non-exhaustive parametrization for qudits with D > 2
Here we provide some more details regarding examples with D > 2 mentioned in the main text.In particular, here we explain how to derive the simplified conditions that v 1 and v 2 in Eq. ( 28) need to satisfy.These simplified conditions help obtain v 1 and v 2 .In some cases, one can guess good candidates, and in other cases, these simplified conditions help in finding examples numerically.However, an explicit construction of all possible v 1 and v 2 is still lacking.
Given that the Clifford group matrix τ r,m is complete, the expression v * 1 τ k,−l v T 1 from the first line of Eq. ( 28) can be decomposed into We can look into the two terms separately.If the first term vanishes for all of {k, l} ̸ = {0, 0}, we just obtain the L 1 condition Eq. ( 18).Thus, a L 2 circuit would necessarily require that for some {k, l} ̸ = {0, 0}, the first term is non-vanishing.The set of these {k, l} can be directly determined with the specific form of θ p,q .The second term must vanish for this set of {k, l}s to satisfy Eq. (91).Since τ r,m is a complete basis in the matrix space, the necessary and sufficient condition for the vanishing of the second term can be expressed as: for every {r, m}, either α r,m = 0 or d θ * p d ,q d θ p d +s,q d +t ω p d m−q d r = 0, ∀{s, t} ̸ = {0, 0}.
(92) The latter one can also be directly calculated by substituting with the specific form of θ p,q .If it does not vanish for some {r, m}s, we obtain the condition that α r,m must vanish for those pairs.Since α r,m is related to v 1 , this serves as a necessary and sufficient condition v 1 should satisfy.
In conclusion, v 1 is implicitly defined by its conjugation action on those τ k,−l , with {k, l} determined from the first term in Eq. (91), such that v * 1 τ k,−l v T 1 has zero overlap with certain τ r,m , with {r, m} determined from Eq. (92).
We mention that this condition is more convenient than the original L 2 condition since we use the property of the Clifford group matrix.To construct v 1 satisfying the condition, one can either guess directly or use some ansatz of v 1 to solve the condition.The same procedure works for v 2 .
For the class given by the Eq.(29) of the main text, the above procedure results in the implicit equations for v i where i ∈ {1, 2} and the equation holds for ∀k, l = even, (k, l) ̸ = (0, 0) and r, v ∈ {0, D/2}.A naive example of v i is just the identity.In D = 6 we can use the ansatz v i = σ x ⊗ κ i where σ x is the Pauli-x operator.We numerically find that v i satisfies the above condition for almost any κ i ∈ SU(3).
For the class given by Eq. ( 30), the condition can be expressed as Tr v * i τ 0,l v T i τ 0,v = 0 for ∀v, l ̸ = 0 where i ∈ {1, 2}.A simple example that satisfies the condition is the generalized Hamamard gate H in higher dimension defined as HσH † = τ .D Proof of Eq. (31) In this subsection, we prove that the left equality in Eq. ( 21) along with unitarity implies the left equality in Eq. (31).The same argument holds for the right one.To prove this implication, we introduce the definitions By proving TrA † A = TrB † B, we establish that if A vanishes, then B also vanishes.This can be proved graphically.To this end, we introduce a four-folded gate given by = .The first term of TrA † A is The contraction represents a contraction between the same leg from the first and second layer as well as the third and fourth layer.Similarly, the contraction represents a contraction between the first and fourth layer as well as second and third layer.It is worth noting that the result remains the same if we exchange the second and fourth layers while simultaneously exchange and .Therefore, we arrive at (96) The R.H.S. just corresponds to the first term in TrB † B. By performing similar calculations, we can establish the agreement between each term of TrA † A and TrB † B.
The fact that a Hierarchical condition along with unitarity implies that its 180 • rotated version of the condition can be straightforwardly extended to higher-level Hierarchical circuits.The derivation presented here applies almost imediately.

E Numerical values of the parameters of the gates
In Fig. 2(b) and Fig. 4, we choose single site operators as v 3 = v 4 = I, v 1 = v 2 = v.The parameters for v = e ir(cos θσz+sin θ cos ϕσx+sin θ sin ϕσy) are listed in Table 1.
which are input and output legs of the gate, resulting in the dual local gate ũ.It is formally introduced by reshuffling the indices ũ = , ⟨j| ⟨l| ũ |i⟩ |k⟩ = ⟨k| ⟨l| u |i⟩ |j⟩ .(7)

Figure 1 : 2
Figure 1: The different Hierarchical conditions lead to different behaviors of 2-point correlation functions.a) L 2 Hierarchical circuits with non-vanishing correlators along three directions.b) General Hierarchical circuits that satisfy conditions for the k L/R level of Hierarchy from the left and right, respectively.Their correlators are non-vanishing in the middle part and on the light rays.The middle part is restricted to the inner light cone given by the velocity v satisfying −v L ≤ v ≤ v R , with v L/R following from k L/R level Hierarchical conditions.c) summary of the property of the two-point correlators.For a circuit without any Hierarchical conditions, we may equivalently consider it as k = ∞.

Figure 2 :
Figure 2: (a) The correlation function for D = 6 dimension, supporting on 2 sites.The gate is a non-ergodic member of L 2 with parametrization given by Eq. (29).(b) The correlation function for the qubit case supported on 3 sites.The gate is an ergodic element of L 2 , with parameters given in the App.E. The asymmetry between the left and right sides results from even/odd effects.In both figures, a i = b j = h, with h a random normalized traceless Hermitian operator.j is fixed at 20.The location of a multi-site operator is defined as the location of its left end.

Figure 3 :
Figure 3: This figure demonstrates the exponential decay of the correlation function along both the time axis (i.e., i = j) and the light cone (i.e., i = t + j).The two-qubit gate from L 2 , and the operators a and b are identical to the one used in Fig. 2(b).The solid line represents a linear fit of the exponential decay.

Figure 4 :
Figure 4: The correlation function of a random two-site observable following a quantum quench starting from the Bell state |Ψ L ⟩ defined in Eq. (63).The dynamics is governed by a L 2 circuit (symbols).We show its linear fit by solid lines.The blue and green data points show an exponential decay of the expectation value of observable, indicative of thermalization dynamics.In contrast, the red points remains constant at long times implying nonthermalization and non-ergodicity.The parameters for these three circuits can be found in App.E.

2 0 1 :
Table The left set of parameters is the one for Fig.2(b) and the blue line in Fig.4.The middle set corresponds to the green line in Fig.4and the right set corresponds to the red line in Fig.4.
clude a prefactor D in the definition to ensure that the autocorrelation function at time 0 is normalized to 1.Alternatively, this can be viewed as a quench from the b j 1 state, i.e. b j applied to the maximally mixed state.The correlations in the folded picture are graphically expressed as )The factor 1 D 2L comes from the normalized infinite temperature state ρ ∞ = I D 2L D 2L .We also in-5Hilbert Schmidt normalized means that Tra † i ai = 1