Chaos from Equivariant Fields on Fuzzy $S^4$

We examine the $5d$ Yang-Mills matrix model in $0+1$-dimensions with $U(4N)$ gauge symmetry and a mass deformation term. We determine the explicit $SU(4)\approx SO(6)$ equivariant parametrizations of the gauge field and the fluctuations about the classical four concentric fuzzy four sphere configuration and obtain the low energy reduced actions(LEAs) by tracing over the $S_F^4$s for the first five lowest matrix levels. The LEA's so obtained have potentials bounded from below indicating that the equivariant fluctuations about the $S_F^4$ do not lead to any instabilities. These reduced systems exhibit chaotic dynamics, which we reveal by computing their Lyapunov exponents.Using our numerical results, we explore various aspects of chaotic dynamics emerging from the LEAs. In particular, we model how the largest Lyapunov exponents change as a function of the energy. We also show that, in the Euclidean signature, the LEAs support the usual kink type soliton solutions, i.e. instantons in $1+0$-dimensions, which may be seen as the imprints of the topological fluxes penetrating the concentric $S_F^4$s due to the equivariance conditions, and preventing them to shrink to zero radius.


Introduction
Matrix models associated to M-theory an string theories [1][2][3] have been under investigation from various perspectives ever since their discovery over twenty years ago. The broad span of interest on the subject is reflected in the literature (for a recent review, see [4] and references therein.) Among these, the BFSS model [1] is a supersymmetric U (N ) gauge theory consisting of nine N × N matrices in its bosonic part, whose entries depend on time only and it also goes by the common name of matrix quantum mechanics in the literature. It is associated to type II-A string theory [4][5][6] and appears as the DLCQ (discrete light-cone quantization) of M-theory on flat backgrounds [4][5][6]. The massive deformation of the BFSS theory, preserving the supersymmetry, is known as the BMN model and describes the DLCQ of M -theory on pp-wave backgrounds [2,3]. These matrix models describe systems of N coincident D0-branes, respectively in flat and spherical backgrounds. The latter is due to the fact that fuzzy 2-spheres appear as vacuum configurations in the BMN model. At large N and strong coupling/low temperature limit, D0-branes form a black brane, i.e. a string theoretic black hole [4,5]; the structure of this gravity dual is discussed in varying amount of detail in several references [4,5]. In [7], a model on how the fuzzy spheres in the BMN model collapse to form a black hole is discussed.
The perspectives gained from the matrix models have recently started to motivate numerous investigations oriented to acquire new information on the properties of black holes, such as their thermalization, evaporation processes as well as their microscopic constituents via the study of BFSS and BMN models at large N , using both analytical insights and ever increasingly numerical techniques [8][9][10][11][12][13][14][15][16][17]. In some of these works, for the BMN and BFSS models in the large N and high temperature limit numerical evidence for the fast thermalization is obtained [11][12][13][14]. From a more general perspective, the latter is an example of the fast scrambling conjecture, which has been proposed by Sekino and Susskind [18] and which may be stated as the fact that in black holes, chaotic dynamics set in faster than in any other physical system and the rate at which this occurs is logarithmic in the number degrees of freedom of the black hole, that is, it is proportional to the logarithm of its Bekenstein-Hawking entropy. Even though the aforementioned limit is distinct from the one in which the gravity dual is obtained, it is the natural limit in which the classical mechanics provides a good approximation to the quantum mechanical matrix model 1 . This limit is free from fermions as the latter contributes to the dynamics of the bosonic matrices only at low temperatures. It has been also noted that [14], since, numerical studies performed so far do not show a phase transition occurring between the low and the high temperature limits of these matrix models [8][9][10]19,20], it is reasonable to expect that features like fast scrambling of blackholes in the gravity dual could survive at the high temperature limit too. In fact, chaotic dynamics in the BFSS model is studied in [14] in this classical limit by calculating the Lyapunov spectrum, where it was also demonstrated that a classical analogue of fast scrambling is valid for this system as the scrambling time is found to be proportional to log N 2 . In [21] simulations of the BFSS model is performed at intermediate temperatures to numerically study the black hole horizon.
Even the matrix models at small values of N appear to be highly non-trivial manybody system whose complete solution evades us to date. Recently, there also has been some interest in examining the chaotic dynamics emerging from such models [15,22] (See, also [23], for earlier attempts). These studies also provided some qualitative implications on the black hole phases in such models; for instance, in [15] it has been argued that the edge of the chaotic region found in a SU (2) YM matrix model with only two matrices, i.e. the smallest matrix model with non-trivial dynamics, corresponds to the end of the black hole phase. In [24] a novel approach has been developed to estimate the ground state energy of this smallest matrix model. Authors of [22], considered simple ansatzes for the BMN model at N = 2, 3 satisfying the Gauss law constraint to probe the chaotic dynamics.
Fuzzy two sphere and its direct sums are not the only compact spherical geometries appearing in M-theory. In fact, it has been known for a quite long time that fuzzy four spheres make their appearance in matrix models as longitudinal five branes [25]. For the purposes of this paper however, presence of fuzzy four sphere solutions in Yang Mills 5matrix model with massive deformation term plays the central role. The latter has also been known in the literature for quite a while [26] and may be obtained as the reduction of the YM theory in 5 + 1-dimensions to 0 + 1-dimensions keeping only the time dependence of the matrix elements. It can also be conceived as a deformation of a subsector of the bosonic part of the BFSS model, as we will explain in more detail in the next section. Contrary to the fuzzy two sphere solutions in the BMN model, fuzzy four spheres in the mass deformed YM-matrix model are classical solutions for negative mass squared (µ 2 = −8), which may be an indication of tachyonic instabilities. Nevertheless, it was recently shown by Steinacker [27] that in pure YM 5-matrix model, one-loop quantum corrections stabilizes the radius of the fuzzy four sphere and prevents its collapse by shrinking to zero radius. In this paper, we will mainly focus on a mass deformed U (4N ) YM 5-matrix model and consider the exact parametrization of SU (4) equivariant fluctuation modes about the four concentric fuzzy four sphere configurations. Using the equivariant parameterization of the gauge field and the fluctuations, we perform the traces over the fuzzy four spheres at first five lowest lying levels and obtain the corresponding low energy reduced actions (LEAs). We demonstrate that the potentials of all of these reduced effective actions are bounded from below, from which we infer that the negativity of µ 2 does not actually cause any instabilities under equivariant fluctuations. As we will briefly discuss in section 3, this feature of our treatment may also be viewed as a consequence of the fact that the equivariant parametrization of the fluctuations introduces topological fluxes through the fuzzy four sphere, preventing it to shrink to zero radius.
Equivariant parametrizations breaks the U (4) symmetry of the concentric S 4 F configuration down to U (1) × U (1) × U (1) and this is further reduced to only U (1) × U (1) in LEAs as one of the gauge fields completely decouple after tracing over S 4 F . The gauge fields in the reduced actions are not dynamical, and their equations of motion lead to the constraints, which may, in fact, be seen as the residue of the Gauss law constraint on the matrix model enforcing the physical states to be the singlets of the gauge symmetry group. In the LEAs the latter condition simply translates to the requirement that the two complex fields appearing in the LEAs be real, that is, uncharged under the abelian gauge fields. This breaks the We utilize the LEAs to explore the chaotic structure emerging from the matrix model with the fuzzy four sphere background. For the reduced action obtained after tracing at the matrix levels 4N = 16, 40, 80, 120 and 224 corresponding to the first five S 4 F levels, n = 1, · · · , 5, we numerically solve the Hamilton's equations of motion and compute the Lyapunov spectrum at several different energies, revealing the chaotic dynamics. Except at n = 1 the phase space of the LEAs are all ten dimensional. At n = 1, however, out of five of the generalized coordinates and corresponding velocities in the LEA, three of them combine to appear only in a particular form thereby reducing the dimension of the phase space to six. Our findings are detailed by the plots given in Section 4, demonstrating the time development of the Lyapunov exponents, which all converge to fixed values well before the computation time used in the numerical calculations. Table 1 summarizes our numerical findings for the Lyapunov exponents.
The paper is organized as follows. In section 2, for completeness, we give a brief review of S 4 F and how it appears in Yang-Mills matrix models. In section 3, we first determine the exact parametrization of the gauge field and the fluctuations, which are restricted to transform as a scalar and vectors, respectively, under the combined adjoint action of SO(5) isometry of S 4 F and the SU (4) gauge symmetry generators in SO (5). In this section we also obtain the LEAs by tracing over the S 4 F at several different matrix levels, and elaborate on their structure. In section 4 we focus on the dynamical structure the LEAs. After discussing the implications of the Gauss law constraint, we present our results exhibiting the chaotic dynamics emerging from the LEAs. In section 5, we examine the properties of the LEAs in Euclidean signature, and make evident through a number of examples that, they possess Z 2 kink solutions, i.e. instantons in 1 + 0-dimensions. These may be seen as the imprints of the non-trivial topological fluxes piercing the S 4 F , which were mentioned above. Motivated by the recent work of Maldacena and Milehkin [28] on relaxing the Gauss law constraint in BFSS and BMN models (see [29] for supporting numerical work), in Section 6, we revisit the gauge symmetry of the LEAs and present a concise treatment on the consequences of not imposing the Gauss law constraint, which leads us to conclude the presence of massless excitations (Goldstone bosons) associated to these LEAs. We close the paper by summarizing our findings. Appendices collect reference formulas, intermediate steps of some of the analytic calculations, explicit form of the LEAs for n ≥ 2, the whole sets of the corresponding minima of the associated potentials, as well as a short general description of the calculation of Lyapunov exponents.

