From dual-unitary to quantum Bernoulli circuits: Role of the entangling power in constructing a quantum ergodic hierarchy

Deterministic classical dynamical systems have an ergodic hierarchy, from ergodic through mixing, to Bernoulli systems that are"as random as a coin-toss". Dual-unitary circuits have been recently introduced as solvable models of many-body nonintegrable quantum chaotic systems having a hierarchy of ergodic properties. We extend this to include the apex of a putative quantum ergodic hierarchy which is Bernoulli, in the sense that correlations of single and two-particle observables vanish at space-time separated points. We derive a condition based on the entangling power $e_p(U)$ of the basic two-particle unitary building block, $U$, of the circuit, that guarantees mixing, and when maximized, corresponds to Bernoulli circuits. Additionally we show, both analytically and numerically, how local-averaging over random realizations of the single-particle unitaries, $u_i$ and $v_i$ such that the building block is $U^\prime = (u_1 \otimes u_2 ) U (v_1 \otimes v_2 )$ leads to an identification of the average mixing rate as being determined predominantly by the entangling power $e_p(U)$. Finally we provide several, both analytical and numerical, ways to construct dual-unitary operators covering the entire possible range of entangling power. We construct a coupled quantum cat map which is dual-unitary for all local dimensions and a 2-unitary or perfect tensor for odd local dimensions, and can be used to build Bernoulli circuits.

