Shape stability of pasta phases: Lasagna case

The stability of periodically placed slabs occurring in neutron stars (lasagna phase) is examined by exact geometrical methods for the first time. It appears that the slabs are stable against any shape perturbation modes for the whole range of volume fraction occupied by the slab. The calculations are done in the framework of the liquid drop model and obtained results are universal - they do not depend on model parameters like surface tension or charge density. The results shows that the transition to other pasta shapes requires crossing the finite energy barrier


I. INTRODUCTION
The very long structures appearing in neutron star matter, called pasta phases, are a commonly accepted phenomenon. Following the seminal work [1] the pasta phases have been studied in many different manners. Many approaches are based on the compressible liquid drop model (CLDM) where the system is described by two homogeneous phases separated by a sharp boundary with non-zero surface tension [2][3][4]6]. The competition between the Coulomb and surface energy leads to different shapes which are usually described in the Wigner-Seitz approximation, where the geometry of the phase is imposed at the beginning. The spherical and cylindrical cells are not the correct unit cells as they are not able to fill the whole space by periodic placement. Structures obeying periodicity, in the form of gyroid and diamondlike shapes, were introduced in [7,8]. However, they must be treated only as an approximation of a true shape, because they do not satisfy the necessary condition for the cell energy extremum. It was shown in [9] that the condition relates the mean curvature H(x) of the cluster surface and electrostatic potential Φ(x) where σ and ∆ρ are surface tension and charge difference between the cluster and its surroundings and C is a gauge dependent constant. It means the surface curvature depends on the potential distribution and is not constant in general. In fact, the phases considered in [7] represent surfaces with H = 0 and as such they cannot be the true solution of Eq. (1). In works [7,8] it was also shown that they have larger energy than the flat slabs. Nevertheless, the energy difference is not large, so it seems that such triple-periodic structures are likely. Recent analysis of pasta phases, including the periodicity, based on quantum or classical Molecular Dynamics [10][11][12][13][14] and Time-Dependent Hartree-Fock [15,16] or Hartree-Fock with twist-averaged boundary conditions [17] have shown the existence of triply-periodic in the * skubis@pk.edu.pl form of gyroid, diamond, twisted spaghetti and waffles. The general motivation of our work is to confirm the existence of such non-trivial periodic structures also in the framework of CLDM. Therefore, we propose to examine the possibility of transition from lasagna phase to other shapes by its shape deformation. The inspection of various perturbation modes could indicate to what kind of structures the lasagna is going to transform. Such consideration corresponds to the stability analysis of the lasagna phase. The answer to the stability question is not obvious. Though the contribution from surface energy for flat slab is always positive, the Coulomb energy is not and could destabilize the lasagna.
In the work [9] an overall discussion of the rigorous treatment of periodically placed proton clusters of any shape was presented. The stable cluster surface should satisfy not only the Eq. (1), which represents the necessary condition for the extremum, but also the condition for the minimum coming from the inspection of second variation of the energy functional.
The second order analysis, in limited sense, was already carried for cylinders and balls in the [5], where particular deformation mode was examined in the isolated Wigner-Seitz cell. Periodically placed cylinders and balls were considered in [6] and later in [19], but we need to be aware that these structures do not represent the proper minimum determined by the Eq. (1), their shape was assumed apriori. Up to now, the only known true periodic solution of Eq. (1) is lasagna (the boundaries of slabs coincide with the equipotential surfaces) and thus the stability analysis of these structure may be carried out.
Partial stability analysis of lasagna phase was done in by Pethick and Potekhin in the context of elastic properties of the phase in the work [19]. In order to determine the elasticity coefficient, they considered one particular deformation mode and moreover the expression for the energy was valid only in the limit of wavelength going to infinity. Such analysis corresponds to the stability consideration, however, being limited to the only one type of deformation. Here, we present general, unconstrained stability analysis for any kind of deformation with finite wavelength, what, to the authors' knowledge, has never been carried out.
FIG. 1. The unperturbed charge distribution and its potential Φ0.