Basics
We launch the developments in this section by considering the Yang-Mills 5-matrix model in Minkowski signature and with U (4N ) 2 gauge symmetry, whose action may be given as [4,26,27] where X a (a : 1, . . . 5) are 4N × 4N Hermitian matrices transforming under the adjoint representation of U (4N ) as and T r stands for the normalized trace T r1 4N = 1. For future reference we write out the potential part of L Y M separately as There are two distinct deformations of S Y M preserving its U (4N ) gauge and the SO(5) global symmetries. One of these is obtained by adding a fifth order Chern-Simons term to S Y M (i.e. a Myers like term) which is given as [26] S CS = 1 g 2 dt T r while the other is a massive deformation term of the form It is useful to write out the potential terms for S 1 = S Y M + S M ass and S 2 = S Y M + S CS explicitly: (2.7) In this paper, our main interest is on the massive deformations and hence we will focus on the dynamics emerging from S 1 , but will, although very briefly, also consider the consequences of some of the developments presented in the paper for S 3 = S Y M + S M ass + S CS .
We may as well note that S 1 , S 2 and S 3 may be thought as deformations of a subsector of the bosonic part of the BFSS [1] matrix quantum mechanics. As is already well-known, BFSS model can be conceived to emerge from the dimensional reduction of the YM theory in 9 + 1 dimensions to 0 + 1-dimensions [4] with the SO(9, 1) of the former yielding to the global SO(9) of BFSS on the nine matrices X I (I : 1 . . . 9), while the latter may be further broken to SO(5) × SO(4), via the addition of S M ass and/or S CS terms. To be more precise, these deformations terms, involving only X a (a : 1, . . . 5), spontaneously breaks SO(9) down to SO(5) × SO(4) and naturally split the X I to an SO(5) vector X a and a SO(4) vector X α (α : 1 . . . 4). Then, S 1 , S 2 and S 3 emerges by focusing on the sector of X a 's only. V 2 is extremized by the matrices fulfilling This equation has two immediate solutions [26], one of which is the diagonal matrices while the other is given by a fuzzy four sphere S 4 F , where λ is forced to take on a fixed value λ = 2 n+2 , which depends on the matrix level of S 4 F . This latter fact requires attention, since it implies that the direct sum of fuzzy S 4 solve (2.8) if and only if they are fuzzy spheres of the same matrix level. For the configuration (2.9) the potential take the value zero as is readily observed from S CS , whereas a simple calculation shows that [26] it, in fact, has a lower (negative) value for the fuzzy S 4 solutions. Thus at the classical level, fuzzy S 4 appears to be a more stable solution than the diagonal commuting matrices. Further numerical studies have, however, revealed that the fuzzy S 4 is not a minimum of V 2 , but instead a saddle point [30]. 3 As for the massive deformation the potential part of S 1 is extremized by the matrices fulfilling Fuzzy four spheres S 4 F and their direct sums (even from different matrix levels) are solutions of this equation for µ 2 = −8. In a recent article Steinacker [27] showed that quantum corrections in the pure YM 5-matrix model stabilizes the radius of the fuzzy four sphere. We will see that superficial instability implied by negativity of µ 2 does not actually lead to a problem when we consider the equivariant fluctuations in S 1 about the S 4 F backgrounds. The reason for this is essentially that the potential of the emergent equivariantly reduced action is bounded from below at any finite matrix level. We may as well interpret this outcome as being due to the fact that the equivariant parametrization of the fluctuations introduces topological fluxes through the S 4 F stabilizing its radius. We will further elucidate on these points in the Sections 3 and 4. Non-trivial topological fluxes leaves its imprints as kink type solutions of the reduced action in Euclidean signature, as we will exhibit in Section 5.
For the action S 3 = S Y M + S CS + S M ass the potential is extremized by the matrices fulfilling [26] [ whose fuzzy four sphere solutions need to satisfy the relation λ = 8+µ 2 4(n+2) , which includes the previous two cases of S 1 and S 2 as µ 2 = −8 and µ 2 = 0, respectively. The equation of motion for S 3 isẌ in the gauge A = 0 and subject to the Gauss law constraint a [X a ,Ẋ a ] = 0. Equations of motion for S 1 and S 2 are readily inferred from (2.12) by the remark following (2.11). Being static, S 4 F configurations satisfying any one of the equations (2.8), (2.10), (2.11) also satisfy the corresponding equation of motions as well as the Gauss law constraint.
For reference it is useful to list the variation in the potentials under the variation of the matrices X a → X a + δX a . These are given as (2.13)

Models in the Euclidean Signature
Wick rotating to the Euclidean signature we make the changes, (2.14) In section 5, we will consider the kink-type solutions, that is to say, the instantons of the low energy effective actions in 1 + 0-dimensions, which we obtain from the equivariant reduction of S E 1 .

Brief Review of Fuzzy S 4
In this subsection we collect some of the main features of the fuzzy S 4 construction [4,25,26,[31][32][33]. To start with we note that S 4 is embedded in R 5 as Construction of fuzzy S 4 proceeds as follows. Let us denote by Γ a , (a : 1, · · · , 5) the Hermitian 4 × 4 gamma matrices associated to SO(5) fulfilling the defining anticommutation relations For concreteness we take them to be given in the form We introduce the n-fold tensor product acting on the n-fold completely symmetrized tensor product space which is the carrier space of the (0, n) IRR 4 of SO (5). Obviously the latter is equivalent to the completely symmetric tensor product Sym n (0, 1) of the fundamental 4-dimensional spinor representation of SO(5) acting on C 4 . Dimension of this representation and hence that of H n is given as X a are then N × N Hermitian matrices satisfying the relations which are the defining relations for the fuzzy four sphere, S 4 F . In fact (2.21a) may be seen as the SO(5) invariant condition giving the radius of S 4 F as r n := n(n + 4). This construction appears to be quite analogous to that of the fuzzy two sphere [34] with SO(3) of the latter replaced with SO (5), nevertheless only to the extent until one recognizes that the commutation relations of X a do not close but instead they are given by where G ab are the ten generators of SO(5) in its (0, n) IRR satisfying the commutation relations, G ab are anti-hermitian by the definition (2.22). Under the SO(5) transformations generated by G ab , X a transform as vectors (i.e. in the (1, 0) IRR of SO (5)) of SO (5) since We have that It is also useful note that, G ab and X a together generate the SO(6) in its (n, 0, 0) IRR. This structure may be compactly expressed by writing the 15 generators of SO(6) as with the commutation relations taking the usual form We now observe that the SO(5) invariant condition (2.21a) can be expressed as the difference of the quadratic Casimir operators: as (n, 0, 0) of SO(6) branches solely to (0, n) of SO(5) and hence they are of the same dimension N . The relation (2.21b) can also be expressed equivalently in the form Another noteworthy feature of S 4 F is that there is a S 2 F attaching to each point of S 4 F . In other words, there is a S 2 F bundle over S 4 F with the total space being CP 3 F . This fact is reflected in the fields defined on S 4 F carrying an intrinsic spin, whose rank can be up to n [26,32].
We may also record that the commutative limit is achieved by taking n → ∞. The intricate structure of S 4 F with S 2 F fibers leads to where x a , with x a x a = 1, are the coordinates of S 4 ⊂ R 5 and ω ab is antisymmetric in its indices, satisfy x a ω ab = 0, as seen by taking the commutative limit of the first equality in (2.29) and generate, in fact, the S 2 in the fibration S 2 → S 4 → CP 3 . Detailed discussion on this may be found in [35].
In the commutative limit, adjoint action of X a and G ab become the differential operators [26] ad X a → ∇ a := 2i (ω ab ∂ where the derivatives with respect to ω ab 's are given via We may note for future use that

Symmetry Constraints and Parametrization of the Fields
We now focus on the YM model with mass term, whose action is given by S 1 The potential V 1 has an extremum given by four concentric S 4 F , satisfying the equation (2.10). We observe that U (4N ) gauge symmetry of the action S 1 is broken down to U (N )×U (4), and the commutant of (3.2) is just U (4). Let us denote with F a the fluctuations about (2.10). We may write where the r.h.s is introduced as a self-evident short hand notation, which will be used in the ensuing developments. We are interested in finding the equivariant fluctuations about the configuration (3.2) 5 . To be more precise, we would like to concentrate on those F a , which are invariant under the SO(5) rotations of S 4 F up to SU (4) ⊂ U (4) gauge transformations. For this purpose we proceed as follows. We introduce the equivariant symmetry generators as are the generators of SO(5) in the 4-dimensional fundamental spinor representation (0, 1). They constitute the subset of generators Σ AB : in the fundamental (1, 0, 0) spinor representation of SO (6). Evidently, W ab satisfies the SO(5) commutation relations and its SO(5) representation content is given by the decomposition of the tensor product (0, n) ⊗ (0, 1) into a direct sum of IRRs as We digress a moment to note that, this structure can be lifted to SO(6) group by writing W AB = (W ab , W a6 ) with the representation content (n, 0, 0) ⊗ (1, 0, 0) ≡ (n + 1, 0, 0) ⊕ (n − 1, 1, 0) , (3.6) whose branching under the SO(5) simply yields (3.5). Same branching under SO(5) also holds for the complex conjugate representation (0, 0, n) ⊗ (0, 0, 1) ≡ (0, 0, n + 1) ⊕ (0, 1, n − 1). Let us examine the adjoint action of W ab . It consists of two terms, first of which generates the infinitesimal SO(5) rotations of the S 4 F , while the second term is responsible for generating the SO(6) gauge transformations in SO (5). Adjoint action of W ab carries the tensor product of the reducible representation (3.5) with itself and using general formulas on tensor product of SO(5) IRRs, it has the following decomposition in terms of SO(5) IRRs with the coefficients in bold indicating the multiplicities of the respective IRRs. For the case of n = 1, the decomposition takes the form To study the equivaraint fluctuations about the configuration (3.2), we impose the symmetry constraints, that is the equivaraince conditions, first of which means that the gauge field A is simply required to transform as a scalar of SO(5) under the adjoint action of W ab , which is quite natural as it does not carry any SO (5) index, while the second requires that the fluctuations F a introduced in (3.3) transform as a vector of SO(5) to comply with the equivariance condition. We infer from the decomposition (3.7) that the space of rotational invariants that may be constructed form the intertwiners of (0, n) and (0, 1) IRRs of SO (5) is three-dimensional, while the space of vectors that may be constructed in term of the intertwiners and X a is of dimension seven. The aforementioned intertwiners may be introduced via the projection operators to the IRRs appearing in the r.h.s. of (3.5). We may express these three projections as , P 2 I = P I , P † I = P I , I : 1, 2, 3 , and will be useful in what follows. Adjoint representation of SO(6) branches under SO(5) as 15 → 5 ⊕ 10, or in the Dynkin notation: (1, 0, 1) ≡ (1, 0) ⊕ (0, 2) . Thus, further insight on how SU (4) ≈ SO(6) generators sits in these intertwiners is gained by observing that Q I contain, ten of these generators as Σ ab , and the remaining five as Γ a , as seen from (3.12).
Using (2.30), commutative limit of (3.12) takes the form lim n→∞ (G · Σ) 2 n 2 = 8(x a Γ a + 1 4 ) . (3.14) Consequently, we find for q I := lim n→∞ Q I : We may argue that, equivariant parametrization of the fluctuations introduces topological fluxes through the concentric S 4 F s, preventing the latter to shrink to zero radius. Without going into any technicalities regarding the S 2 fibre coordinates ω ab over S 4 , using (3.15b), this reasoning is supported by the fact that in the commutative limit the topological flux piercing the S 4 F may be linked to the second Chern number on S 4 : where p 2 = 1 2 (1 − q 2 ) stand for the rank four projectors generating the projective module over the algebra of functions over S 4 .
We may solve the constraints given in (3.9a) and (3.9b) as follows. To satisfy (3.9a), we may choose to parameterize the gauge field A as where we have introduced α µ ≡ α µ (t) are real functions of time only, and eliminated Q 2 in favor of 1 4N using I Q I = −1 4N . From this form of the gauge field, it is readily observed that the U (4) gauge symmetry is broken down to U (1) × U (1) × U (1). However, later on we will see that term proportional to identity matrix in (3.17) does not survive the dimensional reduction and the gauge symmetry of the reduced actions will eventually be U (1) × U (1). For the fluctuations satisfying (3.9b), a convenient parameterization that befits the ensuing developments turns out be n ≡ X a n ,Γ a = 1 N ⊗ Γ a n ≡ Γ a n . (3.19) The 1 n factors appearing in the last three terms of F a via,X a andΓ a are naturally expected for the convergence of F a . Similar analysis performed in [36][37][38][39] on equivariant parametrization of fluctuations over S 2 F and S 2 F × S 2 F , shares the same features. Indeed, as n → ∞, we find Demanding the fluctuations f a to be tangential to S 4 means that x a f a = 0 has to be fulfilled. Since x a ∇ a = 0, as noted after (2.32), the latter condition is satisfied if and only if φ 3 (t), χ 3 (t) and φ 4 (t) all vanish in this limit.

Reduction over S 4 F and the Low Energy Effective Action
We are now in a position to exploit the explicit parameterizations given (3.17) and (3.18) to determine the low energy effective action that emerges from the action S 1 by tracing over the concentric S 4 F 's. After a straightforward but a long calculation, with some intermediate steps relegated to the Appendix B, we may write down the covariant derivatives in the form where we have used for the fields φ i (t) and χ i (t) for (i : 1, 2) only. Trace of the kinetic term is evaluated to be 1 2 T r(D t X a ) 2 = n(n + 4) (n + 1) 2 |D t φ| 2 + n(n + 4) (n + 3) 2 |D t χ| 2 + 2(n + 4)(n 5 + 8n 4 + 18n 3 + 8n 2 − 11n) where we have introduced the complex fields φ(t) = φ 1 (t) + iφ 2 (t) and χ(t) = χ 1 (t) + iχ 2 (t), with the covariant derivatives taking the usual form D t φ = ∂ t φ + iα 1 φ, D t χ = ∂ t χ + iα 3 χ.
Although (3.23) does not appear to be manifestly positive definite, a simple analysis using Mathematica confirms that it is so, as it should be by construction. Thus, it is possible and useful to make a linear field redefinition in the sector spanned by φ 3 , φ 4 and χ 3 , i.e. convert to a basis, in which the kinetic term (3.23) is diagonalized. In the next section, we will naturally work with such linearly redefined fields, which diagonalize the kinetic term for specific values of n from n = 1 to n = 5.

(3.24)
A few remarks are in order. Firstly, we see that there are terms which are linear in the fields φ 3 , φ 4 and χ 3 . For finite values of n, which is essentially going to be our main focus in the next section, these terms cause no harm. In the n → ∞ the coefficients of these terms diverge linearly with n. Thus, alluding to our previous remark following (3.20), for the finiteness of the large n limit, it will suffice to assume that φ 3 , φ 4 and χ 3 vanish faster than 1 n . The mass terms then converges to −2µ 2 (|φ| 2 + |χ| 2 ), while the kinetic term is given by |D t φ| 2 + |D t χ| 2 in this limit. Let us also note that, we have µ 2 = −8, since we are inspecting this term about the S 4 F configuration satisfying (2.10). The exact form of the constant term C(n) is immaterial; in the next section we will adjust the overall constant term in the reduced Lagrangians so as to set the minimum value of the potential to zero.
Analytic calculation of the trace of the interaction term 1 4 T r[X a , X b ] 2 in the action S 1 for the equivariant parametrizations (3.17) and (3.18) turns out to be quite a formidable task as it involves large number of rather complicated traces. As an alternative approach, we have evaluated the traces for this term for the values of n = 1, 2, 3, 4, 5 using Mathematica. These already correspond to reasonably large span of matrix sizes 4N = 40, 80, 140, 224, respectively and gives us ample information for exploring the dynamics of the low energy reduced action. This is what we take up in the next section.

Gauge Symmetry & the Gauss Law Constraint
The equivariantly reduced action obtained from S 1 is invariant under the U (1) × U (1) gauge transformations with the remaining fields φ 3 , φ 4 and χ 3 being real and thus uncharged under U (1) × U (1). We observe this manifestly from (3.23) and the mass term (3.24). The interaction term has the same gauge symmetry too by construction and can be manifestly seen from (6.1) at the level n = 1.
As the time derivatives of U (1) gauge fields α 1 (t) and α 3 (t) appear nowhere in the action, these fields have no dynamics of their own. Their equations of motion give us the Gauss law constraints: 1 2i We make the gauge choice α 1 (t) = 0 and α 3 (t) = 0, which essentially amounts to the reality conditions φ * = φ and χ * = χ. We can be more precise by first noting that the Gauss law constraints does not break the U (1)×U (1) gauge symmetry completely, but a residual Z 2 ×Z 2 remains. Indeed, writing φ ≡ (φ 1 , φ 2 ) = |φ|(cos θ, sin θ) and χ ≡ (χ 1 , χ 2 ) = |χ|(cos σ, sin σ), we may express the constraints in the form Therefore, the remaining Z 2 × Z 2 symmetry is encoded in the gauge functions as Λ 1 (t) = Λ 0 1 + πk and Λ 3 (t) = Λ 0 3 + πk, where Λ 0 1 and Λ 0 3 are constants and k ∈ Z 2 . This indicates that, for both of the gauge functions, Λ 1 and Λ 3 , we have more generally Due to (4.3), we have θ(t) = θ 0 + πk and σ(t) = σ 0 + πk, and (4.4) holds for both θ(t) and σ(t), as well. Having noted these points, in what follows we set φ 2 (t) and χ(t) to zero (i.e., we have both θ 0 and σ 0 set to zero). Then, the Z 2 × Z 2 symmetry is implemented by In Section 5 we consider the structure of the LEAs in the Euclidean time τ . Due to the Z 2 × Z 2 symmetry, we will be able to explore kink type solutions of the LEAs by imposing appropriate boundary conditions on φ 1 (τ ) and χ 1 (τ ) as τ → ±∞. Availability to impose topologically non-trivial boundary conditions on the latter can be attributed to the property (4.4) of the restricted gauge functions, which holds the same in the Euclidean signature. There is also the possibility of not imposing the Gauss law constraint as recently discussed for BFSS and BMN matrix models in [28]. This leads to presence of Goldstone bosons for the LEAs that we have obtained, as we will briefly discuss and demonstrate in Section 6.

Dynamics of the Reduced Action and Chaos
The explicit form of the equivariantly reduced Lagrangian for n = 2 is given below Let us note that in L (n) for n = 2, 3, 4, 5, we have i) performed the linear transformation among the fields φ 3 → φ 3 , φ 4 → φ 4 , χ 3 → χ 3 which diagonalizes the kinetic term (3.23), and dropped the 's in the final form, ii) have set µ 2 = −8 as we have already remarked to do so in the paragraph after (3.24), iii) have imposed the Gauss law constraints as we just discussed, iv) adjusted the constant in the final form of each L (n) , so that the associated potentials, V (n) , take the value zero at their minima, v) introduced an over-dot,(...), to denote the time derivatives and vi) set the coupling constant g to one, as it has no effect on the classical physics save for determining a global normalization in the energy unit.
For n = 1, the reduced action takes a simpler form as compared to the cases for n ≥ 2, which is given as where we have introduced Φ = φ 3 + φ 4 − χ 3 , which is the only combination of the constituent dynamical variables φ 3 , φ 4 and χ 3 upon which L (n=1) depends. This is to be expected in view of (3.8). For L (n=1) too, all items following (4.5) are performed as well, except the item i), which is already taken care of with the introduction of Φ(t). An important feature of the reduced Lagrangians is that their potentials are all bounded from below. Therefore, we infer that at any level n the equivariant fluctuations about the concentric S 4 F solutions, with µ 2 = −8 do not cause any instability. The absolute minima of the potentials are given in (B.7), (5.6), (B.8) and (B.9).
We find that the reduced Lagrangians, L (n) , have chaotic dynamics. One of the basic tools to probe the presence of chaos in a dynamical system is to compute the Lyapunov exponents, which measures the exponential growth in perturbations [44]. If, say, x(t) is a phase space coordinate, in a chaotic system the perturbation in x(t), denoted by δx(t), deviates exponentially from its initial value at t = 0; |δx(t)| = |δx(0)|e λ L t , λ L being the Lyapunov exponent corresponding to the phase space variable x(t).
The phase space corresponding to the LEA are 10-dimensional, except for the n = 1 case, and spanned by where p i are the corresponding conjugate momenta and the Hamiltonians, H (n) , are obtained from L (n) in the usual manner using H = p iqi −L. We have obtained the Lyapunov spectrum for n = 1, 2, 3, 4, 5 at various energies (as determined by the initial conditions) using numerical solutions for the Hamilton's equations of motion. For n = 1, dimension of the phase space reduces to 6 as easily observed from L (1) in (4.6). The table below summarizes our numerical findings for the largest Lyapunov exponent, λ max , for n = 1, 2, 3, 4, 5 at several different values of the energy E. We have shaded the n = 1 column in this table to indicate the noted differences between this case and the rest.  In the figures (Fig.1) through (Fig.20) given below the time evolution of the Lyapunov exponents, λ i (i : 1, · · · , 10) for n = 1, 2, 3, 4, 5 at the energies E = 20, 30, 100, 250 are given. In almost all of these cases, we observe that all the Lyapunov exponents rapidly converge to constant values, well before the computation time used for numerical calculations; convergence to constant values being slower only at lower energies. Our results indicate that chaotic dynamics is present even at energies as low as E = 20, 30 and becomes significantly expressed at energies E = 100 and higher.
Plots in (Fig.21) and (Fig.22) show the evolution of λ max at E = 100 and E = 250 for n = 2, 3, 4, 5. As n increases further a slow decrease in the value of λ max is clearly observed from these figures and also more clearly from (Fig.23), which also indicates that the decrease in λ max is approximately linear for n ≥ 3. If we keep on increasing n, while keeping the energy fixed, it is to be expected that λ max gets smaller and smaller. This is not surprising though since, keeping E fixed, in this case, would require us to choose initial conditions which are ever more closer to zero, and hence would not be of use to probe the chaotic dynamics effectively. Hence, it appears reasonable to expect that λ max will remain significantly larger than zero for increasing values of n, provided that the system has sufficiently high energy.

Kink Solutions
We now consider the matrix model in the Euclidean signature S E 1 given in (2.14) and the corresponding LEAs. The latter have multiple degenerate vacua supporting kink solutions, i.e. instantons in 1 + 0-dimensions, whose features are sketched out in what follows.

Kinks at level n = 1
The Lagrangian is given as where stands for derivatives with respect to the Euclidean time τ . We have also scaled out an unimportant factor of 5 2 compared to (4.6). We can easily see from (5.1) that there are three different pairs of vacua, which are given by the configurations Since either φ 1 or χ 1 vanish in these vacua, we infer that, the kink solutions could be of the type with topological indices (±1, 0) or (0, ±1) ∈ Z 2 ⊕ Z 2 . These are the familiar kink solutions [45,46]. Indeed, we find that the equations of motion are of the form which have the following solutions

Kinks at levels n ≥ 2
For n = 2, 3, 4, 5, the number of degenerate vacua increases. This may be expected due to the larger number degrees of freedom in the LEAs. Nevertheless, a similar structure in vacuum configurations to that of n = 1 is observed, and allow for the kink solutions. At n = 3, for instance, we have eight pairs of degenerate vacua, which are given as The equations of motion for L (n=3) are coupled non-linear differential equations, which are not easily to solve. We may look at the linearized system of equations around one of the minima. For notational simplicity, let us write (φ 1 , χ 1 , φ 3 , φ 4 , χ 3 ) ≡ (S 1 , S 2 , S 3 , S 4 , S 5 ) := S and also introduce S = S 0 + s, where S 0 is one of the vacuum configurations in (5.6) and s stand for the fluctuations. The linearized system of equations is given by where no sum over repeated indices is implied on th l.h.s and C i can be easily read off from (B.4). For, say, S 0 The leading order profiles of the solutions of these equations which are regular as τ → ∞ are given below, while the complete solutions are given in Appendix B. where c i (i : 1 · · · 5) are arbitrary constants. (5.9) provides asymptotic profile of the kink solutions near S 0 + .

Kinks in the Presence of the Chern Simons Term
Here we will confine our discussion only to the level n = 1. In the Euclidean signature, the reduced action obtained from S 3 takes the form Last term in (5.10) is the contribution coming from the reduction of S CS as can be clearly observed from the fifth order terms it contains. Due to presence of this term the potential of (5.10) does not have an absolute minimum, in fact, naturally, it is not bounded from below. Nevertheless, for µ 2 < 0, it is a matter of a simple calculation to see that the potential still have degenerate local minima. This still allows for kink solutions, which are, however, only metastable and will decay under sufficiently large perturbations. As a concrete example, we have for instance, for µ 2 = −1, the local minima occurring at and a kink solution to the equations of (5.10) is given by

Gauge Symmetry Revisited
In [28] Maldacena and Milekhin considered the BFSS model without imposing the Gauss law constraint, i.e. without the SU (N )-singlet condition on the physical states. This is based on the fact that the BFSS model with A = 0 is still well-defined even in the absence of the Gauss law constraint, at the expense that the SU (N ) is no longer a local but a global gauge symmetry group in this situation.
We observe that such a possibility is also valid and applicable to the low energy reduced actions that we have obtained in this paper. The latter have U (1) × U (1) gauge symmetry and the Gauss law constraint was given in (4.2), imposing the U (1) × U (1) singlet condition meaning that the complex fields φ(t) and χ(t) are uncharged, i.e. real, under the gauge fields α 1 (t) and α 3 (t), respectively. If do not impose the constraint we can still set the gauge fields α 1 (t) and α 3 (t) to zero in the LEAs at any level n. In this case, the LEAs have a global U (1) × U (1) symmetry, which is spontaneously broken by several different vacuum configurations, and hence imply the existence of Goldstone bosons in these LEAs. For instance, at the level n = 1, we have the action with three different vacuum configurations, as easily recognized from (5.2) spontaneously breaking the U (1) × U (1) symmetry. Thus, we immediately infer that the fluctuations around each of these vacuum configurations give rise to one Goldstone boson in the usual manner that it arises in an abelian gauge theory with degenerate vacua. As a concrete example, we have, with φ = 2 + σ 1 + iσ 2 , Φ = − 1 2 + ρ and denoting the fluctuation around χ = 0 still with χ, potential part of L (n=1) takes the form showing that σ 2 is massless, i.e. it is the Goldstone boson associated with this particular vacuum configuration. Similar analysis show the existence of Goldstone bosons for the other two vacuum configurations in (6.2). Finally, we note that for n = 2, 3, 4, 5 we infer from (B.7), (5.6), (B.8) and (B.9) that, at each value of n, there are eight distinct vacuum configurations determined by the values of φ, χ and the real fields φ 3 , χ 3 and φ 4 , each of which comes with a Goldstone boson.

Conclusions
In this paper we have studied the 5d mass-deformed Yang-Mills matrix model with U (4N ) gauge symmetry. We have found the exact form of the SU (4) ≈ SO(6) equivariant parametrizations of the gauge field and the fluctuations about the classical four concentric fuzzy four sphere configuration and used them to calculate the LEAs by performing traces over the S 4 F s for the first five lowest matrix levels. The LEA's obtained in this manner have potentials bounded from below, which indicates that the equivariant fluctuations about the S 4 F configurations with a tachyonic mass term (µ 2 = −8)do not lead to any instabilities. We have showed through detailed numerical computations that these reduced systems have chaotic dynamics and calculated their Lyapunov exponents at various energies. In the Euclidean signature, we have demonstrated that the LEAs support usual kink type solutions, i.e. instantons in 1 + 0-dimensions. The latter may be viewed as the residual topologically non-trivial configurations, linked the topological fluxes penetrating the concentric S 4 F s due to the equivariance conditions, and preventing them to shrink to zero radius .

A.1. Branching Rules & Relations Dynkin and Heights Weight Labels
Irreducible representations of SO(2k) and SO(2k − 1)can be given in terms of the highest weight labels [λ] ≡ (λ 1 , λ 2 , · · · , λ k−1 , λ k ) and [µ] ≡ (µ 1 , µ 2 , · · · , µ k−1 ) respectively. Branching of the IRR [λ] of SO(2k) under SO(2k − 1) IRRs follows from the rule The relationship between Dykin labels and highest weight labels for SO(5) IRRs is For the SO (6) IRRs the correspondence is given by Here we present the details of the calculation leading to (3.12). We have Multiplying both sides of (A.6) with abcdm and using ijklm γ i γ j γ k γ l = 24γ m as can be inferred from (2.21b) one obtains To handle the left hand side of (A.7) we make use of the determinant Left hand side of (A.7) then reads 4γ a γ b γ c γ d G ab G cd + 4γ c γ d γ a γ b G ab G cd + 4γ a γ c γ d γ b G ab G cd +4γ a γ d γ b γ c G ab G cd + 4γ c γ a γ b γ d G ab G cd + 4γ d γ a γ c γ b G ab G cd = 4(6γ a γ b γ c γ d G ab G cd − 12G ab G ab + 24γ a γ c G ab G cb ) , (A.9) where the second line follows after making use of {Γ a , Γ b } = 2δ ab for rearrangements. Therefore, (A.7) can be cast into the form (G · Σ) 2 = γ a γ b γ c γ d G ab G cd = abcdm γ m G ab G cd + 2G ab G ab − 4γ a γ c G ab G cb .
Employing G ab = 1 2 [X a , X b ], we can expand the first term in r.h.s of (A.10) as after using (2.21b). Substituting G ab G ab = −4n(n + 4)1 4N for the second term and simplifying the third term as and finally, combining all the terms together we get (G · Σ) 2 = 12 Γ a Γ b G ab + 8(n + 2)X a Γ a + 8n(n + 4)1 4N .  .17) and (3.18) in D t X a = ∂ t X a − i[A, X a ], we find