The complexity of interacting many-body systems, classical or quantum, matches their importance and prevalence that one cannot easily overstate. This is hardly surprising, given that even low-dimensional systems can display chaos that originates from the lack of sufficient constants of motion or nonintegrability. However, lowdimensional chaos is reasonably well-understood in terms of phase-space structures or spectral properties. Seemingly abstract toy models have played a crucial role in this process. The iconic one-dimensional logistic map has a long history with attributions to John von Neumann and Ulam from 1940 through early computer studies in the 1950's [1]. It came to a head with the more modern and definitive studies of May [2], Li, Yorke [3], Feigenbaum [4] and others, culminating in the discovery of universality in one-dimensional maps.
Closer to Hamiltonian systems, two-dimensional area preserving maps were studied as early as 1937 when Hopf introduced [5] what is now called the bakers map as an example of a mixing system [6], and we find the so-called Arnold cat map in a book by Arnold and Avez, published in 1968 [7]. The quantum avatars of the cat and bakers maps were introduced respectively by Berry and Hannay in 1980 [8], and by Balazs and Voros in 1987 [9, 10]. Such quantum maps [11] did provide vast simplifications over say a hydrogen atom in crossed electric and magnetic fields [12]. They were useful in understanding semiclassical methods and identifying classical structures in quantum states [13][14][15][16][17][18]. However, they were arguably not nearly as successful in being simple enough to penetrate the veils of quantum chaos, something which their classical counterparts did with ease.
Random unitary circuit models involving random two-body local interactions have been extensively studied in the recent past [19][20][21][22][23] as minimal models of many-body quantum nonintegrable systems. In such systems, the growth of bipartite enatnglement is found to be linear using a "minimal spacetime membrane" approach and is validated in the limit of large local Hilbert space dimension [20]. It is notable that quantum supremacy was purportedly demonstrated with a realization of such circuits in 2-D with superconducting junctions [24]. In addition, the introduction of a special class of unitary circuits of solvable lattice models that are typically nonintegrable is remarkable and hold potential as toy models of many-body quantum chaos [25][26][27][28]. The key aspect that make these models solvable is an underlying spacetime duality that is realized for certain interactions, hence these are referred to as dual-unitary circuits [26]. Quantities such as the two-point correlation functions, operator entanglements, spectral form factors, have become largely analytically accessible in these models [26,27,[29][30][31][32][33][34][35][36]. Using matrix product states as initial states, time evolution under dual-unitary circuits has been analyzed, including the evaluation of block entanglement entropy [37]. The OTOC are analytically tractable for dual-unitary circuits [38] and for random dual-unitary circuits [39]. The effect of perturbations that break the dual nature of the circuits is studied in [40] where the structural stability of these models is motivated as being analogous to that of classically chaotic systems [6]. Computational power of the circuits built from one and two dimensional dual-unitary operators are also studied recently and shows the quantum computational supremacy of these circuits [41].

A. Correlations and the Ergodic hierarchy
The dual-unitary circuits are defined in discrete time and can be interpreted as Floquet systems with a timeperiodic forcing. In fact the dual kicked chain model is a many-body Floquet system and the circuit model may be thought of as being equivalent to it [27]. They have been shown to be capable of displaying, in the thermodynamic limit, a range of properties that could be non-ergodic, ergodic or mixing. The mixing rates of these circuits can be found from a simple contractive channel that is a completely positive map, such as implemented in open quantum systems. The channel depends only on the two-body interaction and governs the spread of correlations across space-time [26]. These properties are inferred from the behavior of correlations such as C AB (t; x, y) = A y (t)B x − A y B x . Here operators A y (0) and B x are localized initially at y and x sites respectively, and A y (t) is time-evolved under the many-body dynamics.
For comparison, we recall the ergodic hierarchy of classical dynamical systems, which may be stated, largely, in terms of correlations [6,42,43]. Let f (x) and g(x) be two functions on phase space Ω, x = (q, p) ∈ Ω and dµ(x) be an invariant measure under the dynamics x → x t = φ t (x). For Hamiltonian systems the relevant measure is simply the Lebesgue measure or volume on phase space dqdp, in which case dµ(φ −t x) = dµ(x) from Liouville's theorem. Consider the correlation where for any function h(x), h(x) = Ω h(x)dµ(x), is its phase space average. Denote the time-average of any function h (t) as h (t) = lim T→∞ where u i and v i are single-particle operators, U and U are local unitary (LU) equivalent, and e p (U) = e p (U ) [48,49]. The dual-unitary property is preserved under local unitaries, thus U is also dual-unitary. However, the dynamical properties of circuits constructed from U and U are very different. For example, the two-point correlation functions and the mixing rates, depend on the local unitaries. The qualitative nature of the circuit changes with the use of local operators or "fields". This is not surprising as for example the Heisenberg model on the ring is integrable, however adding a nonhomogeneous local magnetic field can drive it to an ergodic phase. Similarly, local unitaries can enhance entangling power to typical values if sandwiched between nonlocal operators [51,52].
Nevertheless we show in this paper that the entangling power plays a central role in two respects: (1) it is shown that there exists a threshold such that if e p (U) > e * p the dual-unitary circuit built from it is guaranteed to be mixing, and (2) the local unitary averaged mixing properties of dual circuits, when the single particle operators u i and v i are sampled uniformly, are to a very good approximation solely determined by e p (U). The threshold e * p = (q 2 − 2)/(q 2 − 1) is a simple function of the local, single-particle, Hilbert space dimension q. For e p (U) ≤ e * p the circuit could be mixing, ergodic or non-ergodic. If e p (U) = 0 then the circuit cannot be mixing, and if and only if e p (U) = 1, the circuit is Bernoulli. Thus e * p signals a strong enough entangling power that no matter what the local fields are, the circuit is mixing. For qubits q = 2, it turns out that e p (U) cannot exceed e * p = 2/3, hence Bernoulli circuits cannot be constructed from qubits. However, for all q > 2 the entangling power can reach e p (U) = 1. The only case where it was not known if this was the case till very recently was q = 6, but an explicit construction of 2-unitary operators in this dimension has settled that question [54]. This is summarized in the following rudimentary phase diagram. Entangling power e p (U) The special value of the entangling power, e * p , appears as a natural border for various families of dual-unitaries, including qubits. We show explicit examples, including a quantum cat map, wherein e p (U) = e * p and the circuit is not mixing, therefore the bound for guaranteed mixing is tight. In fact e * p is the final threshold of a series of entangling powers that are multiples of 1/(q 2 − 1) that when breached ensure a minimum number of modes remain mixing whatever the local unitaries maybe. If e p (U) = 1, the U has maximum possible entangling power, and is 2-unitary by the definition in [47], then so is U and the Bernoulli property of circuits is local unitary invariant.

C. Various constructions of dual-unitary and 2-unitary operators
Using an ensemble with the Fourier gate (a dual-unitary operator) and different local unitaries, the dynamical properties of dual circuits were studied in [55]. Similarly, [56] constructed the dual-unitary operators from diagonal unitaries and the swap operator, to exhibit mixing circuits. However, we show below that the entangling power of such dual-unitaries do not explore the whole possible range. In [46], we constructed a non-linear iterative algorithm to produce ensembles of unitaries which are arbitrarily close to dual-unitaries and in some cases 2-unitaries and these provide examples with large entangling powers. Thus we are able to cover a wide dynamical range and indeed construct various families of dual operators and find their entangling powers.
Except for local dimension 2, the parameterization of dual-unitary operators is unknown. There are examples that have been constructed and studied in the past based on the Fourier transform (its duality has been discussed via its operator entanglement in [52,57]), complex Hadamard matrices, and diagonal unitaries. However, the systematic exploration of these must also consider their entangling power as this has a crucial impact on the nature of the circuit as we have motivated above. In this paper we provide various constructions of families of dual-unitary operators that are explicit and analytical and yield ensembles with a range of possible entangling powers.
We also provide analytical examples of 2-unitary operators that are to our knowledge new, and related to higher dimensional quantum cat maps. Such 2-unitaries imply the existence of corresponding absolutely maximally entangled (AME) states, wherein every bipartite cut is maximally mixed, and these are equivalent to perfect tensors with 4 indices [47,52]. The fact that there are no Bernoulli circuits for qubits is then related to the well-known fact from quantum information that there are no absolutely maximally entangled states of 4 qubits [58]. Thus the present work ties in quantum information theoretic measures with dual-unitary circuits and it is conceivable that it generalizes for other special circuits as well.

D. Structure of the paper
The structure of this paper is as follows. In Sec. (II) dual-unitary circuits are recalled and the roles of the matrix realignment and partial transpose-two central operations in quantum information, are developed. The Bernoulli circuits are constructed in Sec. (II B) and necessary and sufficient conditions for this property is discussed. It is also pointed out that correlations between operators with support on two sites retain the property of being Bernoulli. The explicit relation between entangling power and mixing in dual-unitary circuits is developed in Sec.(III). The mixing rates are defined and the sufficiency condition for the dual-unitary circuit to be mixing is derived in Sec. (III A). Local, single-particle, averaged spectral radius, which is related to the measure of mixing, is calculated and found to be a simple function of entangling power in Sec. (III B). Explicit analytical calculations of mixing rates for the case of two qubit unitary operators is carried in Sec. (III D).
Dual-unitary ensembles are constructed and their corresponding CPTP maps, involved in the calculation of correlation functions, are studied in Sec. (IV). The dual-unitaries constructed from block-diagonal unitary operators are in Sec. (IV A), the "dual-CUE" resulting from an iterative map acting on the circular unitary ensemble (CUE) is in Sec. (IV B). The designs from permutation matrices and the role of Latin squares is discussed in Sec. (IV C) and finally a quantized coupled cat map is constructed in Sec. (IV D). The cat map provides an analytically simple construction of dual-unitary operators in all local dimensions and 2-unitaries or perfect tensors in any odd local dimension. The latter is to our knowledge the first such construction of 2-unitaries valid in all odd-dimensions. Finally we summarize and discuss some possible directions in Sec. (V).

II. DUAL-UNITARY TO BERNOULLI CIRCUITS VIA OPERATOR ENTANGLEMENTS
In this section we review the dual-unitary circuits for being self-contained and motivate how this naturally gives rise to their being cast in information theoretic language. We point out that the two important rearrangements of matrix elements in quantum information theory that are used to detect bipartite separability of states-the realignment [59,60] and the partial transpose [61] show up prominently in the dual-unitary circuits. This allows us to introduce the entangling power that plays an important role subsequently. These considerations lead naturally to the definition of Bernoulli circuits as being at the apex of an quantum ergodic hierarchy.

A. Dual-unitary circuits, Realignment and partial transpose
Consider a periodic lattice of length 2L and each site, located for convenience at half-integer values, is occupied by a particle whose internal degree of freedom is a state in a Hilbert space H q of dimension q. The unitary U is a two-particle operator in H q ⊗ H q and acts on all pairs of nearest neighbor sites. It is convenient to use the diagrammatic approach as explained in [26]. Unitary operators U and U † acting on bipartite Hilbert space are diagrammatically represented as Connecting the legs is equivalent to contracting the indices, hence the unitarity of U, iα|UU † |jβ = δ ij δ αβ is diagrammatically The one-time evolution operator can be thought of as having two stages, in the first the 2L particle state is evolved under the unitary U ⊗L , acting on pairs (i, i + 1 /2), and in the other on (i + 1 /2, i + 1). Here T 2L |k 1 ⊗ |k 2 ⊗ · · · |k 2L = |k 2 ⊗ |k 3 ⊗ · · · |k 2L ⊗ |k 1 is the site translation operator. The evolution operator over time t is simply U(t) = U t . For example with L = 4, corresponding to 8 particles, U 2 is diagrammatically j=0 be an orthonormal operator basis in the single-particle local Hilbert space H q under the Hilbert-Schmidt inner product tr a † j a k = qδ j,k . If we choose a 0 = 1 q , it follows that tr a j = 0, for j > 0. The dynamical correlation function of single-particle operators in the infinite temperature state can be found from the following quantities (where i, j = 0) [26]: By reshaping of the matrix elements of U, its realignment R 1 (denoted byŨ in [26]) is defined as the permutation This can be used to define the evolution in the "spatial direction" by the transfer operator With periodic boundary conditions, these models satisfy tr U t = tr U L R and are in a sense dual to each other [25]. It has been shown in [26] that if U R 1 is also unitary, the correlations in Eq. (6) vanish inside the light cone |y − x| < t and the resultant circuit has been called dual-unitary . These circuits are built of dual-unitary operators.
This is represented diagrammatically as The realignment in Eq. (7) is intimately related to bipartite operator entanglement [62,63]. It has been pointed out [46] that the unitary U is dual iff the operator entanglement of U is maximum. Besides this, the realignment in Eq. (7) has important applications in investigating the complete positivity of quantum maps [64,65], and in characterizing the entanglement of mixed states [59,60]. To visualize the realignment R 1 in Eq. (7, consider a 4 × 4 matrix defined with the elements in the product basis as X iα,jβ = iα|X|jβ . Each 2 × 2 block of 4 × 4 matrix X is vectorized and arranged as a column in the X R 1 matrix. This can be visualized easily as: X 00,00 X 00,01 X 00,10 X 00,11 X 01,00 X 01,01 X 01,10 X 01,11 X 10,00 X 10,01 X 10,10 X 10,11 X 11,00 X 11,01 X 11,10 X 11,11     (10) X 00,00 X 10,00 X 00,10 X 10,10 X 01,00 X 11,00 X 01,10 X 11,10 X 00,01 X 10,01 X 00,11 X 10,11 X 01,01 X 11,01 X 01,11 X 11,11     (11) To see the connection to entanglement consider the four-party pure state on ⊗ 4 H q is a maximally entangled state, and U acts on the 12 subsystem. The partition 12|34 is maximally entangled as the reduced density matrix ρ 12 = UU † /q 2 = 1 q 2 /q 2 , due to the unitarity of U. The 31|42 partition has the state ρ 31 = U R 1 U R 1 † /q 2 , and if U is dual-unitary this partition is also maximally entangled [52]. The operator Schmidt decomposition of U is such that tr X † j X k = X j |X k = tr Y † j Y k = Y j |Y k = δ jk , that is the operators X and Y form orthonormal bases. The vectors |X j are columnwise vectorization of the q × q matrices X k , that is j|X k |i = ij|X k . It is easy to see that (X ⊗ Y) R 1 = |Y X * |, and hence Hence, the eigenvalues of ρ 31 are the Schmidt coefficients γ j ≥ 0 which determines the entanglement of U: iff only one of them is nonvanishing, the operator is of product form. As is clear from the Schmidt decomposition, the γ j are LU invariants under local unitary transformations of the form Eq. (3), and may be used to define a variety of entropies of entanglement. As ∑ q 2 −1 j=0 γ j /q 2 = 1, {γ j /q 2 , 0 ≤ j ≤ q 2 − 1} defines a discrete probability distribution, in terms of which the entropies maybe defined. For our purposes, the linear entropy, the operator entanglement [49] defined as where 0 ≤ E(U) ≤ E(S) = 1 − 1 q 2 , is useful. Here S is the swap operator in Eq. (2). From this, the operator entanglement E(U) is maximum iff U R 1 is also unitary, and all γ j = 1. Thus the maximal operator entanglement is equivalent to the unitary operator being dual-unitary, a property that is LU invariant as E(U ) = E(U). The OTOC for bipartite systems [66,67], averaged over all local unitary operators was shown to be related to the Loschmidt Echo [67]. The averaged OTOC was related to the operator entanglement of bipartite unitary operators [68], and provides another operational meaning to the operator entanglement by relating it to the measures of chaos in bipartite systems.
Returning to the circuit, using the fact that U is invariant under two-site shifts, UT 2 2L = T 2 2L U, to find the dynamical correlation function in Eq. (6), it is sufficient to consider the pair of functions . These two correlations are when the single-particle operators are to the "left" or "right" leg of a two-particle U. The only non-zero correlation for dual-unitary circuits lie along the light cone and was shown, quite remarkably, in [26] to be given in terms of a completely positive trace preserving (CPTP) map constructed only from the two-particle unitary operator U: This can be seen from the diagrammatic evaluation of the correlation functions C ij The CPTP maps M ± are M + (a) = 1 q These linear maps take density matrices to density matrices and are unital (preserve the identity). Let the action M + (ρ) = ρ be denoted as M + [U]|ρ = |ρ , where |ρ is the q 2 dimensional column vector formed from its rowvectorization, and the dependence on the unitary U is made explicit. That is, jl|ρ ≡ j|ρ|l and However, as α|M + (ρ)|β = α|ρ |β = αβ|ρ we identify Similarly, we get The reordering of matrix elements used in the above equations are defined as follows: Partial transpose T 1 : Partial transpose T 2 : To visualize these operations, consider a 4 × 4 matrix with matrix elements X iα,jβ = iα|X|jβ as given in Eq. (10).
The effect of local operators is as follows: First, we recall some basic properties of the channel. M ± is a Liouvillian superoperator representation of the channel M ± . M ± is a non-Hermitian operator whose eigenvalues {λ j } are bounded λ j ≤ 1 [26]. The channel M ± is unital: it fixes the maximally mixed state ρ = 1 q /q, thus there is always one eigenvalue = 1, the "trivial eigenvalue" corresponding to the (unnormalized) eigenvector |ρ = |1 q /q = |Φ + / √ q, where |Φ + is the maximally entangled state from Eq. (12). We note that M + (a) is complementary [69][70][71] to the channel where the subsystem initially in the most mixed state is traced out, which is the traditional process of tracing out a complex environment. The channel E + (a) is obtained by replacing U in M + (a) by US. The Liouvillian superoperator of the channel E + (a) is [62], and if U is dual-unitary, this channel is 1 Here the result from Eq. 26f is used. The spectrum of M + arranged in increasing order of its magnitude is used below as λ 0 = 1 ≥ |λ 1 | ≥ |λ 2 | ≥ · · · ≥ |λ q 2 −1 |. Generically the eigenvalues of both M ± will be referred to as λ k below. The authors of Ref. ([26]), classified the ergodic properties of the dual-unitary circuits from the nontrivial eigenvalues (k > 0) of M ± as Non-interacting: All 2(q 2 − 1) eigenvalues λ k = 1, the correlations remain constant.
Ergodic and non-mixing: All nontrivial eigenvalues λ k = 1 but there exists atleast one eigenvalue with unit modulus |λ k | = 1.
following the discussion of Bernoulli systems in the introduction, we add Bernoulli: All 2(q 2 − 1) nontrivial eigenvalues vanish: λ k = 0.
Bernoulli implies mixing, but not the other way around, and we will motivate below that this is indeed a unique class by itself and lies at the apex of a possible quantum ergodic hierarchy.

B. Bernoulli circuits, entangling power and 2-unitarity
To begin, note that the correlations on the lightcone vanishes if Eq. (18) and (19) reduces as follows as tr(a) = 0. This motivates the introduction of "T-duality" based on the T 2 partial transpose, see Eq. (23). Although we do not imply any evident duality as such, it is convenient to refer to this in parity to the definition of dualunitarity.
Definition 2 (T-dual unitary) If U and U T 2 are both unitary, U is T-dual. Diagrammatically, Using U T 2 , similar to defining E(U) in Eq. (15) it can be shown [52] that In the context of the 4-party state in Eq. (12), ρ 14 = U T 2 U T 2 † /q 2 is the state of the 14|32 partition and E(US) is the linear entropy of entanglement of this bipartite split. The quantity that in a way brings together measures of duality and T-duality is the entangling power e p (U). This is defined as the average entanglement produced by the action of a unitary operator U on an ensemble of pure product states distributed according to the Haar measure, where E (·) can be any entanglement measure and C q is scaling factor. By considering the linear entropy as an entanglement measure, the entangling power e p (U) of U has been related to the operator entanglements E(U) and E(US) [49]. The normalized entangling power 0 ≤ e p (U) ≤ 1, with C q = (q + 1)/(q − 1) is then: where E(S) = q 2 −1 q 2 has already been introduced. Note that this choice of scaling is different from that used in [49], but is convenient to have the maximum value be 1. The entangling power e p (U) is maximum iff both U R 2 and U T 2 of U is unitary. Thus 2-unitary operators have been previously defined [47].

Definition 3 (2-unitary) If U is dual-unitary and T-dual then it is 2-unitary.
Note that if U and U R 1 is unitary then U R 2 is also unitary and vice versa. Similarly, for T-dual property, if U and U T 1 is unitary then U T 2 is also unitary and vice versa. Hence there is no need to distinguish between the two realignments for defining the dual property, or to differentiate the two partial transpositions in defining the T-dual property. The two classes of dual-unitary operators are connected by the swap operator. In particular if U is T-dual (dual), then US and SU are both dual (T-dual) unitary operators. This follows as from Eq. (26d) (SA) R 2 = SA T 1 , hence if A is T-dual, SA is dual-unitary, as the swap is unitary. Similarly, AS is also dual-unitary . As simple examples, the identity is T-dual, and swap S is dual, cnot is T-dual and dcnot= S × cnot is dual. The operator entanglements E(U), E(US) and the entangling power e p (U) are all local unitary invariants (LUI).
If U is 2-unitary, then U T 2 is also unitary and hence from Eq. (21) Hence the nontrivial part of the spectrum of M + [U] vanishes, that is λ k = 0 for k > 0, and the correlations decay instantaneously. Local unitary transformations do not alter the situation as 2-unitarity is LUI and hence . Note that in this case the spectrum of M − [U], from Eq. (22) is also identical with the only non-zero eigenvalue corresponding to the trivial mode. Bernoulli circuits are not "chiral". Now assuming that λ k = 0 for k > 0, we show that U T 2 is unitary. The condition on the eigenvalues of M + [U] implies that it is rank 1, and hence we can assume a form M + [U] = |φ 2 φ 1 |, where |φ 2 is the (right) eigenvector corresponding to the eigenvalue 1 if φ 1 |φ 2 = 1. However, this is the trivial eigenvalue and the corresponding normalized eigenvector from the fact that the map is unital is |1 q / √ q, where |1 q = |Φ + . Considering the left eigenvector corresponding to the eigenvalue 1, we conclude that |φ 1 is the same as well.
, assuming that its only nonzero eigenvalue is the trivial one, leads to U T 2 † U T 2 = 1 q 2 and hence U T 2 is unitary as was to be shown. As U is assumed to be dual, it is in fact 2-unitary.
If U is 2-unitary the 4-party state in Eq. (12) is such that all the bipartite splits are maximally entangled, in other words it is an absolutely maximally entangled state (AME) [47,72]. It is denoted as AME(4, q) to indicate that there are 4 particles each of q dimensions. It is well-known that there are no AME(4, 2) states: 4 qubits cannot form an AME state [58]. This implies that there are no dual-unitary circuits of qubits that are Bernoulli. It is also known that AME(4, q) exists for all q, with the possible exception of q = 6. This has been demonstrated explicitly by constructing such 2-unitaries as permutations based on orthogonal Latin squares [73]. We will return to this when constructing dual operators from permutations further below.
The classification of a circuit based on the maps M ± is valid for single-particle observables, but may not be in the same class for observables with larger support over a few particles, or "extended observables". While a study of such correlations in dual-unitary circuits is, to our knowledge, absent, at least in some limited way Bernoulli circuits seem straightforward examples where correlations of two-particle observables continue to vanish. To illustrate this, the dynamical correlation function between the operators defined on the two sites is shown in Eq. (33). As any two-particle operator can be written as a linear combination of such product operators a k , it is sufficient to consider the correlations from such pairs. Again, as for single particle, we adopt an orthonormal basis with identity as one of them and consider the correlation function with By using the duality conditions, Eq. (9), the correlations vanish inside the lightcone and is shown in the Appendix (A).
For the correlations along the lightcone, x 1 = t − 1 2 , x 2 = t, the repeated application of unitarity from Eq. (4), Applying the T-duality condition Eq. (28) to the highlighted part in Eq. (34) and repeating its application, the vanishing of the correlation function is evident, as The vanishing correlation function for the case x 1 = x 2 = t − 1 2 is similar and shown in Appendix (A). More general cases of extendable observables need to be studied, but we do not expect them to vanish in general. This is also in line with classical Bernoulli systems where the Bernoulli partition is singled out for independence [43]. There exists unitaries even with q = 2 such that M k + [U] = |Φ + Φ + | for some k > 1, so that all single-particle correlations decay to zero beyond t = k. However, unlike the Bernoulli case, local transformations can change the nature of this correlation and it need not continue to hold even for two-site observables. We will construct examples of Bernoulli and such cases using cat maps, among other possibilities. Finally we observe that if U is not a dual-unitary operator, but is T-dual alone, the correlations along the lightcone still vanish for t = 0.

A. Rates, and a sufficiency condition for dual circuits to be mixing
To separate the trivial eigenvalue so that the largest eigenvalue is a nontrivial one, we definẽ The eigenvalues ofM + [U] are the same as that of M + [U], except that the trivial eigenvalue is replaced by 0. This follows as Φ + | is the left eigenvector with eigenvalue 1 as may also be verified directly using the form Eq. (21) for example. Therefore we can peel off the projector on the maximally entangled state |Φ + Φ + | from a spectral or Jordan decomposition of M + [U]. The remaining states are orthogonal or can be chosen to be orthogonal to |Φ + andM + [U]||Φ + = 0.
Focusing on M + [U], the non-trivial eigenvalues determine the decay rates If any |λ k | = 1, these correspond to non-decaying modes, on the other hand there is instant mixing along a mode if the corresponding λ k = 0. Note that |λ k | ≤ 1 and hence µ k ≥ 0. Arranging the decay rates in an increasing order µ 1 ≤ · · · ≤ µ q 2 −1 , ignoring the trivial mode as µ 0 = 0, we may define µ 1 as a degree of mixing. It is derived from the largest magnitude eigenvalue (the spectral radius |λ 1 |) ofM + [U] and vanishes when any one of the non-trivial modes does not decay. This is infinite if all the µ k = ∞ for k > 0, which is the case if the circuit is Bernoulli. The converse is also true: if µ 1 = ∞, then the circuit is Bernoulli. For the most part, we will concentrate on estimating the largest nontrivial eigenvalues |λ 1 | which controls these rates as well as the spectral gap 1 − |λ 1 |.
This proves the earlier observation that the singular values of U T 2 and the eigenvalues of U T 2 U T 2 † , are local unitary invariant. However their realignment, which results in the CPTP map (18), does not have this property, and hence the decay rates depend on the locals used.
Therefore the eigenvalues of M + [U] are not LUI. We note that the eigenmode corresponding to the trivial eigenvalue of 1 remains |Φ + as (u ⊗ u * )|Φ + = |Φ + , where u is any local unitary and u * is its element-wise complex conjugate, which follows from the vectorization of the identity uu † = 1, as |ABC = (A ⊗ C T )|B . As µ 1 is dependent on the local operators, we may define the maximal mixing rate associated with U as which is the maximum possible mixing rate under the action of arbitrary local or single-particle unitaries. A similar maximal mixing rate can also be defined from M − [U]. Another useful measure is the average mixing rate where the u i and v i are independently drawn from the Haar measure on the local single-particle spaces. One may also first perform the average on the spectral radius |λ 1 | before finding the rate. All of these quantities, measuring the mixing induced by the channel, are by construction LUI and we address their connection to LUI measures based on U itself, in particular to the entangling power e p (U). The primary connection is to the norm of the channel. The Eq. (39) implies that this norm is a LUI, as where A = tr(AA † ) is the Hilbert-Schmidt or Frobenius norm. To see this, use Eq. (21), and the easily verified As U is dual-unitary, E(U) = E(S) and from Eq. (30) e p (U) = E(US)/E(S). Using the relation in Eq. (29), we then arrive at the above connection, which is also the norm of M − [U]. LetM + = VTV † be a Schur decomposition [74] ofM + , where V is unitary and T is an upper-triangular matrix whose diagonal entries are the eigenvalues λ i . It then follows from the unitary invariance of the norm that M + 2 = T 2 and hence These consequences are immediate: , it follows that the smallest eigenvalue is bounded from above as 2. An upper-bound for the largest nontrivial eigenvalue of M + or its spectral radius |λ 1 | immediately follows which also gives a lower-bound on the spectral gap 1 − |λ 1 |, and the rate of mixing.
3. The above two are subsumed in a series of inequalities which follow from the observation that The upper-bound in Eq. (45) implies that whenever e p (U) > 0, |λ q 2 −1 | < 1 and there will be at least one mixing mode and hence iff U is the swap or locally equivalent to it, all the |λ k | = 1 as e p (U) = 0. Thus the presence of mixing modes is generic. On the other extreme, if e p (U) = 1, |λ 1 | = 0 and hence all nontrivial λ k = 0 and we are in the Bernoulli circuit case, consistent with our previous observations. We observe that the upper-bound for |λ 1 | in Eq. (46) is trivial, as it will not be less than 1, unless This is then a sufficiency condition for the dual-unitary circuit to be mixing as it ensures that |λ 1 | < 1, and hence |λ k | < 1 for all k > 0, and for both M + and M − . The inequality in Eq. (46) is a tight one, that is there exists dual-unitary gates U such that e p (U) = e * p = (q 2 − 2)/(q 2 − 1) and |λ 1 | = 1, in fact λ 1 = 1, as we show explicitly via examples when discussing specific instances. The qualitative content of the sufficiency condition is plain: if the interaction as measured by the entangling power of the unitary bricks of the quantum circuit is sufficiently large, the resultant circuit is always mixing, whatever be the local fields. For the case of qubits, q = 2, the sufficiency condition e p (U) > 2/3 is never met, as the maximum possible entangling power is 2/3 and 2-unitaries do not exist, as we have already noted. Gates such as dcnot [51,52,75] have the maximum entangling power of 2/3, and circuits built of them need not be mixing depending on the local operations used.
While the smallest eigenvalue |λ q 2 −1 | < 1 for any entangling unitary U, and sufficiently entangling ones ensure mixing, there is a hierarchy of intermediate cases. It follows from Eq. (47) that for 1 ≤ k ≤ q 2 − 1, where e * p,1 ≡ e * p and e * p,q 2 −1 = 0. Thus the influence of the increasing entangling power of the two-particle unitary U on the mixing properties of the many-body circuit is encapsulated by this periodic gradation: if e p (U) > e * p,k there are at least q 2 − k mixing modes. This is illustrated in Fig. (1).
The entangling power does encode at all levels, the degree of mixing in the dual-unitary circuits. Further, as e * p is a tight bound, so also the others maybe, in the sense that there exist U such that e p (U) = e * p,k and λ 1 = · · · = λ k = 1. The rather large upper-bound on the spectral radius Eq. (46) is due to the factor q 2 − 1, which we now go on to whittle down.
As the entangling power of the two-particle unitary U is increased, by multiples of 1/(q 2 − 1), beyond e * p,k , this ensures that q 2 − k modes of the dual-unitary circuit are mixing. At least one mode is mixing if e p (U) > e * p,q 2 −1 = 0 and all the q 2 − 1 nontrivial modes are mixing if e p (U) > e p,1 = e * p .

B. Single-particle averaged spectral radius and rates
We wish to go beyond the bounds such as in Eq. (46) to find approximations. Local unitary operations (37) and can substantially lower it as well. We will motivate, both numerically and analytically that µ + [U], the average mixing rate, is an universal function of the entangling power e p (U), independent of the local dimension for q > 2. Analytical results include the case of q = 2, which is different from the rest, as well as the approximation that arise from which is seen to be qualitatively valid and quantitatively so near e p (U) = 1. Thus the averaging effectively removes the q 2 − 1 factor from the upper bound. It is quite remarkable, that while e p (U) is just one of many local unitary invariants, it plays a central role in the mixing properties of the quantum circuit built of U.
The spectral radius of an arbitrary matrix A is the largest absolute value of its eigenvalues. A formula of Gelfand is well known for the spectral radius and is given as the limit A k 1/k as k → ∞, where A is any norm [74], but we will continue to imply the Hilbert-Schmidt norm. Thus we have |λ 1 | = lim k→∞ M k + 1/k . For simplicity we will writeM + for what we callM + [U]. We have that the norm M + is a LUI and is given by Eq. (42). However the norm of M k + for k > 1 is not LUI. Consider Eq. (39), and identify u ≡ v † 2 and v ≡ u † 1 . We then have that the local-unitary Haar averaged spectral radius ofM + is The invariance of the Haar measure under multiplication by a fixed unitary is used, as well as the limit and the integration is interchanged.
To begin with consider To perform the Haar average (A(u ⊗ u * )B) 2 use the identity (see Appendix B for a derivation) which uses the results from Ref. [76][77][78] for calculating Haar average of monomials over the unitary group With X = A † A and Y = BB † and the identification of and hence an exact evaluation, where we have used Eq. (42). This then implies in the spirit of Eq. (46), a local unitary averaged inequality: which implies For the mixing rate, This gives a sufficiency condition for the circuit to be mixing as in the sense that there will exist locals that will make the circuit mixing. However, these are weak conditions and it appears from numerical results that only e p (U) > 0 is required for suitable local unitaries to make the circuit mixing. More importantly, a comparison of the inequalities in Eqs. (46) and (56), shows how the local unitary averaging reduces the effect of the factor (q 2 − 1) at even k = 2. We now indulge in an approximate evaluation for the norm of the higher powers ofM + . Continuing to denote Haar averages by overbars where B = [A(u ⊗ u * )] k−1 . For k > 2 the Haar average seems hard to perform exactly and we resort to treating the equality in Eq. (53) as remaining valid even in this case although B = (A(u ⊗ u * )) k−1 depends on u. Similar reasoning as before, including the fact that u ⊗ u * leaves |Φ + invariant, implies still that tr Y R 2 = 0. But now we can recursively separate out one more factor of A(u ⊗ u * )B from B to average B 2 . We get the estimate and hence finally, with A =M + and using Eq. (42), Using the Gelfand formula leads to an estimate of the local unitary average of the largest nontrivial eigenvalue as where f q ≈ 1 is a potential "fudge factor". Thus the (q 2 − 1) factor in the norm ofM + [U] is essentially eliminated by the local unitary averaging procedure with the exactly same factor appearing in the Haar average formula of Eq. (52). Shown in Fig. (2) are eigenvalues λ 1 and λ q 2 −1 of M + [U ] when an ensemble of local unitaries have been applied to fixed dual-unitary U, see Eq. (39). Except for the case when e p (U) = 0.75, these are the largest and smallest nontrivial eigenvalues for 3 realizations of U = D 1 S when D 1 are themselves diagonal matrices of random phases uniformly distributed on the unit circle. This is an ensemble of dual-unitary matrices, which we discuss in the section (IV A) as the D 1 S ensemble. Its entangling power does not exceed q/(q + 1). The realizations have been chosen to have a representative range of values of the entangling power, which is also displayed. The case when e p (U) = 3/4 is a dual-unitary obtained from an algorithm described in Sec. IV B. We will presently discuss the reason for this choice.
In the figures a circle of radius 1 − e p (U) is also marked. For small e p (U) the |λ 1 | is large, but so is |λ q 2 −1 | and both nearly coincide with the 1 − e p (U) circle. In all cases it is clear that the inequality in Eq. (45) holds. We notice also that there are substantial number of real eigenvalues, apart from the one at 1: this arises due to a symmetry in the spectra that if λ is an eigenvalue then so is λ * . There is a similar "equatorial" enhancement for the real Gaussian random or Ginibre ensemble [79]. As the entangling power increases, the separation between the eigenvalues is also more prominent, but the ring of λ 1 values is essentially governed by the same circle, and indicates that that Eq. (62) is valid. This is shown in Fig. (3) where quite remarkably the simple formula in Eq. (62) fits reasonably well the average of |λ 1 | for a variety of ensembles of dual-unitaries in a large range of entangling power, with f q = 1.1 for 3 ≤ q ≤ 6. The dual-unitaries used are from two sources: first is the D 1 S ensemble used in the last figure, whose entangling power is limited. The other is a numerical iterative procedure that generates an ensemble dubbed dual-CUE [46] and which we recall in the next section. It is capable of producing high entangling power dual-unitaries including 2-unitaries. Perturbations of such 2-unitaries are used to generate unitaries with entangling power arbitrarily close to 1 to populate the inset of the figure, which shows how well the Eq. (62) works with f q = 1 near e p (U) = 1.
Thus the factor f q appears to be a slow function of the entangling power. Note that we have made use of q = 4 for the inset, as q = 3 is a special case wherein e p (U) = 1 appears to be isolated in the set of dual-unitary matrices.
Numerical results indicate that for q = 3 there is a gap whose left extreme seems larger than e * p , but is otherwise not yet known [52].
The average mixing rate µ + and the maximal rate ν + for dual-unitaries with entangling power close to e p (U) = 0 for q = 2, 3 and 4 are shown in Fig. (4). These dual-unitaries are of the form D 1 S, D 1 is diagonal with entries e i φ where φ is uniform and random in [−π, π) and controls how far from the swap the dual-unitary is. Firstly, it is remarkable that for q = 2, ν + seems to lie exactly on the e p (U)/2 line. We shortly derive the analytical form of ν + for qubits, which indeed is a singular case. For q = 3, 4 the ν + has a slightly smaller slope than q = 2 and both lie on the same line. The µ + line for q = 2 has a slope of ≈ 3/8, while for q = 3, 4 we find a scatter of points. We discuss complementary data and comparisons with analytics later in Figs. (6) and (9).
The scatter of the average mixing rates near e p (U) = 0 signals the involvement of other local-unitary invariants than the entangling power alone. That these are not fluctuations due to sample size was verified by increasing the sample size by several orders of magnitude and observing that the scatter remained practically invariant. To examine this further, we take gates U 1 and U 2 with the same entangling powers: e p (U 1 ) = e p (U 2 ), but which are not local-unitarily connected, that is there exists no u i and v i such that (u  Fig. (2) shows yet another instance of e p (U) = 3/4 but is not locally equivalent to either U 1 or U 2 . For details of these matrices, see the Appendix (C). The fact that the λ 1 is very differently distributed in the unit circle itself indicates that these operators are not locally connected. It is found that E u i ,v i (|λ 1 |) ≈ 0.48655 for U = U 1 whereas it is ≈ 0.55045 for U = U 2 with N s = 10 6 samples of local unitaries, and seem to have converged to an accuracy of at least 10 −3 .
The peculiar structure seen in Fig. (5) corresponding to U 1 warrants some comments, especially as some aspects of the "triangular features" are also evident in two of the cases in Fig. (2). Very similar figures have appeared in the study of so-called unistochastic matrices [80]. The map corresponding to M + [ (39)) is special in this case and apart from a null 6-dimensional space, it is given by the bistochastic 3 × 3 matrix whose elements are |u ij | 2 , where u = u † 1 v † 2 , as shown in Appendix (C). A bistochastic matrix is unistochastic if its elements are derived as the absolute value squared of the elements of a unitary matrix. The eigenvalues of random unistochastic matrices for q = 3 are conjectured to be bound by the deltoid (3-hypocycloid) [80] which are the origin of the structure observed in the first part of Fig. (5). It must be emphasized that the CPTP maps M + [U ] are in general not related in such a manner to unistochastic matrices and hence such structures are generically absent. In particular we have found them for specific subsets of CPTP maps derived from special dual unitaries with special values of the entangling power (such as 3/4 in q = 3).

C. Inhomogeneous circuits
If the two-particle unitary operators are different, but still dual-unitary, within a circuit, this represents an inhomogeneous case and in this brief interlude we use the results above to bound the correlation decay. The CPTP maps after t time steps are Here is the interaction between site k − 1 and k. The correlation function from Eq. (16) becomes in this case, considering traceless basis operators, where we have used the Schwarz inequality and tr a † j a j = a j |a j = q. The replacement of M ± byM ± , a product ofM ± [U k ] is possible as the a j are traceless.
Allowing for single particle unitaries that are also inhomogeneous (u 1k etc. are i.i.d. random) and averaging over them leads to the inequality Here we note that the averaging in Eq. (61) is exact as the local unitaries are independent at every stage of the averaging. As a special case, assuming single particle inhomogeneity, but uniform two-particle interactions, we have where γ = − ln 1 − e p (U) 2 . Thus the decay of correlations in the quenched disordered dual circuits are exponentially fast and the rate is explicitly given in terms of the entangling powers of the two-particle unitaries used.

D. Qubit case analytics
The qubit case is special in many ways including the fact that the dual-unitary family has only one free parameter. This presents a completely solvable case to find the connections to the entangling power when a restricted volume of local unitary operations are used. All local unitary inequivalent classes of two-qubit operators are characterized by their Cartan decomposition that is a three parameter family [81][82][83]. For dual-unitary subclass, this is reduced to the one parameter family where σ j are Pauli matrices and 0 ≤ J ≤ π/4. Note that and the case J = π/4 corresponds to a swap gate while J = 0 corresponds to the dcnot [83]. A general local or single qubit gate u ∈ SU(2) can be parametrized as, where, θ ∈ [0, π], and φ, ψ ∈ [0, 4π]. We wish to find the eigenvalues of the channel corresponding to a locally equivalent U : M + [U (J)] = (u † ⊗ u T )M + [U(J)]. This turns out to be analytically difficult for us and therefore we take single qubit gates from two different subsets. These are found to be useful, as the first gives an approximation to the average mixing rate µ + , while the second leads to, what appears numerically to be the exact analytical expression for the maximum mixing rate ν + .
The first subset is parametrized as, obtained by setting φ = 0 in Eq. (69). The maximal mixing rate is where e max p = 2/3, the maximum possible value for qubits (See appendix E 1). The prime on ν + is to indicate that the maximum rate is found within a subset of all possible single-particle unitary gates. The connection to the entangling power is via Eq. (68) and displays explicitly its monotonic relation to the maximum mixing rate. Thus ν + [U] ≥ ν + [U], and comparisons with numerical results are shown in Fig. (6). The average mixing rate µ + , again with the restricted set of local unitaries drawn from w, with cos θ being uniformly distributed, is after taking into account the non-smooth changes at θ c and π − θ c . Again the monotonic relation between µ + [U] and e p [U] is made explicit, and agrees with numerical results as shown in Fig. (6). We reiterate that the expressions for µ + and ν + are valid only for the given parametrization of the locals (Eq. 70), and in particular do not agree when all the unitaries are sampled using the CUE. The latter sampling that uses general locals u ∈ SU(2) of the form in Eq. (69) is also shown in Fig. (6). It is interesting that the average mixing rates are approximately the same from both sampling. For the case when J = 0, when the entangling power of U is the maximum possible value of 2/3, µ + [U(J = 0)] = 1. The CUE average in this case is also easy to find and µ + [U(J = 0)] = 1 as well. This reflects the fact that there are no 2-unitary gates for qubits, as for q > 2, µ + [U] = ∞ when e p (U) = 1.
For small e p (U), where from Eq. (72) we infer that µ + [U] ≈ 3e p (U)/8. From Fig. (4) we see that with CUE sampling, µ + [U] ∼ 0.36e p (U) which is indeed close to this. However, Eq. (71) gives ν + [U] ∼ 3e p (U)/8 which is quite different from the numerically obtained estimate of e p (U)/2 from Fig. (4). A better bound on ν + [U], which we believe to be the exact value, is attained with the second subset of local gates of the form, obtained by setting θ = π/2 in Eq. (69). The maximal mixing rate, where the maximum mixing rate is found over single or local qubit gates of the form Eq. (73), and is also written in terms of e p (U) = 2 3 cos 2 (2J) (See appendix E 2). Indeed this is a better bound than the one in Eq. (71), however as this is still based on a subset of local unitaries it is not a proof that this is indeed ν + . In Fig. (6), Eq. (74) is plotted and fits exactly the ν + [U] calculated from general single qubits sampled from the CUE. Further, calculations and numerical results shown in appendix ?? strongly suggest that Eq. (74) is a global optimization using single qubit gates of the general form Eq. (69) from the entire SU(2) space. Notice also that Eq. (74) gives for small e p (U) that ν + [U] ∼ e p (U)/2, agreeing with the corresponding numerical result displayed in Fig. (4). Note that the results for qubits are not surprising as there is only one parameter J and the invariants are all function of these. However, the monotonic dependence of the mixing rates with the entangling power does not follow from this and is a forerunner of the predominantly monotonic dependence for q > 2.
We now construct dual-unitary, and 2-unitary, operators in a variety of different ways, keeping in view the resulting entangling power which we have just demonstrated plays a key role in the ergodic life of a dual-unitary circuit.

IV. DUAL-UNITARY CONSTRUCTIONS AND THEIR ASSOCIATED CPTP MAPS
The channel M ± [U] and the influence of the entangling power of U on its mixing properties have just been discussed in some generality. We now step back to describe families of dual-unitary operators which are the building blocks of the circuits in the first place. There have been few constructions in the literature already [55,56], we augment them considerably, especially with the understanding that their entangling power forms a crucial aspect to control. dual-unitary operators are also those with maximum operator entanglement and while there exists no general parametrization of these, there are few concrete families including a numerical algorithm constructed by the current authors, which results in the "dual-CUE" ensemble [46]. However, the dual-CUE have typically large entangling power and it may be good to construct analytically families that may have smaller entangling power but are parameterized well. For large entangling power and 2-unitaries there are also permutation matrices and complex Hadamard matrices that are quantizations of classical cat maps: paradigms of classical chaos.

A. From block-diagonal unitaries
Dual-unitary matrices are those that remain unitary under realignment. However, of the four permutations involving the partial transpose and realignment, the simplest to visualize is T 2 . Therefore we may construct first T-dual matrices and left or right multiply these by the swap operator S to get dual matrices. For the symmetric q × q spaces, A T 2 amounts to partitioning A into q 2 q × q matrices, the first block starting at the first element, and taking the transpose within these blocks. As the transpose of a unitary matrix remains unitary, this implies that any matrix comprising uniformly of only diagonal, unitary, q × q blocks is T-dual. Denoting this as D q , we have where {U j } is any set of q × q unitary matrices. Hence are two families of dual-unitaries, parametrized by q matrices ∈ U (q).

Diagonal unitary constructions
The blocks U j may in turn have a block structure, an extreme case being when they are all diagonal. In this case, we refer to the diagonal matrix as D 1 and the ensemble SD 1 was studied recently in this context in [56], but was studied from the point of view of entanglement content in [85]. The channel M + [U] depends on which side the swap operator is multiplying the block diagonal matrices, but the entangling power does not. The entangling power of e p (U = D 1 S) is essentially the operator entanglement of diagonal gates. If D 1 is the uniform ensemble of diagonal unitary gates (the phase in the diagonal elements is i.i.d. random, uniform in [−π, π)) moments of the singular values of D R 2 1 maybe derived [85]. For example we find from there that the average entangling power is This compensates for an additional factor of (q − 1)/(q + 1) present in [85] due to different scaling used in the definition of entangling power. Recall that in the present convention, the maximum value of the entangling power is 1. However, we show below that the maximum value of e p (U = D 1 S) is only q/(q + 1), see Eq. (82). Hence this is ≤ e * p , for q > 2 and = e * p for q = 2, where e * p is the border for ensuring mixing in Eq. (48). For q = 2, the distribution of the entangling power may also be explicitly found, as it is observed in [85] that the distribution of the squared singular values of U R 2 follows the arcsin law. This implies the distribution as seen in Fig. (7), where for comparison the case q = 3 is shown as well. These show the members of these dual-unitary family can have a fair spread of the entangling power from nearly 0 to q/(q + 1).

General block diagonal constructions
A more general ensemble is constructed from where 1 < K ≤ q and U j are special classes of q j × q j unitary matrices, where q j are multiples of q such that ∑ K j=1 q j = q 2 . For example when q = 3, we could have K = 2 blocks with q 1 = 6 and q 2 = 3. The extreme case of K = 1 simply corresponds to a full T-dual matrix of q 2 dimension, there is no block structure and is avoided here. The case when K = q is the uniform one discussed above.
Only when q j = q, U j is an arbitrary unitary. If q j = m j q, with 1 < m j < q, U j comprise of m 2 j q × q sub-blocks. These U j may be thought to belong to a tensor product space H m j ⊗ H q , of a m j dimensional system with a q dimensional one. If all such U j are T 2 -dual, that is remain unitary under partial transpose of the q−dimensional subsystem then D q 1 ,··· ,q K is T-dual. Of course this is recursive process and the difficulty is now transferred to finding such special unitaries in lower dimensions.
However, we note that even if one such set of T-dual matrices are known then local unitary transformations on the H m j ⊗ H q provide a parameterization of T-dual q 2 × q 2 matrices of differing nonlocal content in the H q ⊗ H q space. The simplest T-dual matrix is the identity and hence we may take are arbitrary m j and q dimensional unitary matrices. This provides one way of constructing the blocks U i of D q 1 ,··· ,q K , such that it is T-dual. Once these are constructed U = SD q 1 ,··· ,q K and U = D q 1 ,··· ,q K S provide ensembles of dual-unitaries generalizing the easy uniform case. As K = 1 (if K = 1, m 1 = q and such operators are locally equivalent to swap S), they result in nontrivial dual-unitary families.
The entangling power of dual-unitaries depends only on E(US) = E(SU) and these in turn follow from the singular values of the partial transpose such as U T 2 . From Eq. (26e) we get that U T 2 = (SD q 1 ,··· ,q K ) T 2 = S(D q 1 ,··· ,q K ) R 1 . The action of the realignment R 1 on a q 2 × q 2 matrix again proceeds by uniformly dividing it into q 2 q × q subblocks, column vectorize the sub-blocks starting from the top left corner, and going along the blocks column-wise. Thus (D q 1 ,··· ,q K ) R 1 contains as many columns as there are q × q blocks in it and we get that as each q j −dimensional sub-block contains m 2 j q− dimensional ones. Along with ∑ K j=1 m j = q, this renders (D q 1 ,··· ,q K ) R 1 rank-deficient (rank < q 2 ). In the extreme case when K = q blocks are present, m j = 1 uniformly and the rank is lowest at q, which includes the case of D being a diagonal unitary.
This implies, as we show in Appendix D, that for these dual-unitary ensembles constructed from block-diagonal T-dual matrices, where e * p is the bound from Eq. (48). The smallest range of the entangling power occurs for the uniform case when K = q, Hence, for q = 2 the maximum value e p (U) = 2/3 is already reached by the block-diagonal cases above and can be achieved even with diagonal 2-qubit gates. The largest range occurs for the case K = 2, which interestingly coincides with the bound e * p derived from the norm of the channel as a sufficiency condition for mixing. Thus none of the dual-unitary constructions based on diagonal blocks, satisfy the sufficiency condition, Eq. (48) for mixing: there are always members that do not mix.
Whenever (D q 1 ,··· ,q K ) R 1 is rank-deficient the entangling power cannot be maximized as K > 1. As all blockdiagonal based constructions are rank-deficient these will not reach the maximum entangling power, with the exception of q = 2. Quite remarkably, the maximum depends only the number of blocks and not on the block sizes m j . For the case of qutrits, the uniform block-diagonal based ensemble can only reach a maximum entangling power of 3/4. If we take the non-uniform case when K = 2, m 1 = 2 and m 2 = 1, the maximum according to Eq. (81) is 7/8 = 0.875. Interestingly, in this case, the bound is not tight, and is not achieved. For ququads with q = 4, if K = 2, m 1 = m 2 = 2 the bound on e p (U) is 14/15 and is achieved. The distinction is that there are 8 dimensional 2-unitary matrices in the bipartite qubit-ququad split: H 2 ⊗ H 4 [86], however there are no 6 dimensional 2-unitary matrices with the qubit-qutrit bipartition [87].

Channels from block-diagonal based constructions
We will study the CPTP maps M + [U] when U = D q 1 ,··· ,q K S and U = SD q 1 ,··· ,q K and denote them simply as . Recall that D q 1 ,··· ,q K consists of K unitary blocks, denoted U k , on the diagonal of sizes q j = m j q, with ∑ j m j = q.
From Eq. (21) and using Eq. (26f) we get It is not hard to see that the CPTP map There are in all K 2 such blocks with sizes in the set {m k m l |1 ≤ k, l ≤ K}, that is In the simplest uniform case, there are K = q blocks, all the q j = q and m j = 1, M + is diagonal. In this case U R 2 k,l are row vectors of length q 2 and hence the product U R 2 k U R 2 † l is simply an inner product, which in terms of the original matrices is the Hilbert-Schmidt inner product. Hence the (unordered) spectrum of M + [DS] is Hence there are q eigenvalues λ kk = 1, and local unitary operations are needed to render the circuit mixing. Also interestingly, if the blocks contain orthonormal unitary matrices, then all the other eigenvalues are 0, equivalently λ kl = δ kl , and represent the most mixing that can be obtained from this type of dual-unitaries. These continue to hold for the special case when U k are themselves diagonal matrices, and was considered recently in [56]. In the uniform blocks case, the entangling power cannot exceed q/(q + 1) as shown in Eq. (82). The maximum possible entangling power is when K = 2, when e p (U) ≤ e * p . It was mentioned that there are gates U such that e p (U) = e * p = (q 2 − 2)/(q 2 − 1) and λ 1 = 1. It is shown in App. (D) that if U k are 2-unitaries in mixed dimensions m k × q, that is they satisfy the conditions of Eq. (D1), then the entangling power is e * p . However, these conditions imply that there are K blocks in M + [DS] which are where |φ + m k is the maximally entangled state ∑ m k i=1 |ii / √ m k . Thus in general there are K eigenvalues that are 1 and ∑ k (m k − 1) = q − K zero eigenvalues. Hence |λ 1 | = 1, also in the particular case of K = 2. As we have previously noted, there exists no matrices that satisfy these conditions in q = 3 and K = 2, as there exists no 2 × 3 2-unitary matrix. However there do exist such block matrices for q = 4 when we could have two 8 × 8 blocks, each satisfying the conditions of Eq. (D1) as there exists 2-unitary matrices in H 2 ⊗ H 4 . This follows from the existence of a 4-party AME state of two qubits and 2 ququads, which is a corollary of the fact that 6 qubit AME states exist [87]. We show below that the quantum cat map is another example wherein e p (U) = e * p but can be non-mixing. Turning to M + [SD], which is also M − [DS], we get where an identity in Eq. (26e) is used. The difference with the previous case that involved R 2 , Eq. (83), is that here there no "interference" terms between the blocks as a simple schematic of multiplying such realigned matrices will convince a reader. This results in This simplifies for K = q, when U j are uniformly q × q unitary matrices as Thus unlike M + [DS] this is not diagonal and generically there are no other eigenvalue that is 1 apart from the trivial one. Thus the same circuit could be chiral in the nature of their mixing. Also, note that in this case, using Eq.
Thus the channels of locally transformed operators are equivalent to the channels of other members of dual-unitary matrices of the same form.
To construct dual-unitary operator U with an entangling power in excess of the bound e * p , which assures mixing, and to construct Bernoulli circuits from 2-unitaries, we need to look beyond the block-diagonal constructions. A few possibilities are outlined below.

B. Dual-CUE
As shown above, the entangling power e p (U) of block diagonal dual-unitaries is bounded by the number of blocks, and cannot exceed the sufficiency bound in Eq. (48) for a mixing circuit. They certainly cannot result in Bernoulli circuits. Permutations and decorations thereof may exceed it and even reach the maximum possible, and for odd dimensions there exist cat maps that can do this. However, as these are special constructions, an ensemble of dual-unitary matrices with large entangling power still seems inaccessible. However, in Ref. [46], we constructed a nonlinear dynamical map on the space of unitary operators, which acting iteratively, limit to dual operators and which have large entangling power.
We recall this procedure here for completeness. Define the map M R : U (q 2 ) → U (q 2 ), which is schematically represented as follows: R is a linear map, R(U) = U R 2 , whereas the map P D is nonlinear, such that P D(X) = V, where V is the closest unitary operator to X (under any unitarily invariant norm), which is given by the polar decomposition X = V √ XX † . It was observed in [46] that, starting from almost any unitary operator U 0 , an iterative application of M R , U n = M n R (U 0 ), produces unitary U n with an operator entanglement E(U) monotonically approaching the maximum value 1 − 1/q 2 as n → ∞, so that M ∞ R (U 0 ) is dual-unitary . Furthermore, defining the Tsallis entropies, where p j = γ j /q 2 (see Eq. (13)), the α = 2 case corresponds to the operator entanglement E(U) defined in Eq. (15). Under the M R map, it was proved in [46] that showing a nondecreasing trend towards an entropy increase and making plausible that the limiting process results in a maximum entropy state. From a dynamical systems perspective, the dual-unitaries are fixed points of the M 2 R map. To see this, let U be a dual-unitary operator to begin with, then action of M 2 R is as follows: Since U is dual-unitary, so that U R 2 is unitary and it remains invariant under the P D map. As shown in [46] there is overwhelming numerical evidence that almost any unitary is in the basin of attraction of the dual-unitary fixed point set, and is reached under the M R map. Thus if the starting unitaries are sampled from the Haar measure, the resultant ensemble of dual-unitary matrices, M ∞ R (CUE) was dubbed "dual-CUE". The distribution of entangling power e p (U) for dual-CUE is given in Fig. (8) in which the average entangling power of CUE,ē p = (q 2 − 1)/(q 2 + 1) [49,52] and e * p are both indicated for comparison. It can be seen that for q = 3 the M R map also produces a fraction of 2-unitaries whose e p (U) = 1. If e p (U) > e * p , the dual-unitary circuit is guaranteed to be mixing. For comparison, the average entangling power of CUEē p is also indicated.
Additionally, it is possible to define other maps wherein the partial transpose is added as a step in the protocol so that the resultant matrices are 2-unitary. It was shown in [46] that especially for q = 3 and q = 4, this was successful in producing ensembles of 2-unitary matrices, with e p (U) = 1, from which Bernoulli circuits may be constructed. Thus, although these maps result typically in numerically close to dual and 2-unitary matrices, they are for practical purposes as good as analytical ones and crucially enable us to explore gates with large entangling powers. In Appendix (F), we provide examples of a dual-unitary operator with entangling power e p (U) = 8/9 > e * p = 7/8 and 2-unitary operator obtained via this numerical procedure for q = 3. Note that this procedure is capable of identifying exact operators if the elements are simple enough to be recognized, as in these examples. FIG. 9. The maximal mixing rate ν + and the average mixing rate µ + for q = 3 and K = 1, 3 cases of dual-unitary operators of the form D 3 S, and D 1 S. In the former case the 3 blocks are random unitaries drawn from the CUE, while in the latter the diagonal entries are random phases drawn uniformly from the unit circle. The maximum possible entangling power 3/4 is indicated by dashed vertical line. The inset shows the case of dual-CUE, again with q = 3, which can explore a larger range of entangling power.
The maximum and the average mixing rates ν + and µ + of Eqs. (40) and (41) are shown for q = 3 in Fig. (9), using three types of dual-unitary constructions. Two of these are based on block diagonal constructions, while the D 1 S uses diagonal unitary matrices, the D 3 S uses random 3 × 3 CUE blocks on the diagonals. While both of these have the same range of possible entangling power, namely [0, q/(q + 1)], it is seen that D 3 S typically has larger entangling power. The third ensemble used is the dual-CUE, which can have as large as e p (U) = 1, but q = 3 is a peculiar case, where there is a gap at e p (U) = 1. Nevertheless we can see that this realizes gates of larger entangling power than the diagonal blocks based ensembles. For larger entangling power ν + does not appear to be a function of the entangling power alone and these numerical results indicate that it is a more complex quantity. For smaller values of e p (U) it has been noted in Fig. (4) that it is simple function of the entangling power. On the other hand the single-particle averaged µ + seems to be determined by the entangling power and is close to the estimation from |λ 1 |.

C. From permutations
Permutation matrices form a source of unitary gates with a wide range of entangling powers [73]. They have been used to construct gates with e p (U) = 1 for all dimensions except q = 2 and q = 6. The swap and identity operators are two permutations with vanishing entangling power but have different nonlocal contents. The characterization of permutation matrices that are dual and T-dual does not seem to have been systematically exposed, although the ideas are present in [73].
If [q] = {1, · · · , q}, a permutation of q 2 elements [q] × [q], written as P(i, j) = (k ij , l ij ), is represented by the pair of matrices K = (k ij ) and L = (l ij ). Corresponding permutation operator P acts on the product basis as P(|i |j ) = k ij l ij , and represented as For example, for q = 3, if P permutes the elements of As P R 2 = ∑ i,j k ij |i l ij j|, if k ij k in = δ jn , then P R 2 † P R 2 = I. It is also easy to see that the converse is true. Similarly considering conditions under which P R 2 P R 2 † = I, one concludes that the permutation operator P is dual-unitary, that is P R 2 is unitary, if and only if k ij k in = δ jn and l ij l mj = δ im . This imposes a constraint on the K and L matrices. The elements along the each row of the K matrix and the elements along the column of the L matrix must be distinct. Thus the example above does represent a dual-unitary permutation.
Similarly, the permutation operator P is T-dual if and only if k ij k mj = δ im and l ij l in = δ jn . The elements in each column of the K matrix and the elements in each row of the L matrix must be distinct. The conditions for dual and T-dual together defines the conditions for P to be 2-unitary. This requires that both K and L be Latin squares, with unique entries along any row or column. In addition for the combined KL to be a permutation, we must require that these Latin squares be orthogonal: overlaying the squares result in no repetition of the order pairs in  While permutations provide a discrete set of dual, T-dual, and 2-unitaries, one may use a kind of gauge-freedom they enjoy to define continuous families that are not local unitarily connected. Multiplying by diagonal unitaries on either sides does not change the dual nature of the permutations, we define generically locally inequivalent ensembles of unitaries that are "enphased" permutations as where −π ≤ θ ij < π. It is easy see that the same condition for duality, T-duality and 2-unitary given above for P continues to hold. However, if we take a single dual permutation, the members of its enphased family will have different entangling powers. The exceptions are 2-unitaries, which remain 2-unitary under any enphasing. This follows as 2-unitaries are both dual and T-dual, and enphasing does not change either of these properties.

D. From cat maps
It seems natural to look at quantized classically chaotic systems for dual and 2-unitary operators. Quantum maps have long been studied whose classical counterparts are dynamical systems with a high degree of chaoticity [11], in particular the baker [10,88,89] and the cat maps [8,90]. The cat maps are linear automorphisms of the 2-torus that are hyperbolic, Anasov or C-systems. Their quantization was first studied by Berry and Hannay [8] and provides a rich source of models for which semiclassical methods are exact. But they also have peculiar properties such as periodicity in time, whose origins lie in number-theoretic aspects of the classical cat itself [8,14,91]. As our concentration in this work is not primarily with these interesting models, but rather use them as examples for the class of unitary operators to use in quantum circuits, we restrict ourselves to a minimal description. Coupled cat maps from the point of view of duality have been studied in [92] and coupled quantum cat maps have also been studied in [93].
The Arnold cat map is area-preserving and defined on the torus [x, p] ∈ [0, 1) 2 by As we are interested in bipartite unitary operators, two coupled cat maps are considered in which case the automorphism is of the 4-dimensional torus and the map on it is symplectic. There have been previous studies involving such coupled maps, such as in [94]. However, we define the quantum propagator via the classical generating function, a route taken in the absence of a Hamiltonian for example in the construction of the quantum bakers map [10]. The linearity of the map implies a quadratic generating function. Considering the generating function on R 4 a symplectic transformation is obtained from p i = ∂F 1 /∂x i and p i = −∂F 1 /∂x i [95]. This linear map is a differentiable one on the compactified unit torus T 4 if a, b, c are integers and c|a, c|b. We take a = −1, b = −2 and c = 1 to get the 4-dimensional symplectic map which we recognize as Arnold cat maps on the (x 1 , p 2 ) and (x 2 , p 1 ) subspaces. The corresponding quantum propagator is given by the prescription of Jordan and Dirac x 1 x 2 |U|x 1 x 2 = e iF 1 /h [96,97]. On the unit torus with q states in each degree of freedom h = 1/q and the position and momentum take discrete values k/q with 0 ≤ k ≤ q − 1. This leads to the quantum cat map corresponding to the classical map in Eq. (93) acting on the bipartite system in H q ⊗ H q : where α, β, j, k ∈ {0, 1, · · · , q − 1} label position eigenkets.
The operator U C is dual-unitary for all q and it is 2-unitary operator for all odd local dimension q (See appendix G). To our knowledge, this is the first time a simple non-permutation class of unitaries have been shown to be 2-unitary. Note that this is a complex Hadamard matrix. Examples of complex Hadamard matrices are known to be 2-unitary for q = 3 [47]. Besides, being related to dynamical systems, and their similarity to Fourier transforms, they may be implementable in efficient quantum circuits or other physical setups [98,99]. Thus while e p (U C ) = 1 for q odd, it is amusing that for q even, is at the border of the sufficiency condition for mixing. While the unitary U C does not provide a family or ensemble of dual operators, a one parameter family may be found by setting a = 0 and c = 1 for arbitrary b in the generating function Eq. (92). The corresponding quantum maps are This forms a dual-unitary family for all b ∈ R. However its entangling power is no larger than q/(q + 1) (reached when b = 1), which one recognizes from Eq. (82) as that of a uniform block-diagonal case. Indeed, it is not hard to see that where F q is the q− dimensional discrete Fourier transform and kα|D q 2 (b)|jβ = e 2πib(jβ)/q δ kj δ αβ is a diagonal unitary. Hence this family of cats is locally equivalent to a swap-diagonal family such as we considered earlier, irrespective of the parity of q. We note that U C in Eq. (94) also has a decomposition of the form SD q 2 (F q ⊗ F q )D q 2 where D are diagonal matrices, but the fact that the local product is sandwiched between two nonlocal ones makes the structure special and require that the parameters a, b be fine-tuned. Complex Hadmard matrices have been used in similar constructions by [55] and it will be interesting to compare them and find their entangling power.

Cat map channels
The channel M + for the quantum cat map in Eq. (94) is given by For q odd, δ 2(α−k) mod q,0 = δ αk and M + [U C ] = |Φ + Φ + | as behoves a 2-unitary cat, and the circuit is Bernoulli. For q even, the cat map U C is dual but not 2-unitary. The M + in Eq. (97) takes a simple form as where |Φ + is the maximally entangled state in Eq. (12) and are two other maximally entangled states. The operator |Φ + Φ + | as always corresponds to the trivial eigenvalue, and Ψ Ψ is non-diagonalizable, as |Ψ and Ψ are orthogonal. However, it is easy to see that, and hence the cat map for q is such that correlations vanish after the second time step. We also note that M − [U C ] = M + [U C ] as U C S = SU C , the cat is parity or exchange symmetric. Under local single-particle unitaries though the cat map channel changes as Eq. (39). It is sufficient to consider one local unitary to get M + [U C ] = (u ⊗ u * )M + [U C ], where u ∈ U (q). Then in general M 2 Therefore any state that is orthogonal to both |Φ + and |Ψ is a (right) eigenstate of M + [U C ] with eigenvalue 0, and hence the channel is rank-2. Ignoring the trivial eigenvalue 1, we get If u is the discrete Fourier transform on q sites, with 0 ≤ ϕ 1,2 < 1 being phases, then λ 1 (U C ) = cos(π ϕ 2 ). Thus the dual circuit can be rendered non-ergodic if ϕ 2 = 0 or ergodic if ϕ 2 = 0.This shows that the cat map for q even, with suitable local unitaries, is an example of a case when e p (U) = e * p , but the circuit is non-mixing. A generic choice of u will however render the circuit mixing. Indeed, if u is chosen from the Haar measure, it is a straightforward application of Eq. (B2) to get the local-unitary averaged This quantum cat map is therefore a solvable case wherein the general relation between the largest nontrivial local single-particle averaged eigenvalue and the entangling power, Eq. (62) takes this exact and simple form.

V. SUMMARY AND DISCUSSIONS
This paper has considered the recently introduced dual-unitary circuits that are models of many-body quantum chaotic systems [26]. We showed that the partial-transpose and realignment operations central to the problem of separability of quantum states [59][60][61], play a central role in these circuits as well. Defining T-dual gates then, we showed how the CPTP maps that govern the decay of correlations along the lightcone is related to the partialtranspose operation. This leads naturally to Bernoulli circuits that are at the apex of a putative quantum ergodic hierarchy. These are made of 2-unitaries or perfect tensors which remain unitary under both realignment and partial transpose. Determining necessary and sufficient conditions for the circuit to be Bernoulli, we showed how two-site observables retain the property of having Bernoulli correlations.
A simple observation that the norm of the CPTP maps determining the decay of correlations is directly related to the entangling power of the two-particle unitary building the circuit, enables us to derive a sufficiency condition for mixing. There are unitaries U such that e p (U) = e * p and the circuit is not mixing, and not even ergodic. For the case of qubits it is also true that if e p (U) ≤ e * p , there exists local unitary transformations such that the circuit can be made non-ergodic. This follows from the results in the Supplementary material of [26]. It is unclear if this generalizes to q > 2, but we believe that this true. However, increasing entangling power does ensure that an increasing number of modes are rendered mixing, whatever the local unitaries maybe.
While the entangling power are local (single-site or single-particle) unitary invariants represent local fields, the mixing properties of the circuit are not. We show however that local unitary averaged mixing rates depend simply and universally on the entangling power. That these averages are in fact typical values is shown by numerical results. This paper has also systematically enumerated several methods of constructing dual-unitary operators. Importantly, we also calculate the entangling power in these constructions as indeed these play a central part in the ergodic properties of the quantum circuit.
There are very many natural directions that this work throws up. The use of T-duality by itself or in conjunction with realignment, such that the unitary remains unitary under serial applications of these rearangements, in quantum circuits maybe interesting. This is related to using operators other than dual-unitary ones. The rudimentary phase-diagram where the ergodicity properties depend on the entangling power alone, will likely expand to include at least two LUI quantities such as the operator entanglements E(U) and E(US) or the entangling power and the so-called gate typicality [52]. Another direction is the use of dual-unitary operators with particles with different local dimensions, such as a mixture of qubits and qutrits [86]. A closely related possibility is to have multi-particle interactions, rather than only nearest neighbor ones. For example, a tripartite 3-body interaction can enable Bernoulli or perfect tensor circuits even with qubits as there are AME states of 6 qubits [86,87].
Other directions that seem intriguing are comparisons of the ergodic properties as discussed via correlation functions and those from standard quantum chaos literature, such as spectral statistics. Finite systems, with recurrences are natural in this context, and the nature of eigenstates of the Floquet operator and the development of complexity in nonequilibrium states are of interest. Also, we have used a cat operator that is capable of being a 2-unitary. It will be interesting to study other physical systems whose Floquet operators are indeed 2-unitary or in general dual operators. It is unclear to us what quantity would play the role of entangling power in organizing the ergodic hierarchy of possible classical limits of such dual-unitary circuits.
The vanishing correlation function along the lightcone, x 1 , x 2 = t − 1 2 , due to the application of T-duality condition (28), and C ijkl (x 1 , x 2 , t) = 1 q 4t a i a j a k a l . By substituting Eq. (B2) into Eq. (B1), we get the following identity: (B3)

Appendix C: Unistochastic matrices and CPTP map M + [U ]
The nontrivial largest and smallest magnitude eigenvalues λ 1 and λ q 2 −1 of M + [U ] is plotted in the Fig. (2), where U = (u 1 ⊗ u 2 )U(v 1 ⊗ v 2 ). The following dual-unitary U with e p (U) = 3/4 is used to generate λ 1 and λ q 2 −1 of M + [U ] as shown in the lower right figure in Fig. (2). This unitary is obtained by iterative application of the nonlinear map M R (U) to a permutation operator. Note the elegant form the dual-unitary operator obtained from a numerical algorithm, and is given as : The structure of the real and imaginary parts of the nontrivial largest eigenvalue λ 1 of M + [U i ] in Fig. (5) deserves an explanatory note. First note that the dual-unitary operator U 1 = D 3 S and U 2 = D 2 S are having same entangling power e p (U i ) = 3/4. The local unitary inequivalence between U 1 and U 2 is evident from the different λ 1 distribution seen in the Fig. (5). The rank of M + [U 1 ] is 3 and M + [U 2 ] is full. The block diagonal unitary operators used to construct these dual-unitary operators are, where a ≈ 0.315167.
The geometric structure in the Fig. (5) corresponding to the statistics of λ 1 of M + [U 1 ] is called deltoid, a 3hypocycloid. Here we briefly explain the reason behind this structure as these structures were also reported in Ref. ([80]) in the context of random unistochastic matrices. By using Eq. (21), the operator M + [U 1 ]is given by its matrix elements as Here M + [U 1 ] |jβ = δ jβ |jβ is used. The subspace contributing the three eigenvalues are spanned by the vectors {|ii }. The corresponding matrix elements are ii|M + [U 1 ]|jj = i|u|j i|u * |j = u ij 2 . Thus the eigenvalue of M + [U 1 ] is same as that of the unistochastic matrix B = u ij 2 of dimension 3 × 3 generated from the local unitary operators u = u ij . In ref. ([80]), the authors constructed an ensemble of unistochastic matrices of dimension 3 × 3 starting from the nontrivial permutation operators and showed that the spectrum is bounded by a 3-hypocycloid. They conjecture that the spectrum of any random unistochastic matrices of dimension 3 × 3 forms a 3-hypocycloid.
They also construct random unistochastic matrices for the dimension N × N and show that the spectrum forms N-hypocycloid for such constructions. But in general the spectrum may deviate for N > 3.
The dual-unitary U 1 = D 4 S and U 2 with e p (U i ) = 4/5 are given below. Here U 2 is obtained from the iterative application of M R map to a permutation operator. Again note the elegant form of the dual-unitary obtained from the M R map. U 1 and U 2 are not local unitarily equivalent and the rank of M + [U 1 ] is 4 and M + [U 2 ] is of full rank. The distribution of real and imaginary values of nontrivial largest eigenvalue λ 1 of M + [U i ] is plotted in Fig.(11). There is a small fraction of eigenvalues that deviates from the 4-hypocycloid. The single-particle averages E u i ,v i (|λ 1 |) = 0.46170 for U 1 and E u i ,v i (|λ 1 |) = 0.4940 for U 2 .
FIG. 11. The real and imaginary parts of the nontrivial largest eigenvalue λ 1 of M + [U i ] for two dual-unitary U i with e p (U i ) = 4/5. A finite fraction of spectrum of M + [U i ] is deviated from the 4-hypocycloid. The matrix U 2 converges from iterating the map M for a large number times n to a particular permutation matrix P.

Appendix D: Maximum entangling power of block diagonal unitaries
The bound in Eq. (81) is attained when the blocks U j of Eq. (79) are not only such that their partial transposes are unitary, but that their realignments are also unitary; that is they are 2-unitaries in the bipartite H m j ⊗ H q space. As m j < q these spaces are not symmetric under exchange and the conditions of 2-unitarity are modified accordingly [86]. They are while U R 1 j U R 1 † j is not proportional to identity. The operator Schmidt decomposition of U j is in general where {a j,k } and {b j,k } are orthonormal operator bases in H m j and H q respectively, tr a j,k a † j,l = δ lk and tr b j,k b † j,l = δ lk . The Schmidt values satisfy ∑ m 2 j k=1λ j,k = m j q, and U j is maximally entangled ifλ j,k = q/m j for all k, from which the second relation in Eq. (D1) follows. If all the U j are maximally entangled, as (A ⊗ B) R 2 = |A B * |. It follows that (D q 1 ,··· ,q K ) R 2 (D q 1 ,··· ,q K ) R 2 † = K j=1 q m j m 2 j ∑ k=1 |a j,k a j,k |, as |a j,k a j ,l | = 0 for j = j , their supports not overlapping. The singular values of (D q 1 ,··· ,q K ) R 2 are therefore q/m j with a m 2 j fold degeneracy. Using Eq. (15) we get Hence as E(SU) = E(D q 1 ,··· ,q K ), the bound in Eq. (81) arises. Despite the restriction, explicit forms of the roots λ 1,2,3 , and their ordered versions, are cumbersome. They satisfy the following relations, λ 1 λ 2 λ 3 = sin 2 2J, λ 1 + λ 2 + λ 3 = cos φ sin 2J, and are independent of ψ. As the roots depend on cos φ it is sufficient to restrict attention to the interval φ ∈ [0, π]. The absolute values of λ i are shown in Fig. (13) as a function of φ for a representative value of J = π/16. The phenomenology of the roots shown is independent of the value of J used. For φ = 0, it is easy to see that λ 1,2 = ±i √ sin 2J, while the real root is λ 3 = sin 2J. At φ = π, all three roots are real: λ 3 = √ sin 2J while λ 1 = − √ sin 2J and λ 2 = − sin 2J. The root with the largest magnitude, |λ 1 |, decreases from |λ 1 | = |λ 2 | = √ sin 2J at φ = 0 and then increases again to the same value at φ = π. We are interested in its smallest value, and it turns out that this is obtained when φ = π/2, independent of the value of J used. as a function of (θ, φ) for a fixed value of J. Pairs of (θ, φ) angles for which the bound is reached are seen as lighter in the plot Fig. 15 where the function takes large values (diverges in principle). For J = 0, bound is reached at θ = π/2 for all values of φ. The dependence of the minimum of |λ 1 (J, θ, φ)|, obtained by searching over different values of (θ, φ) pairs, on J or equivalently, e p (U) is shown quantitatively for different values of J in Fig. (16).
The T-dual property is restricted to only odd dimensions. This can be shown by considering the matrix elements For q odd, the second Kronecker delta product is nonzero iff k = j and hence there are only diagonal entries (= 1) in this case and T-duality follows. If q is even, there are nonzero elements (= (−1) β ) when α = β and (k − j) = ±q/2. As the unitary cat map U in Eq. (94) is dual for all dimensions q, it is 2-unitary for q odd.