Analytic SU(N) Skyrmions at finite Baryon density

We construct analytic (3+1)-dimensional Skyrmions living at finite Baryon density in the SU(N) Skyrme model that are not trivial embeddings of SU(2) into SU(N). We used Euler angles decomposition for arbitrary N and the generalized hedgehog Ansatz at finite Baryon density. The Skyrmions of high topological charge that we find represent smooth Baryonic layers whose properties can be computed explicitly. In particular, we determine the energy to Baryon charge ratio for any N showing the smoothness of the large N limit. The closeness to the BPS bound of these configurations can also be analyzed. The energy density profiles of these finite density Skyrmions have \textit{lasagna-like} shape in agreement with recent experimental findings. The shear modulus can be precisely estimated as well and our analytical result is close to recent numerical studies in the literature.


I. INTRODUCTION
The characterization of the phase diagram of the low energy limit of QCD at finite Baryon density and low temperatures has motivated intense research in the last two decades, see [1] and references therein. Analytic models are scarce and new exact results are hard to obtain. A well known example is the (3+1)dimensional Nambu-Jona-Lasinio (NJL) model that shares some of the analytical difficulties of the low energy limit of QCD (see [2] for a review). Together with the uselessness of perturbation theory at low energy, this means that the complicated phase diagram of low energy QCD cannot be easily analyzed with the available analytic techniques (see [3][4][5] and references therein).
A remarkable feature of low energy QCD at finite Baryon density is that at low temperature very complex structures appear. When the Baryon density is increased, a phase that is commonly defined as nuclear pasta appears. In [6][7][8][9][10][11][12][13], the presence of "baryonic layers" was disclosed, which will be the main focus of the present paper. Such a name arises from the fact that most of the baryonic charge and energy density is concentrated within lasagna-shaped regions in three dimensions 1 . Many physical properties of these configurations are currently under investigation, such as the elasticity of nuclear pasta and their transport properties [10][11][12][13]. The high topological charge of nuclear pasta makes it hard to study analytically.
As powerful numerical techniques are available to analyze these configurations (see, for instance, [10][11][12][13] and references therein), why should one insist so much in finding analytic solutions? There are many reasons to strive for analytic solutions even when numerical techniques are available. Firstly, it could be enough to remind all the fundamental concepts that we have understood thanks to the availability of the Kerr solutions in General Relativity and of the non-Abelian monopoles and instantons in Yang-Mills-Higgs theory. Secondly, as in the present case, analytic solutions can disclose relevant physical properties of very complex structures which are difficult to analyze even numerically.
Until recently, these types of non-homogeneous condensates in the low energy limit of QCD in (3+1)dimensions could not be properly understood analytically. A further problem is that, computationally, the large N f and large N c limits must be addressed carefully [14,15]. One of the goals of the present paper is to shed light on the large N f behavior of these complex structures.
A simplified version of the low energy limit of QCD that encodes many relevant features is the (1+1)dimensional version of the NJL model, also known as chiral Gross-Neveu model [16][17][18][19]. Such a model possesses a crystalline phase at low temperature and finite Baryon density [20][21][22][23]. These results suggest that ordered structures must also appear in the low energy limit of QCD. At leading order in the´t Hooft expansions [24][25][26], the low energy limit of QCD is described by the Skyrme theory [27] (see [28,29] for reviews). Despite the Bosonic nature of the Skyrmion field U , its solitons represent Baryons (see [26,[30][31][32][33]).
Here, we will analyze the appearance of complex structures at finite Baryon density in the SU (N ) Skyrme model in (3+1)-dimensions. We will focus on the analytic computations of relevant physical properties, such as the energy density, the energy per Baryon and the shear modulus of nuclear-lasagna like structures living at finite density 2 . We will compute their corresponding scaling with N . We will combine the use of Euler angles for SU (N ) developed in [39][40][41] together with the use of nonspherical hedgehog Ansatz introduced in (see [42][43][44][45][46][47][48][49][50]).

II. SKYRME ACTION
The action of the Skyrme model in four dimensions is where R µ = U −1 ∂ µ U = R j µ t j with U ∈ SU (N ), t i the SU (N ) generators, and where K and λ are the Skyrme couplings, g is the the metric determinant, 3 and F µν = [R µ , R ν ]. The field equations are We will construct topologically non-trivial solutions at finite Baryon density. Our main goal is to determine the scaling with N of relevant physical quantities. As we want to analyze Skyrmions of high topological charge living in flat spaces at finite Baryon density, we will consider the following metric: ds 2 = −dt 2 + L 2 r dr 2 + L 2 γ dγ 2 + L 2 ϕ dϕ 2 , (II.3) 2 Pioneering works on the Skyrme model at finite density are [34][35][36][37][38] and references therein. 3 We remind the reader that the N of the SU (N ) of the Skyrme model corresponds to N f . while the range of coordinates is with the caveat that, despite the chosen values, they are not periodic! The parameters L r , L γ and L ϕ represent the size of the box within which the Skyrmion is confined.
A. Quantities of high physical interest.
Firstly, the main goal of the paper is to compute the energy per Baryon and its large N behavior. Therefore, only solutions with non-vanishing Baryon charge have been considered. The usual definition of Baryon charge in the Skyrme model (see [26,27,31,32]) is so that a necessary condition in order to have nontrivial topological charge is From the geometrical point of view, the above condition can be interpreted as saying that the Skyrmion "fills a three-dimensional spatial volume", at least locally. On the other hand, such a condition is not sufficient in general. One also has to require that the spatial integral of ρ B be a non-vanishing integer : (II. 8) Usually, this second requirement allows to fix some of the parameters and integration constants of the Ansatz, as we will see in the following. However, there are more global conditions to be satisfied, as it will be explained below. Hence, in the following we will only consider solutions satisfying both the condition in Eq. (II.7) and the one in Eq. (II.8).
Secondly, the energy density (the 0 − 0 component of the energy-momentum tensor) reads where F µν = [R µ , R ν ]. Thus, the total energy E of the Skyrmion is the spatial integral of the above quantity We define a Skyrmion U to be static if its energy density defined above is static. In other words, a Skyrmion is static, if it corresponds to a static distribution of energy density. It is worth to note that this definition is more general than the naive definition of a static Skyrmion as a static SU (N )-valued configuration U which does not depend on time. In particular an elegant approach to avoid Derrick's famous nogo theorem on the existence of solitons corresponds to search for a time-periodic Ansatz such that the energy density of the configuration is still static, as it happens for boson stars [51] (in the simpler case of U (1)-charged scalar field: see [52] and references therein). The Ansatz to be defined in the next sections will have exactly this property. Moreover, unlike what happens for the usual Bosons star Ansatz for U (1)-charged scalar fields, the present Ansatz for SU (N )-valued scalar field also possesses a non-trivial topological charge. Thus, we are interested in solutions in which the energy density has non-trivial local maxima, which could be identified with the position of the Skyrmions. Given a solution of SU (N ) with Baryonic charge B and energy E living in the metric (II.3) we have already mentioned, it is very interesting to analyze the following quantity (which is nothing but the energy per Baryon of the configuration g (N, a)) E B def = g (N, a) , (II. 10) where a is any set of integration constants which characterizes the given solution. It is especially interesting to understand the behavior of g (N, a) defined above when N is large (the´t Hooft limit). Here and in the following we will call the function g (N, a) the "gfactor". The very deep question is whether or not, in the given family of solutions one is considering, one can define and if this limit is well defined. In particular, one would like to know whether or not "the closeness to the BPS bound" improves when N is large. Indeed, it is worth to remind that in the SU (2) case all the known solutions with non-vanishing topological charge exceed the bound by at least the 20%. Hence, one would like to know whether, in the´t Hooft limit, the "closeness of Skyrmions to the BPS bound" is finite or whether it grows without bound. This issue is deeply related with the so-called Veneziano limit [25], which is a variant of the´t Hooft limit in which also the flavour number N f goes to infinity in such a way that N c /N f stays finite. The Veneziano limit allows to take into account the effects of quarks while keeping the advantages of the´t Hooft topological expansion. Since, to arrive at the Skyrme model as an effective low energy limit of QCD, N c must be already large, the large N limit which we are considering here (in which N is the one of the SU (N ) Skyrme model), can be considered as a sort of Veneziano limit applied to the Skyrme model itself. The fact that such a limit is smooth is a very non-trivial result which would be very difficult to prove directly on the QCD Lagrangian. The above discussion clearly shows that in order to declare a solution of the Skyrme field equations as "physically interesting" two criteria must be satisfied: 1) The topological charge of the solution must be non-vanishing 2) The energy density T 00 as function of the coordinates must have an interesting pattern.

III. LOCAL SOLUTIONS
Using the Euler angles for SU (N ) determined in [39,40] together with the Ansatz for non-spherical Skyrmions living at finite Baryon density in [42][43][44][45][46][47][48][49][50], one arrives at the following Ansatz for the SU (N ) Skyrmion: with a suitable choice of k in su(N ) and h(r) in the Cartan subalgebra H to be specified below, m a nonvanishing integer number, and where we recall that the metric is given by (II.3). When necessary to expand with respect to the basis of su(N ), we will also write with (see app. A) In general we will use the simplifying notations As for k, for c j arbitrary complex numbers, forming the components of the vector c ∈ C N −1 , we choose where the prime indicates derivation w.r.t. r.
The proof is given in appendix B. Exploiting (III.3) and (III.6) we can further simplify the equations of motion, which can be put in the following form.
where h.c. stays for Hermitian conjugate, and α j are a suitable choice of simple roots of SU (N ), defined in App. A 1. Indeed, using (B.3) and (B.5), we can rewrite (III.7) as (III.10) Now we use general properties of simple roots. Since λ j are eigenmatrices relative to simple roots, it happens that or [λ j , λ k ] = 0 or it is an eigenmatrix relative to a positive root. 4 Similar considerations follow for λ † j w.r.t. negative roots. It follows that none of these terms can lie in H, so that projecting equation (III.10) on H we get (h ′′ belongs in H by definition) we get (III.8), while projecting on the complement we get (III.9). These equations could be expressed even more explicitly in components, by exploiting (III.3) and us-ing that α j (J k ) = C AN−1 j,k are the components of the Cartan matrix of SU (N ), as defined in App. A 2, so that α j (h (n) ) = 2y j−1 . However, such explicit expression is not necessary in order to get the general solution.

A. Explicit solutions
We want now to find all the solutions of the equations (III.8) and (III.9). To these hand we will make use of some technical fact explained in appendix A 2. Let us first consider (III.9). Using (A. 16) (III.11) 4 That is a linear combination of simple roots with non negative integer coefficients We will assume c to be generic, with this meaning that all the c j are non zero. Since E j,j+2 , including their conjugates, are all linearly independent, this gives Since α j are real valued, we get also The firsts of these give We have two possibilities. Or h ′′ = 0, or not. We will now show that the second case leads to a contradiction. First, notice that if h ′′ = 0 then it must be α j (h ′′ ) = 0 for at least one j (since the α j are linearly independent), so that all α j (h ′′ ) are equal and different from zero. From the second of (III.12) we have that there must exist signs ε j such that (III.14) Deriving it w.r.t. r must give (III.13), so that ε j = 1 for all j, and we are left with the linear system of equations Since the α j are linearly independent (of rank N − 1) this is a set of N − 2 linearly independent equations for h ′ ∈ H. Since H is N − 1 dimensional, the space of solutions is one dimensional and the general solution of it is where f is an arbitrary function and v ∈ H is the unique matrix satisfying α j (v) = 1 for all j (which we will compute later, for now it is sufficient to know it exists). We now replace this solution in (III.8). We immediately get Since we have assumed h ′′ = 0, we have f ′ = 0 and, therefore, After applying α k to this equality, using that α k (v) = 1 and noticing that α k (J j ) = C AN−1 k,j are the components of the Cartan matrix, we get This relation can be inverted easily: if we consider 1 at varying j to be the components of a vector in R N −1 , we can apply the inverse Cartan matrix to both members, thus getting Since λ is positive and the same is true for the elements of the inverse Cartan matrix (A.15), we see that this led us to a contradiction. Therefore, the only possibility is that f ′ (r) = 0, which is equivalent to h ′′ (r) = 0. Hence, we proceed in investigating the first possibility, h ′′ = 0. In this case (III.8) is automatically satisfied and (III.9) reduces to (III.14). Its solution is where a is a constant and v ∈ H is the unique matrix solving α j (v) = ε j , j = 1, . . . , N − 1 where ε j ∈ {0, 1} (and ε 1 = 1). Since ε 1 is fixed, this gives 2 N −2 solution for every choice of c j in k. As we will see in the explicit example of SU (4), however, not all of these are really distinct solutions. There is a convenient way to express v explicitly. Indeed, let us write 5 h = arv ε , where a is a constant and v ε ∈ H is a matrix  5 We omit an irrelevant additive integration constant We have thus proved where a is a real constant and ε j are signs, with ε 1 = 1.
This solutions are only local solutions, which means that they solve the differential equations. They do not extend automatically to global solutions, that are solution with a well defined Baryon number. Looking for global solutions is the task of the next section.

IV. GLOBAL SOLUTIONS
Up to now we have found the most general solution of the differential Skyrme equation. Nevertheless, it is not sufficient to determine a Skyrmion, since global conditions have to be imposed in order to get a solution with a well defined topological charge. This condition is not simply equivalent to impose that the topological charge must be integer (this is just a consequence of the right topological condition) but that it has to wrap a homological cycle an entire number of times (mathematically, it has to cover a cycle, that means to be a surjective map with a well defined degree). We will normalize the parametrizations so to have all ranges in [0, 2π].

A. Statement of the problem
The difficulty in passing from local solutions to global solutions is twofold. In order to illustrate it, let us consider the specific example of SU (4) when k is given by c i = 1. For getting a well defined global solution, the function is expected to provide a good coordinate of the image, of the solution. Since the target space of the map is compact, this requires that if we extend the range of γ to the whole R, g(γ) must result to be a periodic function. Now, a simple calculation show that the eigenvalues of k are ±µ + , ±µ − , with This means that, for a suitable unitary constant matrix U , we have In particular, its elements have periodicities T ± with is not rational, they have not a common period and the orbit never close, so is not a periodic function but, rather, its orbit describes a curve which densely covers a bi-torus in SU (4). In particular, it is not possible to use g(γ) as a good factor to get a finite covering of a cycle, despite it gives a solution of the equations of motion. It doesn't provide a solution with a well defined topological number and must be discarded. One has to tackle the problem of looking for acceptable matrices k, t hat are matrices generating a well defined period.
Assuming we have solved the periodicity problem, there is a second subtlety to be tackled: how to determine the right range of the coordinates in order to correctly cover a cycle. First notice that π 3 (SU (N )) = Z. This suggests that homotopically we have just one representative for any given topological (Baryonic) charge. Moreover, since π 2 (SU (N )) = 0, we have also H 3 (SU (N ), Z) = Z, so we have also a unique homological representative. Nevertheless, the solutions have not to be identified under deformation, but at most under gauge equivalence. But since the action is not gauge invariant, in our case all different representatives in a given equivalence class must to be considered as different solutions.
We will distinguish three different classes of solutions. The first two classes have canonical representatives: the ones of SU (2)-type, which belong in every class, and the ones of SO(3)-type, which belong in even classes only. They can be simply understood as follows. For any given N we can embed the representations of su(2) into su(N ). Exponentiating, they will give realisations of SU (2) or SO(3), depending on the specific representation. This give rise to pure SU (2)type or SO(3)-type solutions. However, they can be continuously deformed, by varying the corresponding c when allowed, giving rise to solutions that are not embeddings, so we can consider them as true SU (N ) solutions. But there exist a third class of solutions that cannot be obtained as continuous deformations of embeddings. Their existence is due to the fact that SU (N ) has a center isomorphic to Z N , which acts continuously on SU (N ), see appendix A. In particular, if Γ is a normal subgroup of the center, then one can construct the group SU (N ) Γ := SU (N )/Γ. The new class of solutions are generated by cycles in SU (N ) that reduce to cycles of SU (N ) Γ after the quotient. We will call them genuine SU (N ) solutions. We will consider them carefully in the explicit examples of SU (3) and SU (4), where everything is exactly computable, but now we shortly describe the SU (2)type and SO(3)-type, where some details are a priory known, see App. C. An SU (2)-type cycle has the form where h ′ is constant and the coordinate must run as follows. The range of r must be T /4, where T is the period of e h ′ r . The range of γ must be T k , the period of e γk (with m = 1!), and the range of φ must be T k /2. Therefore, the convenient choice for the coordinates is corresponding to the Baryon number where B 0 is the fundamental charge of the given Skyrmion. For SO(3)-type cycles the interval for φ must cover an integer period, so that the ranges must be and the corresponding Baryon number is The SO(3)-type can be defined as "di-Baryon class" after the seminal works [31] [32]. These results were extended, keeping spherical symmetry, to the SU (N ) case in [53] [54] [55] [56] leading to numerical non-embedded configurations in the SU (N ) Skyrme model. In the present paper we will generalize those findings to the non-spherical case at finite Baryon density achieving, moreover, analytic solutions.

B. SU (3) Skyrmions
Let us apply the above formalism to the case N = 3. In this case we will see that the problem of periodicity will not arise.

SO(3)-type solutions and genuine SU (3) solutions
The matrix k is We put c 2 = |c 1 | 2 + |c 2 | 2 . Then, the characteristic equation is The eigenvalues are λ 0 = 0 and λ ± = ±i c , so that is periodic with period Now, we pass to determine the Cartan element. We have two possibilities according to the two possible choices for ε: The inverse Cartan matrix for SU (3) is Thus, we find the two solutions The period of exp h + (r) is while the one of exp h + (r) is (IV.15) Now, we have to discuss the global properties in order to fix the ranges of the parameters. To this end, accordingly to appendix C, we have to look for the intersection between the orbit of h ± and the one of γk c .
Using the characteristic equation we immediately see that, see App. E 1, so that the intersection we are looking for is just the unit matrix I. However, we can notice that the orbit of exp h − (r) contains the elements In order to correctly define the solution we thus have to identify the ranges as follows. First, it is convenient to normalise c so that c = 1. This is just equivalent to re-scale the coordinates Φ and γ. Therefore, we fix once for all the metric to be with the caveat that, despite the chosen values, none of the coordinates is periodic! Our Skyrmions are living in a rectangular box. SO(3) type solutions. We already know that r must cover 1/2 of the period of the Cartan torus, which implies that we have to fix a = 1 2 . Hence, our solutions are We can now compute the energy end the factor g + = E 2m . We omit details here, since are particular cases of the general one for generic N considered below. We get where |c 1 | 2 + |c 2 | 2 = 1. In particular, for each value of m, |g + (m, c)| takes its minimum at |c 1 | = |c 2 |, which is Some comments are in order now. The reason for which the solution we have just described are of SO(3)-type can be understood remembering that we are working with 3 × 3 matrices, which carry naturally a representation of spin 1 of the rotation group.
Indeed, the minimum energy case just discussed, in which |c j | = 1/ √ 2, corresponds exactly to the case when the matrices h + and k c are the generators of the group SO(3) in the representation of spin 1. The other solution, for every fixed m, are continuous defor-mations obtained varying c, which does not changes their topological nature, and in particular the Baryon number, but it changes the energy. One can easily check that for generic c the matrices h + and k c do not generate a subgroup. One may wander if this is related to the fact that their energy is not a minimum. The present remark suggests how to look for SU (2)type solutions. Genuine SU (3) type solutions. Since this case does not enter in the canonical classes, we have to manage separately the determination of the correct ranges (then normalised to 2π as specified above). As for r, we will prove in proposition 3 that in order to have r ranging in [0, 2π], one has always to fix a = 1 2 . For what concerns the other coordinates, let us notice that h − (r) does not commute with k c but it commutes with k 2 c . Therefore, for g(γ) = e γkc , we see that g(T k /2) commutes with e h−(r) . This means that we can write shows that we are covering it twice unless we restrict one of the two ranges, of Φ and of γ, to one half the period of g. We choose to reduce Φ, where B has been computed as in App. F. Explicitly, For U − , g results to be independent from c:

SU (2)-type solutions
It is now clear that in order to find SU (2)-type solutions we have to consider deformations of spin 1 2 representations. This can be obtained by "reducing matrices" down to 2×2, and can be achieved by choosing where c is a phase. This is not the same thing as simply putting c 2 = 0 in k c in the sense that we have to choose k = k c before solving equation (III.9). Indeed, in (III.9) we assumed that all simple roots enter the game. This fixes the set of possible choices of h(r), and if in the above solutions we deform smoothly c to (c, 0), we cannot move away from our topological classes. This is confirmed by the fact that if we put c 2 = 0, the matrix k reduces to a 2 × 2 matrix, but the k ± do not allow to reduce the representation down to C 2 . We have to make a discontinuous deformation. The point is that for c 2 = 0 the root α 2 does not enter into equation (III.9) that, indeed, for N = 3 becomes just an identity. This means that when c 2 = 0 we can choose for h(r) any combination with the only caveat that e h(r) must be periodic, so that a and b must be in rational ratio. We can set For q = ±1 we fall down to the previous SO(3)solutions, while, of course, q = 0 provides a canonical embedding of SU (2) into SU (3), thus identifying an SU (2)-type solution. It is worth to mention that, since q ∈ Q it cannot be deformed continuously among the three values, compatibly with the fact that the case q = 0 is not in the same topological class of the other ones and, indeed, we may wander what happens for all the other values of q, since they would generate new genuine SU (3) solutions. However, it results that they have vanishing Baryon number, so that we will not consider them further. Thus we get the solutions (IV.34) The 1/2 factor in the first exponent has been added to ensure that when Φ varies in [0, 2π] it covers half of the period. Finally, we can compute the factor g: We will now consider the class of Skyrmions associated to the matrix k given by We will limit ourselves to the case when all the c j are different from zero. Here, we have to face the problem of establish for which choices of c j the matrix e γkc is periodic. By now, let us assume to have solved it and write down the corresponding solution: where σ = 1 for SO(3)-type solutions and σ = 1/2 for SU (2)-type solutions, and v ε is given by (III.21). For general genuine solutions the value of σ must be computed case by case. For any admissible c these are 2 N −2 solutions (since ε 1 = 1). In principle a could depend on N and ε. However, we will now show that this is not the case and the value of a is completely fixed by requiring that the normalized interval [0, 2π] for r must have the extension necessary to cover once a cycle: and the corresponding map U c ε [t, r, ϕ, γ] has not to cover a cycle more than once, then necessarily a = 1 2 . Proof. The proof is simply based on the same strategy used for example in [39]: one first constructs the invariant measure restricted to the hypothetical cycle; the resulting measure will depend explicitly on some of the coordinates and will vanish at specific value of that coordinate. The good range for such a coordinate to cover just once a cycle is any range between two vanishing points. The nice fact is now that the Haar measure restricted to a cycle, a part from an eventual normalization constant, is just ρ B , which is computed in App.F. Since it results to depend on r via sin(ar), we see that a suitable good interval for r is [0, π/a]. Since we want it to be [0, 2π], it must be a = 1 2 . Therefore, we definitely have a = 1 2 (IV. 41) in any case. Now, we can compute the g factor for our solutions. To this end, first note that according to Appendix B, and we used According to (B.6), (B.7), (B.8), with a = 1 2 , we have Replacing in the expression for T 00 and using that the energy is In a similar way one can compute the baryon number. This is done in appendix F, with the result From these results we immediately get the g-factor: Up to now, we have assumed c to be normalised so that g(γ) = e γkc has period 2π. However, we will not have really found a solution until we will be able to specify for which c the function g is periodic. Therefore, we cannot further postpone to tackle this prob-lem. However, before considering it in general, we want now concentrate on a very particular case, when ε j = 1 for all j. In this case It is clear that, among all possible choices for ε j , this minimises the energy, apart from possible effects due to v . We want also minimise with respect to the c j , assuming the normalisation of c fixed. Introducing a Lagrange multiplicator Λ, we have to extremize the function Deriving with respect to |c j | 2 we get the system 1 being the vector in R N −1 having all elements equal to 1. This gives the solution Interestingly this also solves automatically the periodicity problem. It is easy to see (App. D) that where ζ j are arbitrary phases, give a matrix e γkc that is periodic in γ with period 2π.
Moreover, we have Proof. The first result follows immediately by the well known formulas For the second expression notice that, by using (A. 15), where f , given by (IV.55). From this we get Now, f is the sum of two homogeneous pieces, one of degree 4 end the other of degree 2. Therefore, we can use the Euler theorem 6 to rewrite the last as which completes the proof. ✷ Using this results and noticing that σ 2 Λ = 1/2, we find for the energy per Baryon We can also notice that is the Dynkin index of the given representation of the principal representation of sl(2) in sl(N ), so, the fundamental Baryonic charge associated is Notice that for N odd I N is even, so B is always integer.
Finally, we are also interested in minimizing expression (IV.69) with respects to L a , a = ϕ, r, γ. This is done in general in appendix G. By using the formulas therein and the ones in the last proposition, we get that the minimum is reached at , with corresponding minimal value for an homogeneous function f : R N → R of degree L one has x · gradf = Lf.
Using normalised units (corresponding to λ = 1 and K = (6π 2 ) −1 ) we get Notice that this is independent from N and it is expected to be the absolute minimum with respect to any choice of ε j . We will not try to prove this conjecture here, we will limit ourselves to check it for N = 4 here below. The comparison with [57] is very interesting. The present results are slightly above the bound in [57] due to the time-dependence in the Ansatz. Note however that the present time-dependence cannot be undone as the present solutions wrap in a topologically non-trivial way also around the time direction. To the best of our knowledge, this is the first analytic computations showing explicitly how the closeness to the BPS bound "evolves" with N in the SU (N ) Skyrme model.
To be more specific, as it has been already emphasized, we are interested in topologically non-trivial solutions. In the present context this means that we only consider SU (N ) Ansatz such that As it has been discussed in the previous sections, ρ B represents the Baryon density when it is nonvanishing along three-dimensional space-like hypersurfaces Σ t=const .
In these cases, the integral of ρ B over Σ t=const represents the Baryon charge. While, mathematically, these integrals represent how many time the SU (N )valued Skyrmions wrap around Σ t=const . On the other hand, ρ B can be topologically trivial also along timelike hypersurfaces. In this case, one can also consider the wrapping of the SU (N )-valued configurations along three-dimensional time-like hypersurfaces. The configurations which have been constructed here are, as a direct check easily reveals, topologically non-trivial in two ways. Not only they possess nonvanishing Baryonic charge, they are also wrapped nontrivially along time-like hypersurfaces. Indeed, if one considers U c ε [t, r, ϕ, γ] = e σΦkc e avεr e mγkc , then the corresponding topological density has one space-like component and one time-like component: In particular, it implies that these SU (N ) Skyrmions wrap non-trivially around the three-dimensional timelike {ϕ = const} hypersurfaces. The consequence of this fact is that the time-dependence of the present configuration "cannot be undone" otherwise the winding number corresponding to the {ϕ = const} hypersurfaces would change.

D. Solving the periodicity problem
The solution of this problem is provided in App.E. 7 We discuss here the main results. The vectors c ∈ C N −1 having all components different from zero and allowing for a periodic function g(γ) = e γkc , with period 2π, form a family Following App.E 3, for any two coprime positive integers p and q such that p > q, for N = 4 we can find four families of solutions, each one parametrized by three real phases α 1 , α 2 , α 3 and a real modulus τ ∈ [q, p]. Each of these families is specified by one of the four possible inequivalent choices for the discrete vector ε. Recall that in this case the inverse Cartan matrix is We also have, see App.E 3, which is exactly the matrix representing the diagonal generator of SU (2) in the spin 3/2 representation. However, this is not true in general and we will see that in this series only the one with (p, q) = (3, 1) is deformable to an SU (2) embedding. Let us first look at the coordinate ranges. Regarding the range of r, it is completely fixed by proposition 3. As for the remaining ranges, they must correspond to the period of g 4 , unless there are (finite discrete) subgroups of the U (1) group generated by g 4 , which commute with v a . Since v a does not commute with k j c , j = 1, 2, 3 7 S.L.C. is particularly grateful to Laurent Lafforgue for suggesting him how to tackle this problem in full generality. for some integer ℓ. For the x satisfying this condition, call them x p m, one has for f 3 which is zero for x = jπ for some integer j. Since our coordinates are forced to vary in [0, 2π], the only non trivial possibility is x = π. Putting this back into the previous condition, we must also have p − q = 2ℓ, (IV.89) which means that, since p and q are coprime, this happens only when both p and q are odd. In this case g 4 (π) = −I. Therefore, we see that for p − q odd, there are no discrete symmetries, and the ranges of Φ and γ must coincide with the whole period, so σ = 1. Instead, for p − q even we have g 4 (Φ + π)e so we see that to any point on the image there correspond two different coordinates (Φ, γ) and (Φ + π, γ + π), unless we restrict one of the two ranges to half a period. We choose to do it with Φ, and, in order to keep its range to be [0, 2π], we fix σ = 1/2.
The field is The Baryon number is B a = 2σ p−a m(p 2 + q 2 ), while for the g-factor we get The corresponding minimal energy, expressed in normalized units, is given by (G.16), which in this case becomes: We can further minimise w.r.t. τ . Setting x = τ 2 , we have to find the stationary points in (IV.96) Deriving the expression in the square root w.r.t. x and multiplying by (p 2 + q 2 ) 2 x 3 /6, we get the equation This gives the admissible solutions (x is positive) (IV.98) x 0 is always present, while x ± are stationary points only when the square root is real, that is when (IV.99) Setting z = p/q this means x 4 − 14x 2 + 1 > 0 so (since p/q > 1) and, finally, Taking the second derivative of the above expression and evaluating it in x 0 , we get that x 0 is the absolute minimum (at fixed p and q) if that is (recalling p ≥ q), for 1 ≤ p q < 1 2 (9 + √ 77), (IV.103) otherwise the minimum is placed in x ± . In conclusion The absolute minimum in the family is the minimum of the first row. Setting x = pq/(p 2 + q 2 ), we see that 1 + 10x 2 − 6x has a minimum for x = 3/10, which correspond to p = 3, q = 1. The corresponding absolute minimal energy is exactly (IV.73). This is not surprising at all, since the (p, q) = (3, 1), ε = (1, 1, 1) corresponds to solution (IV.69) for N = 4 (use (IV.57) with Λ = 2 in (E.59)). This corresponds to the undeformed SU (2) embedding, as anticipated.

The case
In this case we get Regarding the ranges, we can do de same exact reasoning as for the previous case, so we get for p − q even . (IV.108) The Baryonic charge is For the g-factor we get (IV.110) The corresponding minimal energy, given by (G. 16), in this case becomes: (IV.111) We can further minimise w.r.t. τ . Setting x = τ 2 , it is immediate to see that in this case the minimum is reached for to which it corresponds the value For fixed q, this is a monotonic decreasing function of p, so there is no an absolute minimum in this family. However, notice that the lower bound is  Reasoning as before, we see that the field is now The Baryonic charge is For the g-factor we get (IV.119) The corresponding minimal energy, given by (G. 16), in this case becomes: (IV.120) This is independent on τ and for fixed q it is a monotonic increasing function of p. It follows that the lower bound is reached for p = q = 1 (the value 1 is enforced by the request that p and q are coprime, but the result depends only on p/q) g c,bound = g c (1, 1) = π 6 (1 + 2 √ 2), (IV.121) 8 We are using the convention that s s s indicates the representation of spin s which, again, is (IV.73). However, this is not allowed, since for p = q = 1 the functions f j are not periodic and the solution of the equations does not yield a well defined global solution! In this particular family the absolute minimum is instead Under this map the inverse Cartan matrix is invariant and ε b → −ε d ≡ ε d , where the last equivalence is by a global rescaling. This sort of duality makes the two families perfectly equivalent and giving the same minima.
Remark: We see that of the four predicted sequences of families the true inequivalent ones are the first three ones, while the d case is not really new. It is natural to expect that such duality extends to any N , but this would require a deeper understanding of the global properties of the relevant moduli space W . To this aim, it would be interesting to investigate the explicit cases N = 5 and N = 6. This, however, goes beyond the scope of the present work.

V. SHEAR MODULUS FOR LASAGNA STATES
On the crust of ultra compact objects, like neutron stars, nucleons form large structures called pasta states. Knowing the elasticity properties of the crust may be very important to understand the structure of the gravitational waves emitted in a collision with a black hole. An important recent result has been found in [12] where, using numerical simulations based on the phenomenological nucleon-nucleon potential, the authors showed that the shear modulus for nuclear lasagna can have a value much larger than previous estimates. Here we give a first principle explanation of it as an application of the Skyrmionic model. To compute the shear modulus associated to lasagna states, our strategy will be to first compute it for the SU (2) case for the solutions determined in [46,49], by employing its relation with the 1 + 1 computations presented in [58].
Let us begin with a review [49]. We will consider the symmetric case 9 in Equations (13) and (16) of [49], namely This means that we are considering configurations in which the SU (2) Skyrmions live in a box of volume V tot , where l 1 is the length along the r direction (which is the coordinate of the profile H in Eq. (13) of [49]). The Baryonic charge corresponding to the Ansatz in Equations (12), (13) and (14) of [49] is (see below Eq. (24) page 5 of [49]). Then, the SU(2) field equations for the Ansatz in Equations (12), (13), (14) and (16) of [49] with a static profile H = H(r), reduce to where B A can be interpreted as the Baryon density per unit of area of the Lasagna configuration (up to π factors). In order to compare directly the present results with the ones in [58], it is convenient to define the rescaled coordinate y as follows Now, equations (V.5), (V.6) and (V.7) (which are equivalent to the results in [49]) can be compared directly with equations (2.4), (2.7) and (2.9) of [58]. In particular, the dictionary between the results of [58] and the present ones is where the left hand side (with respect to the "→") is from [58] while the right hand side comes from the above equations. Equations (2.9) and (2.10) of [58] read L = 2τ I −1/2 (τ ) , which depends on the Baryon charge as well as the size of the box in which the configuration lives. Now, with the above dictionary, we can write the speed of sound of the phonons using Eq. (3.15) of [58]: where the T 00 is given in Eq. (28) of [49]. Thus we have the following expression for the shear modulus G SU (2) in the SU (2) case We can then estimate it as follows. In place of T 00 we use its mean value computed as where E 0 min is the minimal energy corresponding to B = 1. From Table 1 of [58] we see that B/A is independent from B for the minimal energy configuration.
Using the values in the table 10 we get T 00 ≃ 1.26 10 34 erg/cm 3 . 10 Notice that with these values the baryon density is n ≃ 0.0468f m −3 ≈ 0.05f m −3 , the same value used in the simulations of [12] With the same values, from (V.4) we obtain Γ ≃ 0.371, πΓ ≃ 1.166.
Therefore condition (V.12), which is easily solved numerically after noticing that I − 1 2 (τ ) = K(τ 2 ), the first complete elliptic integral, gives Notice that the present value is expected to an approximation from above, since we are using a Skyrmionic effective model. From the above analysis, taking into account (IV.73), we can infer that in any case the true value should be G SU(2) 10 31 erg/cm 3 . The comparison with [12] is very good especially taking into account that we only used the Skyrme model. At this point we can use the new solutions found in the present work to relate the shear modulus for SU (N ) case to the one for SU (2). Let us consider the minimal energy per nucleon (IV.72). After multiplying by B and dividing by the volume, which, because of (IV C) is proportional to λ 3 2 , we getT On the other hand, the Baryon density is which solved for λ and replaced inT 00 gives Assuming the speed of sound to be essentially independent from N , as suggested by the fact that all the component of T µν scale in the same way with N , we get that the dependence of the shear modulus from N is so that we get the final estimate for the value of the shear modulus G SU(N ) of the SU (N ) Skyrme model as

VI. CONCLUSION AND PERSPECTIVES
In conclusion, we constructed the first examples of analytic (3+1)-dimensional Skyrmions living at finite Baryon density in the SU(N) Skyrme model (which are not trivial embeddings of SU (2) into SU (N )) for any N. These results allow to compute explicitly the energy to Baryon charge ratio for any N and to discuss its smooth large N limit as well as the closeness to the BPS bound. The energy density profiles of these finite density Skyrmions have lasagna-like shape. A quite remarkable by-product of the present analysis is that we have been able to estimate analytically the shear modulus of lasagna-shaped configurations which appear at finite Baryon density. Our estimate agrees with recent results [12] based on many body simulations in nuclear physics using phenomenological nucleon-nucleon interaction potentials. In this section we collect some general facts we applied for getting the solutions. Let V, (|) the Ndimensional complex vector space isomorphic to C N , endowed with the canonical hermitian product We see from the definition that su(N ) is a real vector space of dimension N 2 − 1. A basis can be easily obtained as follows. For any j, k = 1, . . . , N we define the matrix E j,k with elements They are called the elementary matrices. With these notations, a basis of su(N ) is given by In particular, the matrices J h span the Cartan subalgebra H.

Roots and simple roots
A concept that is particularly helpful for most of the calculations we need is the one of roots. These are related to the above observation regarding the diagonalizability of ad X for any X ∈ H. The diagonalizability must be checked on su(N ) ⊗ C, which is generated by the complex span of the basis given above, in place of the real span. Notice that the complex span contains the matrices E i,j , i = j. This is sufficient to determine all the eigenvectors and eigenvalues of ad X for any given X ∈ H. To this aim, let us specify H as follows: where with · · · R we mean the span over R of · · · . Thus, we immediately see that so that E j,k and J h are eigenmatrices of the adjoint action of X, with eigenvalues i(c j − c k ) and 0 respectively. The point is that the eigenvalues depend linearly on X. Let us consider the linear operators L j , j = 1, . . . , N defined by Then, we can write ic j = L j (X) so that The linear operators are said the non vanishing roots of su(N ). The corresponding eigenspaces are one dimensional. Beyond these, there is the vanishing root defining the 0 eigenvalue, which eigenspace is H, so that has degeneration equal to the rank r = N − 1.
In particular, the set of non vanishing roots contains a set of r linearly independent roots, having the property that all the remaining roots are are combination of them with all non positive or all non negative integer coefficients. These are called the simple roots and are Finally, for convenience, we introduce the less conventional concept of real valued roots α j,k = −iβ j,k , α j = −iβ j , which we will simply call again roots and simple roots. With this convention, for the simple roots α j , we can also write useful for practical purposes. This also shows that the α j are linearly independent. We name the corresponding eigenvectors λ αj ≡ λ j = E j,j+1 , so that

Some further technical facts
There is a canonical way to introduce a scalar product on the real space spanned by the simple roots. We however bypass the historical construction and employ (A.10) to define the scalar product (α j |α k ) := −Tr(J j J k ). (A.13) On H * R := α 1 , . . . , α N −1 R it is an euclidean scalar product. One then defines the r × r Cartan matrix 13 C AN−1 with components The Cartan matrix is strictly positive definite. Indeed, for any vector (x 1 , . . . , x r ) ∈ R r we have j,k which is strictly positive and vanishes only for x j = 0 for all j. In particular, the Cartan matrix is invertible and, indeed, one easily checks that Another important fact to notice is that for j, k one has Appendix B: Proof of Proposition 1 In order to prove the proposition, it is convenient to work with the coordinates Φ and T = t + L ϕ ϕ. The metric takes the form ds 2 = −L ϕ dT dΦ + L 2 r dr 2 + L 2 γ dγ 2 . With the given Ansatz, after replacing Φ with σΦ for constant σ (for convenience), for the L µ we get For F µν = [R µ , R ν ], with x = e −h(r) ke h(r) , we get the non vanishing components Setting L µ := [L ν , F µν ], the equations of motion are Using that nothing depends on T and that there are no lower T components, these reduce to which proves the proposition.

Further details
Making use of (III.6) and (A.11) we can write Repeating the same calculation: Finally, The last two terms cancel after summation, while the first terms vanish for j = k, so we get where we have changed the order of commutators in the first sum and exchanged the name of variable in the second sum. Therefore Similarly, The first two terms can be treated as above, giving the contribution while the last two terms, after renaming the labels in the first of the sums, give the contribution which in components is We finally get

A further proposition
We want now to state another technical proposition: where ε j is a sign, j = 1, . . . , N − 1, and x := e −h ′ r k c e h ′ r . Then and Proof. First, we have where we used the notation λ j = E j,j+1 . Since λ j is upper diagonal, so is λ j λ k , hence Tr(λ j λ k ) = 0. Similarly, Tr(λ † j λ † k ) = 0. On the other hand This proves (B.6). Now, notice that where we used that α j (h ′ ) 2 = (ε j a) 2 = a 2 . This proves (B.7).
Despite these being very well known facts, in this appendix we want to discuss the difference between SU (2) and SO(3), since it is crucial to identify our solutions. Locally, the two groups coincide, they have the same Lie algebra. However, SU (2) is simply con- a, b, c are the Euler angles. Each of the exponentials has a period (depending on the normalisation of the matrices), say T 3 for a and c, and T 2 for b. The strategy to correctly cover G exactly one time is explained in [40] and works as follows. To be sure to cover G one integer number of times one first allow the coordinates to run each one in the respective period. This number, in general, is larger than one because of redundancies, due to two reasons. The first reason is that the central element, parametrised by b, is chosen in the maximal torus (the exponential of the Cartan matrix). The redundancies correspond to the action of the Weyl group to the torus. This action is determined by the algebra and is the same for both SU (2) and SO (3). It shows that indeed moving b along a period quadruplicates the determination of the points for SU (2) and duplicates for SO (3), and one can reduce the range of b down to T 2 /4 or T 2 /2 respectively. At this point, the difference between SU (2) and SO(3) appears. Indeed, for SO(3) this is the end of the story, it is already covered just one time, while for SU (2) it remains a redundancy and we covered it twice. This redundancy is due to the fact that . (C.4) Therefore, since ∆ 2 = I g(a, b, c) = e aτ3 e bτ2 e cτ3 = e aτ3 e bτ2 ∆ −2 e cτ3 This redundancy is eliminated by reducing the range of a down to T 3 /2 for SU (2). This is the way, relevant to our case, to distinguish the two kind of solutions: if the above intersection id ∆, then the ranges of the variables a, b, c are T 3 /2, T 2 /4, T 3 respectively, and the group is SU (2), otherwise the ranges are T 3 , T 2 /2, T 3 , and the group is SO (3). Finally, we want to add a final remark relevant for recognising genuine solutions: for SO(3) generator τ it happens of course that the orbit exp(xτ ) never meets the center, while if τ is an SU (2) generator, then exp(x/2τ ) is the only non trivial generator of the center of SU (2). No other elements of the center of SU (N ) can be met these kind of orbits.

Appendix D: Representations of SU (2) and periodicity
It is well known from representation theory that spin J representation of SU (2) has generators T 1 , T 2 , T 3 given by the N × N matrices, with N = 2J + 1 Each of these matrices is diagonalizable with eigenvalues given by the ones of T 3 . Since it follows that the periodicity of depends only on the eigenvalues and so all f j have the same periodicity, which is obviously 2π fo odd N and 4π for even N .
On the other hand, let us consider the matrices k c and g(x) = exp(xk c ). The possible periodicity of g depends on the eigenvalues of k c . It is easy to see that the coefficients of the characteristic polynomial of k c depend only on the |c j | 2 so the phases of the c j are irrelevant for the periodicity. In particular, this means that the matrix exp(xT 2 ) with has the same periodicity of f 2 (x).

Appendix E: Solving the periodicity problem
In section IV we have seen that for N higher than 3 there is a further difficulty to overcome in order to find a global solution: generically the matrix g(x) = e xk is not periodic and its orbit densely fills a torus of dimension strictly higher than one. This phenomenon corresponds to the fact that the one parameter subgroup g(x) is not indeed a Lie subgroup but only an imbedded subgroup. Therefore, for arbitrary choices of the coefficients c j , the matrix cannot be used to generate global solutions unless the corresponding g(x) is periodic. We will now tackle this problem in general. For the sake of completeness we will first show that no problems arise in the case N = 3.

The case N = 3
In this simple case we have The corresponding characteristic polynomial is The eigenvalues are therefore 0, ±i c , which are in rational ratios so g(γ) = exp(γk c ) is periodic, in particular, with period 2π/ c . For other purposes, we compute explicitly g(γ). To this aim, let us first notice that, by Cayley-Hamilton theorem, k c satisfies where I and O are the identity and the null 3×3 matrices. This implies k 3 c = − c 2 k c so that any power of k c can be reduced to a power lower than three. Hence for three functions satisfying g 1 (0) = 1, g 2 (0) = g 3 (0) = 0, since e O = I. Deriving (E.5) w.r.t. γ and using the characteristic equation, we get with the Cauchy conditions g 1 (0) = 1, g 2 (0) = g 3 (0) = 0 (so that g ′ 2 (0) = 1). From the first equation we immediately get g 1 (γ) = 1, while deriving the second one and replacing from the third, we get g ′′ 2 (γ) = − c 2 g 2 (γ), g 2 (0) = 0, g ′ 2 (0) = 1, (E.10) which has solution Finally, from the third equation we get (E.12) Therefore e γkc = I + sin( c γ) c k c + 2 sin 2 ( c 2 γ) c 2 k 2 c . (E.13)

The general case
One can in principle solve this problem as follows. Since k is anti hermitian, it can be diagonalized in C, with pure imaginary eigenvalues. Moreover, if λ is an eigenvalue, also −λ = λ * is. Therefore, if S is the integer part of N/2 (so that N = 2S or N = 2S +1 for N even and odd respectively), generically we have S distinct non vanishing eigenvalues. Let U be a unitary matrix such that where σ is the diagonal form of k, say for all j = 1, . . . , S, that means T λ j = n j 2π, with n j a positive integer (obviously, we assume T > 0) for any j = 1, . . . , S. Therefore, so that all pairs of eigenvalues must have rational quotients. Of course, this condition is satisfied for N ≤ 3, and any choice of c is allowed. But for N ≥ 4 we cannot choose the c j arbitrarily: only those values, such that k admits eigenvalues with rational ratios are allowed. Notice that c remains defined up to a real multiplicative constant: if t ∈ R, then k tc = tk c . The eigenvalues are the solutions of the characteristic polynomial 19) of degree N in x. Since k c is antihermitian, its eigenvalues are purely imaginary and, moreover, if µ is a nonvanishing eigenvalue, then also µ * = −µ is an eigenvalue. So the non vanishing eigenvalues are in pairs and, if N is odd, there is at least one zero eigenvalue. Moreover, since in the factorization of the polynomial the nonvanishing eigenvalues µ must appear in the factors (x − µ)(x + µ) = x 2 − µ 2 , we see that the general form of the polynomial must be P N (x) = (x 2 ) n + a 1 (x 2 ) n−1 + . . . + a n , for N = 2n, x[(x 2 ) n + a 1 (x 2 ) n−1 + . . . + a n ], for N = 2n + 1. (E.20) The coefficients a j are not the same for N odd and for N even, but it is convenient to keep the same name so that we can generically write the equation for the non vanishing eigenvalues as y n + a 1 y n−1 + . . . + a n−1 y + a n = 0, We can be more precise: Proposition 6. Using the notation j ≪ k for k − j ≥ 2, we have (E.23) Proof. It can be easily proven by induction. We have already seen it for N = 3. A direct computation shows that it is true also for N = 4, since P 4 (x) = x 4 + x 2 (|c 1 | 2 + |c 2 | 2 + |c 3 | 2 ) + |c 1 | 2 |c 3 | 2 . Now, assume it to be true for N and N − 1. Let k n be the matrix n × n defined as k c with components c 1 , . . . , c n . This way, we se k n as a submatrix of k n+1 obtained erasing the last row and column. Let Developing the determinant with the Laplace rule applied to the last row, we easily find The first addendum contains all the monomials of the stated form but the terms containing |c N | 2 . The second addendum contains all the terms of the stated form containing |c N | 2 . The proposition is proved.
So, for example, Notice that, assuming that all c j are different from zero, we have always a n = 0, so these are truly non zero eigenvalues. Now, condition (E.18) is equivalent to require that there must exist a positive real number z and n positive integers m j , j = 1, . . . , n such that the non vanishing eigenvalues of k c must have the form λ ± j = ±im j z. This happen if the solutions of (E.21) are At this point, we can notice that the coefficient of the above polynomial can be written in terms of the solutions as: . . . (E.32) a n = j1<···<jn (−y j1 ) · · · (−y jn ). (E.33) Comparing with the last proposition, we get the following set of equations for the |c 2 j | =: ζ j : This is a set of n equations in N − 1 real positive variables. We will now show that it has generically an (N − 1 − n)-dimensional space of solutions in the interesting region, which is for ζ j positive. To this end, we assume the generic situation where all m a are different, and order them in an increasing sequence m 1 < m 2 < · · · < m n . We will show later that the condition on the m a cannot be weakened in order to get periodic solutions. Then, we show that there is a simple solution on the boundary of the region of interest, which is (if N is odd we assume the null eigenvalue to be the last one, λ 2n+1 = 0) whereM is a (n − 1) × (n − 1) matrix whose first row has all elements equal to 1 and Since m 2 a < m 2 b for a < b, we see that this determinant is different from zero. The proof of the claim then is an immediate consequence of the implicit function theorem.
Going back to the c j , we then see that in general c j = ξ j ζ j (m, t), (E.52) for arbitrary phases ξ j , j = 1, . . . , N −1, with m ∈ N n > , t ∈ W ⊂ R N −n−1 . The parameters t j parametrize the above family of solutions. We can always assume that the integer m j are coprime. Indeed, if m is a common divisor of m j so that m j = ms j , then we can write m = ms and m can be reabsorbed in z.
Having assumed this, we can now fix z in such the way that e xkc has period 2π. Indeed, since the non vanishing eigenvalues of k c are λ ± j = ±izm j , since the m j are coprime, the common period of the associated exponential is 2π/z. This fixes z = 1. Notice in particular that in this case ζ j1 · · · ζ j k = a1<...<a k ≤n m 2 a1 · · · m 2 a k , k = 2, . . . , n.
One says that these matrices have a moduli space 57) where T N −1 is the torus generated y the phases and W ⊂ R N −n−1 is the moduli space of the system. It is difficult to say something of general about the global properties of W . We will study in general the case N = 4 where all computations are feasible explicitly. Remark: for N = 3 we have n = 1 and, therefore, only one integer m that must be equal to 1 (to be "coprime"). So c must have norm 1 and the fundamental Baryon number is B = 2σ.