II. ENERGY VARIATION FOR SINGLE SLAB
Let's consider proton cluster P in the shape of one slab placed in the center of unit cell with size a, b, c and volume V C = abc. The slab is perpendicular to x-axis and occupies w fraction of the total cell volume, 0 < w < 1. The charge density contrast between phases is ∆ρ = ρ + −ρ − , so the slab is positively charged with the density ρ + = (1 − w) ∆ρ and immersed in negatively charged neutron gas ρ − = −w ∆ρ (we assume the homogeneous electron background). For such charge distribution, the unperturbed potential Φ 0 , Fig. 1, is function of x only and takes the form Such configuration fulfill the necessary conditions for minimum, i.e. equations coming from vanishing variation of total energy in the first order. The sufficient condition for minimum is expressed by the positive value of total energy variation in the second order δ 2ε where the constraints for baryon number and charge conservation are imposed. Such energy variation for any deformation of the proton cluster is expressed by the integral over proton cluster surface ∂P [9] where B 2 = κ 2 1 + κ 2 2 = 0 is the sum of squared principal curvatures, ∂ n Φ 0 is the normal derivative of unperturbed potential, δ Φ is the first order perturbation of the potential caused by the surface deformation . One must remember that here we consider only the deformation which preserves the volume of the cluster, that means Deformations of this kind may be expressed in terms of Fourier series on the y, z-plane. The deformations on each slab face are independent, so we get two series for each face located at x 1 and x 2 . For further calculations it is convenient to introduce expansion based on complex amplitudes α j mk , where m and k are the mode indices and the superscript j corresponds to the face number where we introduce 3-dimensional discrete wave vector and x = (x, y, z). The m and k indices take integer values from −∞ to +∞ except the case when both of them equal zero, which is indicated by an apostrophe in the sum sign. The lack of (0, 0)-mode comes from the volume conservation condition, Eq. (3). Since the deformation function j (y, z) must be real it means that the complex amplitudes should fulfill the following relations It is more natural to introduce the cosine mk,C cos( my b + kz c ) and sine modes mk,S sin( my b + kz c ) in the expansion Eq. 4. Then the relation between complex and real amplitudes is In the Fig. 2 an example of the slab deformation with real amplitudes taking the following values 1 10,S = 2 10,C = 0 and 1 10,C = 2 10,S = 0 is shown. Below we present and discuss the subsequent contribution to the total energy variation δ 2ε , Eq. (2). In terms of complex amplitudes the surface energy contribution to the 2nd order variation takes simple form It is always positive, as B 2 = 0 for the slab and due to relations (6) we get α j mk α j −m,−k = |α j mk | > 0. That is an important fact, because in many cases the surface energy acts as a destabilizer, like, for example, δ 2ε S < 0 in the case of Rayleigh-Plateau instability.
The contribution coming from the electrostatic interaction includes two terms. The first one, determined by the normal derivative of potential Φ 0 , is given by The normal derivatives at the slab faces are ∂ n Φ 0 ( which is always negative. The second term of Coulomb interaction part is associated with perturbation of the potential δ Φ where δ Φ is the first order potential perturbation calculated thanks to the periodic Green function For the three-dimensional unit cell with sizes a, b, c the periodic Green function can be expressed as the sum over discrete modes numbered by three indices n, m, k The prime sign in the summation means that n, m, k cannot vanish simultaneously, similarly as in the Eq. (4). Such Green function fulfills the Poisson equation in the unit cell with periodic boundary conditions [20] (we use CGS units) and the absence of (0,0,0) mode in the expansion (14) is a consequence of the unit cell neutrality. By joining The slab deformation in the unit cell (the grey rectangle). The first face is only with sine mode whereas the second face with cosine mode.
Eqs. (4,13,14) we obtain the change in the energy with respect to the potential variation where ξ ij corresponds to the difference between the faces location and χ mk is the norm of the dimensionless wave vector The vector χ mk describes the manner in which the face is vibrating -the indices m and k determine the number of wavelengths being placed in the face sizes b and c. The function F (ξ, χ) comes from the summation over the n index. As the n-numbered modes do not depend on x, the summation over n can be done separately and defines a function F The above series can be evaluated [21] to a closed form .
The function is positive for any ξ and χ. In some special cases the function simplifies to Taking together all energy terms we get the total energy variation expressed in terms of mode amplitudes α mk As the δ 2ε is the quadratic form of the amplitudes α mk the competition between the terms in curly brackets decides about stability of the face surface. There are three characteristic terms: one from the surface energy being positive and two from Coulomb interactions: first is negative and the second is always positive. It seems the stability consideration depends on the values of surface tension σ or charge contrast ∆ρ but it appears that these parameters may be removed from our analysis. One of the conditions for the minimum of the energy for unit cell is the virial theorem [9]. The theorem takes the form of relation between the surface and Coulomb energy of the cell. For the cell with high symmetry, it has simple form from which we may get the relation between σ and ∆ρ Finally the total energy variation may be written as the quadratic form of α i with its coefficients A ij mk given by The obtained general expression, Eq. (26) for the second order variation of the total energy allows us to determine whether the slab is stable with respect to any deformations preserving its volume. It is worth noting that A ij mk coefficients, which decide about stability, does not depend on the strength of interactions being determined by surface tension σ and charge contrast ∆ρ. In this way, we obtained an interesting result that the stability of pasta depends only on the geometry of phase and mode under consideration and not on the details of strong or electromagnetic interactions.
Before general discussion of the stability of a single slab we show how the above results work in the stability analysis for particular class modes.

III. AN EXAMPLE OF STABILITY ANALYSIS
Let's consider the simplest surface perturbation consisting of combination of sine and cosine modes going along y axis for each face, which corresponds to m = ±1 and k = 0. Then, for the i-th face, the only non-vanishing complex amplitudes α i mk are: where we introduced the real amplitudes i S and i C corresponding the functions sin( 2πy b ) and cos( 2πy b ). The deformation i for i-th face is then a function of y only It is convenient to introduce the vector built of deformation amplitudes Then the total energy variation for such deformation may be written down in the matrix form where the dimensionless matrixM iŝ and its elements are: where w is the volume fraction and χ determines the wavelength of mode in comparison to the cell size χ = 2πa/b. The inspection of eigenvalues and eigenvectors of M allows for a complete stability analysis for the deformation we have chosen. The matrixM given by Eq. (33) possesses the two-fold degenerated two eigenvalues. Further, we call theM -eigenvalues λ l , l = 1 . . . 4 as stability function. For our concrete form ofM the eigenvalues and their eigenvectors are The stability functions λ i do not depend on the details of interactions σ, ∆ρ but only on the geometry of our system which is described by volume fraction w and mode wavelength χ. In Fig. 3 the stability functions are plotted in the w, χ parameter space for their eigenvectors. These vectors represent two classes of modes. We may call them as snaky and hourglass-shaped modes. The snaky modes occur when the deformations on both faces are in phase ( 1 , 3 ) and hourglass modes occur when the deformations on faces are in antiphase ( 2 , 4 ). As we may notice, for all values of volume fraction and mode wavelengths the stability functions are positive, that means, that for both classes of modes the slab is stable. However the stability is not too strong. Careful inspection of λ 1 and λ 2 shows that those functions go to 0 when χ and w approach to some specific values and λ 2 | χ→∞ → 0 for w → 0 or 1 . It means that snaky modes become unstable in the limit of very long waves regardless slab thickness, whereas the hourglass modes become unstable for very thin slab and very short waves. To sum up, we may say the slab becomes asymptotically unstable for very long mode or for very short mode when the slab becomes very thin in comparison to unit cell width. One should note that, the case when w ≈ 0 or w ≈ 1 must be treated with caution because in reality the cluster surface has finite thickness and for very thin slab the validity of liquid drop model could be questioned.

IV. ELASTIC PROPERTIES
Nuclear pastas share their elastic properties with liquid crystals. If the wavelength of the snaky mode is very large in comparison to the cell size, it corresponds to the so-called splay deformation of the liquid crystal. It allows to determine elastic constant K 1 for the lasagna phase what was shown by Pethick et al. in [19]. By definition, the constant K 1 , relates the deformation energy with the transverse derivative of the deformation field when the mode wavelength becomes very large. In our notation the relation takes the form where the bracket .. means the average over the slab surface. Taking the definition of K 1 in limit of the very long mode, we get (here denote the mode amplitude only). The deformation energy for the snaky mode is Applying the limit in Eq.(38) we get what exactly corresponds to the result of [19] if the ∆ρ is replaced by the unperturbed Coulomb energy of the cell ε C,0 = 1 6 πa 2 ∆ρ 2 (1 − w) 2 w 2 . In comparison to [19] our approach represents an improvement. Pethick et al. used a summation for the deformation energy which is in fact divergent. Then, the authors treated the divergent series as an asymptotic expansion with respect to χ −1 and finally they obtained correct result. In our approach we avoid any divergences.

V. GENERAL DISCUSSION OF STABILITY
In the Section III the stability analysis for the simplest slab deformation was carried out. The full stability analysis would require the inclusion of modes for all multiplicities m and k. Writing down the matrixM for all modes it appears to take the block-diagonal form where in mk-th position we get the same matrix as in Eq. (33) with A and B elements taking the same form like in Eq. (35) but with the replacement χ → χ km . Such a form ofM means that the modes with different multiplicity m, k and the cosine-like and sine-like modes do not couple. So, taking fixed m, k one may repeat the discussion from the Section III and finally conclude that the slab is stable for any mode keeping in mind the asymptotic cases, described by the Eqs. (36,37).

VI. CONCLUSIONS
In this work, by use of analytical methods, we have shown that proton clusters having the form of slabs placed periodically in space are stable for all values of volume fraction occupied by the cluster. It is a compelling result. All works based on CLDM (see for example [1,4,6,7], show that different shapes of pasta are preferred for different values of w. Here, we have shown that lasagna phase is stable in the whole range of w. One must remember that our analysis means the lasagna phase represents merely a local minimum. The transition to another geometry is not totally blocked, but requires finite size deformation in order to exceed the energy barrier. It may be interpreted that, at least for some range of volume fraction the lasagna phase is metastable. That could be quite interesting for the dynamics of pasta appearance during the neutron star formation. Our analysis is based on small deformations so it cannot state what range of w the lasagna represents global minimum of the cell energy. First, the global minimum statement requires the knowledge of exact solutions of Eq. (1) for other kind of shape than flat slab. So far, we have not known such solutions in the CLDM approach. The seeking of them mark out the direction of further research of pasta by differential geometry methods.
We are also conscious that the approach, based on the CLDM, has its limitations and inclusion of such effects like finite thickness of the cluster surface or the temperature fluctuations could change the final conclusion concerning lasagna phase stability.