The N = 4 case
Let us apply the above results to the case of SU (4). We have n = 2, so we expect the dimension of W to be 1. The eigenvalues equation for k is 0 = λ 4 + λ 2 c 2 + |c 1 | 2 |c 3 | 2 . (E.58) The four solutions are iλ + , iλ − , −iλ + , −iλ − , with Let q ≤ p a pair of positive coprime integer numbers. Then, we have to solve the system ζ 1 + ζ 2 + ζ 3 = p 2 + q 2 , (E.60) Notice that this gives Now, let us replace in the first equation, so that Since we have to require ζ 2 > 0, we see that it must be ζ 2 1 − (p 2 + q 2 )ζ 1 + p 2 q 2 < 0. (E.65) This is equivalent to say q 2 < ζ 1 < p 2 . (E.66) So we can use τ = √ ζ 1 as a modulus to represent W . The moduli space, including the boundary, is therefore This is sufficient to show that the case when p = q must be excluded, since the solution is no more periodic.

Appendix F: The Baryonic number
The Baryon number is defined by the integral where we used the explicit expressions for the R a . After using (B.25), we get |c j | 2 ε j sin(ar) Tr(h ′ J j ), and using that −ε j Tr(h ′ J j ) = a, we finally get ǫ ijk Tr(R i R j R k ) = 6σm √ g c 2 a sin(ar).
Replacing in the integral and integrating we get Remark: the form is nothing but the pull back on the rectangular box of the volume form Tr(R ∧ R ∧ R) over the cycle, see for example [39].

Appendix G: Minimal energy per Baryon
Let us minimise expression (IV.53) w.r.t. the L a , a = ϕ, r, γ. Let us rewrite it in the form g(L ϕ , L r , L γ ) = DL ϕ L r L γ Deriving w.r.t. L j and setting we get the equations for the stationary points: A 2 x + B 2 y + C 2 xy − z(1 + α 2 x + β 2 y) = 0, (G.4) A 2 x − B 2 y + C 2 xy − z(1 − α 2 x + β 2 y) = 0, (G.5) Solving the first equation w.r.t. z and replacing in the remaining equations, we get z = A 2 x + B 2 y + C 2 xy 1 + α 2 x + β 2 y , (G.7) 0 = B 2 y + B 2 β 2 y 2 − α 2 x 2 (A 2 + C 2 y), (G.8) From the third equation we get which replaced in the second term of the second equation gives Since we are looking for positive x, y, z, the second factor is strictly positive and the only allowed solution is x = B αC . Replacing in (G.10) and then in (G.7), we get