Magnetar QPOs and neutron star crust elasticity

The crust region is a tiny fraction of neutron stars, but it has a variety of physical properties and plays an important role in astronomical observations. One of the properties characterizing the crust is the elasticity. In this review, with the approach of asteroseismology, we systematically examine neutron star oscillations excited by crust elasticity, adopting the Cowling approximation. In particular, by identifying the quasi-periodic oscillations observed in magnetar flares with the torsional oscillations, we make a constraint on the nuclear saturation parameters. In addition, we also discuss how the shear and interface modes depend on the neutron star properties. Once one detects an additional signal associated with neutron star oscillations, one can get a more severe constraint on the saturation parameters and/or neutron star properties, which must be a qualitatively different constraint obtained from terrestrial experiments, and help us to complementarily understand astrophysics and nuclear physics.


Introduction
Neutron stars produced through the core-collapse supernova at the last moment of the massive star's life can become in extreme states, which is quite difficult to realize on Earth.For example, their density easily exceeds the standard nuclear density, ∼ 0.16 fm −3 , and the gravitational and magnetic fields around/inside the star become much stronger than those observed in our solar system [1].Therefore, physics in such extreme states may inversely be revealed through the observations of neutron stars and/or their phenomena.
The stellar structure depends on the equation of state (EOS) for neutron star matter, which is under beta-equilibrium and charge neutrality.However, the determination of EOS for neutron star matter from the experiments is still difficult.This is because the density of nuclei is around the saturation density almost independently of nuclear species, due to the nature of nuclear saturation properties.Namely, the experimental data concentrates around/below the saturation density, while the neutron star properties are mainly determined from the EOS for a higher-density region.Anyway, the saturation parameters characterizing the EOS in a lower-density region are gradually constrained through terrestrial experiments, e.g., [2][3][4][5][6][7][8][9][10], which are strongly associated with a low-mass (or low central density) neutron star properties.
Meanwhile, it is better to note that experiments cannot directly measure the saturation parameters, where one has to estimate them from a kind of empirical relations (or strong correlations) between the saturation parameters and experimental data, suggested by theoretical studies, e.g., [11][12][13][14][15]. Thus, even if the experimental accuracy would be improved well, because of the theoretical uncertainties in these relations, it is not always true that the constraints on the saturation parameters become better [16][17][18].In practice, the estimation of the density-dependence of the nuclear symmetry energy, the so-called slope parameter L (see Eq. ( 1) for the definition), from the parity asymmetry of the polarized electron scattering cross section of 208 Pb strongly depends on the model [19,20].This means that the Several constraints on the neutron star mass and radius obtained from the astronomical observations and terrestrial experiments.MSP J0740+6620 is one of the massive neutron stars discovered up to now, whose mass is M/M ⊙ = 2.08 ± 0.07 [24], while the radius (and mass) of this object are also constrained via NICER as well as PSR J0030-0451.The tidal deformability restricted from the GW170817 event tells us the upper limit of the radius of 1.4M ⊙ neutron star, i.e., R 1.4 ≤ 13.6 km [26].The observation of the X-ray bursts from neutron stars also gives us the constraint on the neutron star mass and radius, although they depend on the theoretical model [40].In addition, one can theoretically exclude the top-left shaded region due to the causality [41].On the other hand, the terrestrial experiments, e.g., by PREX-II [19], by SπRIT [42], and at RCNP [43], give us the constraint on the stellar models constructed with a lower central density (bottom-right region), adopting the mass formula for low-mass neutron stars [44], where the constraints expected with the fiducial value of the nuclear saturation parameters, i.e., L = 60 ± 20 MeV [3,7,10] and K 0 = 240 ± 20 MeV [2], is also shown (see Eq. (1) for the definition of L and K 0 ).Furthermore, for reference, the mass and radius relations constructed with several EOSs, i.e., SLy4 [45,46], SKa [47], SkI3 [48], SkMp [49], DD2 [50], Shen [51], and Togashi [52], are plotted.Taken from [37].
theoretical study must be vitally significant for deriving a technique to directly estimate the saturation parameters from the experimental data, although this is out of scope in this review.
On the other hand, astronomical observations must be important for constraining the EOS in a higher-density region.For example, the discovery of the massive neutron stars, whose mass is ∼ 2M ⊙ , excluded the soft EOSs, with which the maximum mass does not reach the observed mass [21][22][23][24].The observation of gravitational waves from GW170817 [25] could restrict the dimensionless tidal deformability of neutron stars, which leads to the constraint that the 1.4M ⊙ neutron star radius should be less than 13.6 km [26].In addition, through the observation of the pulsar light curve, one may constrain the neutron star properties, especially the stellar compactness, M/R, with stellar mass M and radius R, e.g., [27][28][29][30][31][32].This is because the trajectory of the photon radiating from the neutron star's surface can bend due to the strong gravitational field induced by the neutron star, which is a relativistic effect.In practice, the Neutron star Interior Composition Explorer (NICER) is now operating on an International Space Station, which successfully constrained the neutron star mass and radius for PSR J0030-0451 [33,34] and for PSR J0740-6620 [35,36].
In such a way, the neutron star mass and radius would be constrained more and more in a higher density region with astronomical observations and in a lower density region with the terrestrial experiments, e.g., [37][38][39] (also see Fig. 1).Then, as these constraints from observations and experiments become more and more severe, one may eventually constrain the EOS for neutron star matter.
In addition to the observations (or estimations) of the neutron star mass and radius, the oscillation signals (and also the gravitational waves) from neutron stars, if observed, are another important information to extract the neutron star properties.Since the objects have their own specific oscillation frequencies, one may know the interior information via the observation of their frequencies as an inverse problem.This technique is known as asteroseismology, which is similar to seismology on Earth and helioseismology on the Sun.In fact, the neutron star has a lot of oscillation modes, where each mode can be excited due to the corresponding physics [53].For example, it is well-known that the frequency of fundamental oscillations in a neutron star is strongly associated with the square root of the stellar average density, because the fundamental oscillations are acoustic oscillations.The gravitational waves from the neutron stars may be suitable astronomical information to adopt asteroseismology, but still, they have never been observed from an isolated neutron star.Even so, it has already been suggested that neutron star mass, radius, and EOS by observing the gravitational waves from the (cold) neutron stars, e.g., [54][55][56][57][58][59][60][61][62][63][64].Moreover, this technique is also adopted for extracting the information of protoneutron stars from the supernova gravitational waves, e.g., [65][66][67][68][69][70][71].
Instead of the gravitational waves, the quasi-periodic oscillations (QPOs) discovered in the afterglow following the magnetar giant flares [72][73][74][75][76][77][78], are also valuable information for adopting asteroseismology, where magnetars are strongly magnetized neutron stars.The observed QPO frequencies are in the range of tens of Hz up to kHz (see Sec. 2 for details).Although the emission mechanism of the flare activity in magnetar is still under debate, the observed QPOs are considered to be strongly associated with neutron star oscillations.Taking into account the dynamical time of neutron stars, which is around 0.1 msec, it may be more difficult to theoretically explain the lower frequencies among the observed frequencies.The possible candidates may be either crustal torsional oscillations, magnetic oscillations, or magneto-elastic oscillations, if one assumes that the observed QPO frequencies come from the neutron star oscillations.However, the magnetic oscillations strongly depend on the magnetic field strength and its geometry inside the star, e.g., [79][80][81], while their understanding is quite poor.Thus, here we discuss the crustal properties by identifying the observed QPO frequencies with the crustal torsional oscillations, neglecting magnetic effects (see Sec. 2 for the magnetic effects on the crustal torsional oscillations).In this way, the crustal properties have been constrained indeed, e.g., [82][83][84][85][86][87][88].
Through the approach with asteroseismology, it is shown that the observed QPO frequencies can be identified with the neutron star oscillations.On the other hand, there are still open issues in the magnetar QPOs.For instance, all the observed QPOs are not simultaneously excited after the magnetar flares, where some of the QPOs are excited first, then others.This phenomenon cannot be explained within the linear analysis.To understand the time dependence of the QPO frequency excitation, one may have to discuss beyond a linear analysis, such as a nonlinear coupling.Moreover, to discuss the crustal oscillations, one has to estimate the crustal elasticity, which is characterized by shear modulus.The details are shown in Sec. 4.
The oscillation frequencies in a spherically symmetric neutron star model can be classified into two families with parity, i.e., axial and polar type oscillations, where the axial type oscillations can be excited independently of the polar type oscillations.The crustal torsional oscillations belong to the axial type oscillations.Since the axial type oscillations do not involve the density variation, they are relatively easily excited but less important in the gravitational wave observations.On the other hand, polar-type oscillations involve density variation and stellar deformations, where various oscillation modes can exist.
Owing to the crust elasticity, the shear and interface modes, which are polar type oscillations, are also excited as well as the torsional oscillations, even though they are not so discussed well compared to the torsional oscillations.The interface modes, whose frequencies are less than 100 Hz, may be excited in binary neutron stars by resonating with the orbital frequency.In fact, the precursors have been observed just before the main flare activity of gamma-ray burst at a binary neutron star merger [89], which may be a result of the resonant shattering of neutron star crusts induced by the binary orbital motion [90,91].If so, one could extract the neutron star properties by carefully observing the precursors and identifying them with neutron star oscillations [92][93][94][95].Anyway, to extract the properties from the resonant shattering of neutron star crust, one has to examine the shear and interface modes systematically.
In this review, mainly focusing on the torsional oscillations excited in the neutron star, we extract the neutron star properties by identifying the observed QPO frequencies with the torsional oscillations.For this purpose, we simply assume the Cowling approximation in this review, i.e., we only consider the fluid oscillations with the (unperturbed) background metric.The perturbation equations are derived from the linearized energy-momentum conservation law.By imposing appropriate boundary conditions on the numerical boundary, the problem to solve becomes an eigenvalue problem.The neutron star crust is mainly composed of spherical nuclei, but the non-spherical nuclei may appear in the vicinity of the boundary of crust and core, depending on the nuclear parameters (or EOS).So, first, we start to discuss the observed QPO frequencies with the crustal torsional oscillations excited in the region only composed of spherical nuclei.Then, the discussion is extended to the more realistic situation with the non-spherical nuclei.Through these attempts, we successfully identify all the observed QPOs with the crustal torsional oscillations and finally derive the constraint on the slope parameter as L ≃ 58 − 73 MeV.Furthermore, we also discuss the identification of the high-frequency QPOs observed in GRB 200415A with the overtones of crustal torsional oscillations.Through this identification, we make a constraint on the neutron star mass and radius for GRB 200415A, adopting the fiducial values of nuclear saturation parameters.In addition, we systematically examine the shear and interface modes for various neutron star models, with which we find the relation between the frequencies and the neutron star properties independently of the stiffness of EOS for core region.
This manuscript is organized as follows.In Sec. 2, we shortly mention the observed QPO frequencies in magnetar giant flares, which will be considered as evidence of neutron star oscillations.In Sec. 3, we show the equilibrium neutron star models and the EOS considered in this study.In Sec. 4, we briefly mention the shear modulus inside the neutron star crust, which is an important property for considering the crustal oscillations.In Sec. 5, we extract the crustal information in practice by identifying the observed QPOs with the crustal torsional oscillations.In Sec. 6, we also show the behavior of the shear and interface modes, especially focusing on the association with the neutron star properties.Finally, we conclude in Sec. 7. Unless otherwise mentioned, we adopt geometric units in the following, c = G = 1, where c and G denote the speed of light and the gravitational constant, respectively.

Magnetar QPOs and magneitc effects on the crustal oscillations
Compared to usual pulsars whose surface magnetic field strength is ∼ 10 12 − 10 13 G, the existence of neutron stars with strong magnetic fields, such as ∼ 10 14 − 10 15 G, is observationally known.This strongly magnetized neutron star is known as a magnetar.In a rotating neutron star, the magnetic stress induced by the stellar rotation gradually accumulates in the crust, which is usually supported by crustal elasticity.But, once the magnetic stress becomes significantly strong and the crust elasticity cannot support it, the crust eventually breaks out, and the magnetic energy is released through the flare activities [96,97], where star quakes may also happen.This scenario could be the origin of the QPOs observed in the afterglow following the magnetar giant flares.
Assuming that these QPOs come from neutron star oscillations, one may identify the higher QPO frequencies with polar-type oscillations of neutron stars, e.g., Fig. 25, but the lower frequencies are more difficult to identify with neutron star oscillations.Maybe the possible candidates are only the crustal torsional oscillations or magnetic oscillations (or magneto-elastic oscillations).Meanwhile, the magnetic oscillations definitely depend on the geometry and strength of magnetic fields, which have not been understood well yet.In addition, it is pointed out that crustal torsional oscillations can be excited in the vicinity of the neutron star surface, while magnetic oscillations are confined only to the core region, if magnetic fields are not so strong.This is because the shear velocity is higher than the Alfvén velocity, v A ≡ B/ 4πρ, with the magnetic field strength, B, and density, ρ, at the basis of the crust, which leads to that the oscillations are controlled by shear properties.In general, the critical field strength, above which the magnetic oscillations become dominant, depends on the EOS (or crust properties), which is considered around a few times 10 14 G [79,81,103].So, to avoid the uncertainty in the magnetic geometry and profile, assuming that the field strength is not so strong, we focus on only the crustal oscillations without magnetic effects in this review.
As an advantage without the magnetic effect, the crustal torsional oscillations are completely confined inside the crust, since the axial-type oscillations do not involve radial variations.That is, one can avoid the large uncertainties in the EOSs of the core region, and directly discuss the crust properties by comparing the observed frequencies to the crustal torsional oscillations.On the other hand, if one takes into account the magnetic effects, even torsional oscillations can be coupled with the core region through the magnetic fields, where oscillations easily damped out via magnetic coupling if the magnetic fields are relatively larger [104,105].We note that one has to consider the polar-type oscillations in the whole region of a neutron star, even if the magnetic effects are neglected, because the polar-type oscillations involve radial variations.
We only focus on the neutron star oscillations associated with the crust elasticity in this review, while we briefly mention the magnetic oscillations in the neutron star here.Since the existence of magnetic fields inside the star breaks a spherically symmetric configuration, one cannot examine the magnetic oscillations (and also magneto-elastic oscillations) in the same way as the case on the spherically symmetric objects (as mentioned in the following sections).Assuming a specific distribution of magnetic fields, one may study them by decomposing the perturbative variables with spherical harmonic function, where ℓ-th order perturbation equations are generally coupled with ℓ-th order equations of other values [79,106].Or, one may solve a time evolution of two-dimensional perturbation equations and determine the frequencies via FFT.
With the later approach, it has been shown that the axial-type magnetic oscillations are continuous spectrum instead of the discrete one [107][108][109].This is because the frequencies of magnetic oscillations are characterized by the propagation time determined with the magnetic field length inside the star and Alfvén velocity, while the magnetic field length inside the star continuously changes from the axis to the position related to the closed field line close to the equatorial plane [108].This phenomenon may also be understood with a toy model suggested in [110,111].In any case, taking into account the coupling with the crust, i.e., magneto-elastic oscillations, one may identify the observed QPO frequencies [80,105,112,113].Furthermore, one may constrain the magnetic geometry from the QPO frequencies and estimated field strength of magnetic fields, i.e., the case that the magnetic fields are confined in the crust due to type I superconductor in the core region can be excluded [114].
Compared to the extensive studies about axial-type magnetic oscillations, the study of polar-type magnetic oscillations is relatively poor.Even so, unlike the axial-type oscillations, the polar-type magnetic oscillations seem to be discrete spectrum [115], and they definitely play an important role in the gravitational wave radiations [116,117].1).The bottom panel is an enlarged view of the above panel.Taken from [124].

Crustal equilibrium
Although the details of a neutron star structure depend on the EOS for neutron star matter, its conceptual structure is generally accepted.Under the thin envelope of the atmosphere and/or ocean, the matter forms a Coulomb lattice, which behaves as a solid.This region corresponds to the neutron star crust.Most of the crust region is composed of spherical nuclei with a body-centered cubic (bcc) lattice, while the shape of stable nuclei changes into cylindrical, slablike, cylindrical-hole, and spherical-hole, and eventually becomes uniform matter, as density increases [118,119].The region composed of nonspherical nuclei is the so-called pasta phase, which behaves as a liquid crystal.Meanwhile, the region whose density is higher than the transition density from the non-uniform to the uniform matter is the neutron star core.This transition density depends on the neutron star EOS, which becomes ∼ (1/3 − 1/2) of the nuclear saturation density.The transition density where the nuclear shape changes inside the crust also depends on the EOS, especially on the density dependence of the symmetry energy, L, in Eq. (1) [120][121][122].In addition, the crust thickness is at most 10% of the stellar radius, which depends on the stellar compactness, M/R, and the nuclear saturation parameters [119,123] (see Fig. 2 for an example of the concrete stellar model).In general, the ratio of crust thickness to the stellar radius decreases as the stellar compactness and L increases.In particular, for L ≳ 100 MeV, the pasta phase disappears, where the crust is composed of only the spherical nuclei, and directly transitions into uniform matter [120][121][122] (also see Table 1).
The EOS for neutron star matter generally depends on the adopted nuclear interaction, nuclear models, and the compositions.But for any EOSs, the bulk energy per baryon for the zero-temperature nuclear matter can be expanded in the vicinity of the saturation density, n 0 ≃ 0.16 fm −3 , of symmetric nuclear matter as a function of the baryon number density, n b , and neutron excess, α, via [125] where u ≡ (n b − n 0 )/(3n 0 ), while n b and α are given by n b = n n + n p and α = (n n − n p )/n b , using the neutron and proton number density, n n and n p .We note that the matter Table 1.The EOS parameters for the OI-EOS family adopted in this review, where y is defined as y ≡ −K 0 S 0 /(3n 0 L).SP-C, C-S, S-CH, CH-SH, and SH-U denote the transition densities for each EOS parameter set.The asterisk at the value of K 0 denotes the EOS model with which some pasta phases are not predicted to appear, i.e., the values with * 1, * 2, and * 3 denote the transition densities from cylindrical nuclei to uniform matter, from cylindrical-hole nuclei to uniform matter, and from spherical nuclei to uniform matter, respectively.
K 0 (MeV) −y (MeV fm 3 ) L (MeV) SP-C (fm −3 ) C-S (fm −3 ) S-CH (fm −3 ) CH-SH (fm −3 ) SH-U (fm with α = 0 corresponds to a symmetric nuclear matter, i.e. n n = n p , while that with α = 1 is a pure neutron matter, i.e., n b = n n .The coefficients in this expansion are the nuclear saturation parameters, which characterize each EOS, and the term inside the brackets in front of α 2 corresponds to the nuclear symmetry energy.In particular, w 0 and K 0 are the binding energy and incompressibility of symmetric nuclear matter at n b = n 0 , while S 0 is the symmetry energy at n b = n 0 .
Owing to the nature of the nuclear saturation properties, n 0 , w 0 , and S 0 are well determined via experiments, i.e., n 0 ≈ 0.15 − 0.16 fm −3 , w 0 ≈ −15.8 MeV [9], and S 0 ≈ 31.6 ± 2.7 MeV [10].In contrast, K 0 and L are more difficult to determine experimentally, because they are a density derivative at n b = n 0 , i.e., one needs to know the nuclear information in a somewhat wide density range around the saturation point.Nevertheless, the constraints on these two parameters are improved, where their fiducial values are K 0 = 230 ± 40 MeV [2,4] and L ≃ 60 ± 20 MeV [3,7,10].We note that the constraint on L still seems to have large uncertainties, e.g., L = 106 ± 37 MeV with PREX-II done by the Thomas Jefferson National Accelerator Facility in Virginia [19], or 42 ≤ L ≤ 117 MeV with SπRIT by Radioactive Isotope Beam Factory at RIKEN in Japan [42].
Anyway, in this review, to systematically discuss the dependence of the crustal oscillations on the nuclear saturation parameters, we especially adopt a phenomenological EOS proposed by Oyamatsu and Iida [120,121] (hereafter referred to as OI-EOS).For given values of K 0 and L, the OI-EOS family is constructed by optimizing the values of n 0 , w 0 , and S 0 to reproduce the experimental data for masses and charge radii of stable nuclei in such a way that the bulk energy reduces to Eq. ( 1) in the limit of n b → n 0 and α → 0, adopting the extended Thomas-Fermi theory.The EOS parameters adopted in this review are listed in Table 1.We note that some of the EOS parameter sets are out of the range of fiducial values for K 0 and L, but we adopt such a wide range to systematically discuss the dependence on the saturation parameters.
As an equilibrium crust model (or a neutron star model), we simply consider a nonrotating, strain-free, and spherically symmetric neutron star model in this review.The metric describing such a stellar model is given by where Φ and Λ are the metric functions depending on only the radial coordinate, r.Λ is directly associated with the mass function, m(r), which corresponds to the enclosed gravi-tational mass inside the position r, i.e., e −2Λ = 1 − 2m/r.Adopting an appropriate EOS, e.g., the OI-EOS family in this review, one can construct the stellar model by integrating the Tolman-Oppenheimer-Volkoff (TOV) equation.In general, to construct a stellar model, one integrates the TOV equation outward from the center to the surface, assuming a central density, i.e., the neutron star models become one parameter family for each EOS.In fact, we construct the stellar model in this way for discussing the shear and interface modes in Sec.6, because these modes are polar-type oscillations, which are excited not only in the crust region but also inside the core.On the other hand, the torsional oscillations discussed in Sec. 5, which are axial-type oscillations, are confined inside the crust (elastic) region.So, to remove the uncertainties of the core EOS, we construct only the crust equilibrium for discussing the torsional oscillations, where we integrate the TOV equations inward from the stellar surface to the base of the crust, assuming the stellar mass and radius as a boundary condition at the stellar surface.

Shear modulus
The elasticity inside the crust can be characterized by the shear modulus, µ.The shear modulus, µ sp , in a phase of spherical nuclei with a bcc lattice is approximately described in the limit of zero-temperature as where n i , Z, and a are the ion number density, charge number of the ion, and Wigner-Seitz cell radius, i.e., 4πa 3 /3 = 1/n i [126].This expression has been derived by averaging the overall direction, assuming that each nuclei is only a point-particle [127].The shear modulus in the phase of spherical nuclei has been discussed more by including the phonon contribution [128], the electron-screening effect [129][130][131][132], the effect of polycrystal properties [133], and effects of finite size of atomic nuclei [134], which can modify the shear modulus more or less.But, in this review, we simply adopt µ sp given by Eq. ( 3) as shear modulus in the phase of spherical nuclei.
Compared to the shear modulus in the phase of spherical nuclei, the discussion for the shear modulus in the pasta phases is very poor.Nevertheless, according to Ref. [135], the shear modulus, µ cy , in the phase of cylindrical nuclei can be estimated as where E Coul and w 2 denote the Coulomb energy per volume of a Wigner-Seitz cell and the volume fraction of cylindrical nuclei, while the shear modulus, µ sl , in the phase of slablike nuclei is i.e., the matter composed of slablike nuclei behaves as a fluid within the linear response.This is because the deformation energy in slablike nuclei becomes of higher order with respect to the displacement [135].So, in this review we simply adopt Eqs. ( 4) and ( 5) as the shear modulus in the phases of cylindrical and slablike nuclei.But, on the other hand, it has also been suggested that the elastic constant in the poly-crystalline lasagna (which corresponds to slablike nuclei) may become a (nonzero) tiny value [136,137].If one considers the nonzero shear modulus inside the slablike nuclei, the results shown in this review may be changed.Furthermore, the shear modulus in the phase of cylindrical-hole, µ ch , and sphericalhole nuclei, µ sh , may be considered in a similar way to that in the phase of cylindrical and spherical nuclei.This is because the liquid crystalline structures in the phase of cylindricalhole and spherical-hole nuclei are basically the same as in the phase of cylindrical and spherical nuclei.Thus, we adopt Eqs. ( 4) and (3) for µ ch and µ sh by appropriately replacing the quantities in the formulae [88,138].That is, w 2 in Eq. ( 4) is replaced by the volume  2, where Sp, Cy, CH, and SH respectively denote the phases composed of spherical, cylindrical, cylindrical-hole, and spherical-hole nuclei.The right panel is an enlarged view of the shaded region in the left panel.Taken from [124].fraction of gas of dripped neutron for µ ch , while n i and Z in Eq. ( 3) are replaced by the number density of spherical-holes (bubbles) and the effective charge number, Z sh , of spherical-hole for µ sh , where Z sh is estimated as Z sh = n Q V sh with the effective charge number density of the spherical-hole, n Q , and the volume of the spherical-hole, V sh .n Q is given by the difference of the charge number density inside the spherical-hole from that outside the spherical-hole, i.e., n Q = −n e − (n p − n e ) = −n p , using the proton number density outside the spherical-hole, n p , and the number density of uniform electron gas, n e .
In Fig. 3, we show an example of the distribution of shear modulus inside the crust for the stellar model shown in Fig. 2, where the right panel is an enlarged view of the shaded region in the left panel.

Crustal torsional oscillations
On the crustal equilibrium models mentioned in the previous section, we make a linear perturbation analysis to determine the specific frequencies of crustal torsional oscillations.In this review, we examine the oscillation frequencies within the relativistic framework, while it has also been discussed within the Newtonian framework, e.g., [126,139,140].To derive the perturbation equations, we simply adopt the Cowling approximation, where the metric perturbation is neglected during the fluid oscillations.We note that one can expect to accurately determine the frequencies of the axial-type oscillations even with the Cowling approximation, because the axial-type oscillations do not involve the density variations, i.e., the metric perturbations (corresponding to the perturbations of gravitational potential) should be too small.The perturbation equation can be derived from the linearized energy-momentum conservation law.In practice, introducing the ϕ-component of the Lagrangian displacement, ξ ϕ = Y (r)e iωt ∂ θ P ℓ (cos θ)/ sin θ, where P ℓ (cos θ) denotes the ℓ-th order Legendre polynomial, the perturbation equation becomes [79,141] where the dash denotes the derivative with r and H denotes the effective enthalpy contributing the oscillations, given by for the spherical or cylindrical nuclei [85,142], while for the cylindrical-hole or spherical-hole nuclei [88].
In the expression of the effective enthalpy, H denotes the local enthalpy given by H = ε + p with the energy density, ε, and pressure, p, while A is the baryon number density inside a Wingner-Seitz cell, N s is the number of neutrons inside a Wingner-Seitz cell that do not comove with protons in the nuclei, N i is the number of neutrons inside a cylindrical-hole or spherical-hole, and R is a parameter characterizing a participant ratio in the oscillations, i.e., how much ratio of nucleons outside the cylindrical-hole or spherical-hole comove with it non-dissipatively.That is, all the nucleons inside a Wingner-Seitz cell contribute to the effective enthalpy with R = 1 (maximum enthalpy), while no nucleons outside the cylindrical-hole or spherical-hole do so with R = 0 (minimum enthalpy).We note that R in a phase of spherical-hole nuclei is predicted as ∼ 0.34 − 0.38 at n b = 0.08 fm −3 from the band calculations [143].Meanwhile, assuming that N s comes only from a part of the dripped neutrons inside a Wingner-Seitz cell, N d , the parameter N s /N d can control the fraction of neutrons contributing to the oscillations, i.e., N s /N d = 0 corresponds to the situation that all the dripped neutrons comove with the protons (maximum enthalpy), while N s /N d = 1 is that all the dripped neutrons behave as a superfluid and do not contribute to the oscillations (minimum enthalpy).
To determine the frequencies, one has to impose appropriate boundary conditions in the perturbation equation ( 6) at the basis of the crust (the boundary between the core and crust) and the surface of the crust (or the boundary between the crust and envelope).Both boundaries are essentially equivalent to the boundary between the elastic and fluid region, where we impose the condition that the traction force vanishes [79,141], which reduces to Y ′ = 0. On the other hand, at the interface between the elastic regions with different shapes of nuclei, e.g., between the spherical and cylindrical nuclei, we impose a continuous traction condition, i.e., µY ′ should be continuous at such an interface.Then, the problem to solve becomes an eigenvalue problem with respect to the eigenvalue, ω.The resultant value of ω gives us the frequency, f , via f = ω/2π for each value of ℓ.In this review, we use the notation of n t ℓ for expressing the eigenfrequencies of torsional oscillations with the angular index ℓ and the nodal number in the eigenfunction (although the notation is different from the usual one, such as ℓ t n ).

Torsional oscillations excited in the spherical nuclei
Even though the crust is not only composed of the spherical nuclei as mentioned above, to see the dependence of the frequencies of the torsional oscillations on the EOS parameters and stellar properties, first we consider the torsional oscillations excited in the elastic region composed of only the spherical nuclei.As a first step, we examine the frequencies without neutron superfluidity, i.e., N s /N d = 0 (or H = H).As mentioned in Sec. 3, to construct the equilibrium model, we have to select two parameters for EOS, i.e., K 0 and L, and two for stellar models, M and R. In the left panel of Fig. 4, we show the fundamental frequencies of ℓ = 2 torsional oscillations for a neutron star model with 1.4M ⊙ and 12 km, using the EOSs with some parameter sets listed in Table 1 (see Ref. [86] for the concrete parameter sets adopted here), as a function of L. From this figure, one can observe that the frequencies are less sensitive to K 0 .Thus, we will discuss the dependence of the fundamental frequency as a function of L through the fitting formula by assuming the polynomial form as where c (i) 2 with i = 0, 1, 2 are the adjustable parameters depending on M and R (and N s /N d ), while L 100 is the value of L normalized by 100 MeV.Using this fitting formula, one can estimate the frequency of 0 t 2 less than ∼ 5% accuracy for a neutron star model with 1.4M ⊙ and 12 km (Table 2 in [86]).In the left panel of Fig. 4, we also show the expected frequencies with the fitting formula given by Eq. ( 9) by the thick-solid line.
One may be able to understand the dependence of 0 t 2 on L as follows.As the value of L increases, the nuclear symmetry energy conversely decreases at the sub-nuclear density, considering the density dependence of the bulk energy given by Eq. ( 1).As a result, In the left panel, the fundamental frequencies of ℓ = 2 torsional oscillations for a neutron star model with 1.4M ⊙ and 12 km constructed with several EOS parameter sets are shown as a function of L. The thick-solid line is the expected frequency using the fitting formula given by Eq. ( 9).In the right panel, the expected fundamental frequency of ℓ = 2 torsional oscillations as a function of L, assuming that the neutron star mass and radius are in the range of 1.4 ≤ M/M ⊙ ≤ 1.8 and 10 km ≤ R ≤ 14 km, where for reference the lowest QPO frequency observed in SGR 1806-20, i.e., 18 Hz, is also shown.The vertical solid and broken lines correspond to L = 47.4 and 76.2 MeV, respectively.Taken from [86].
protons easily turn into neutrons, which leads to the situation that the charge number of nuclei gets smaller.This means that µ sp becomes smaller from Eq. ( 3).On the other hand, the fundamental frequency of the ℓ-th torsional oscillations, 0 t ℓ , are estimated as [139].Thus, one can expect that 0 t ℓ decreases with L. Now, assuming the neutron star mass and radius would be in the range of 1.4 ≤ M/M ⊙ ≤ 1.8 and 10 km ≤ R ≤ 14 km, which are reasonable assumptions as a neutron star model, the resultant fundamental frequency of the ℓ = 2 torsional oscillations is shifted, depending on the stellar models.In the right panel of Fig. 4, we show the expected region of 0 t 2 for such neutron star models as a function of L. We note that 0 t 2 decreases with R and M. So, the lower bound of the shaded region corresponds to the stellar model with 1.8M ⊙ and 14 km.On the other hand, the lowest QPO frequency observed in the SGR 1806-20, i.e., 18 Hz, is also shown in the right panel of Fig. 4. Since the ℓ = 2 fundamental torsional oscillations are theoretically the lowest frequency among a lot of eigenfrequencies of the torsional oscillations, 0 t 2 should be equal to or even lower than the observed lowest QPO frequency on the assumption that the magnetar QPOs come from the crustal torsional oscillations.Then, from the right panel of Fig. 4, one can get the constraint on L as L ≳ 47. 4 MeV, if the central object of SGR 1806-20 is a neutron star with M ≤ 1.8M ⊙ and R ≤ 14 km [84,86].Similarly, if the central object is a typical neutron star with 1.4M ⊙ and 10 km, L is constrained as L ≃ 76.2 MeV.
Next, we take into account neutron superfluidity in the crustal torsional oscillations, i.e., N s /N d ̸ = 0. To see how the frequencies depend on N s /N d , we calculate 0 t 2 , assuming that N s /N d is constant inside the inner crust.We confirm that 0 t 2 for a given neutron star model can be well expressed as a function of L using Eq. ( 9) even for N s /N d ̸ = 0.In Fig. 5, we show the L dependence of 0 t 2 for a neutron star model with 1.8M ⊙ and 14 km as varying the ratio of N s /N d with each 0.2 from 0 to 1 (the solid lines from bottom to top).The tendency of why the frequency increases with N s /N d can also be understood with the dependence of frequency on the shear velocity.As mentioned above, the effective enthalpy, H, decreases with N s /N d , which leads to the increase of the shear velocity.So, as a result, one can expect that the frequency also increases since 0 t ℓ ∼ v s .We note that the neutron star model discussed for the L dependence of 0 t 2 in Fig. 5 corresponds to the stellar model, which gives the lower boundary of the shaded region in Fig. 4.That is, assuming that the central object in SGR 1806-20 is a neutron star with M ≤ 1.8M ⊙ and R ≤ 14 km, and also assuming that the lowest QPO frequency observed in SGR 1806-20 comes from For a neutron star model with 1.8M ⊙ and 14 km, the fundamental frequency of ℓ = 2 torsional oscillations is plotted as a function of L, as changing the ratio of N s /N d with each 0.2 from 0 to 1, where N s /N d is assumed to be constant inside the inner crust.Meanwhile, the marks denote the frequencies, using the density-dependent N s /N d derived in [143], and the dashed line is their fitting with Eq. ( 9).Taken from [85].
0 t 2 , one can get the constraint on L as L ≳ L min with the lower limit of L min , depending on N s /N d , where L min is determined as the intersection of the L dependence of 0 t 2 and 18 Hz (dot-dashed line) in Fig. 5, i.e., L min increases with N s /N d [85].
In a realistic stellar model, N s /N d should depend on the density.According to the band calculations in [143], the ratio of N s /N d is only of the order of 10 − 30 percent at n b ∼ 0.01 − 0.4n 0 , even though the ratio of N s /N d is still under debate, e.g., [144][145][146].The fundamental frequencies determined with this density dependence of N s /N d are also shown in Fig. 5 with marks (and the dashed line for their fitting), which are close to the results obtained with N s /N d = 0.2.Anyway, since the frequencies of torsional oscillations strongly depend on the ratio of N s /N d , it is quite important to understand the ratio of N s /N d in a realistic situation.
In a similar way to the fundamental frequencies of the ℓ = 2 torsional oscillations, one can discuss the L dependence of the ℓ-th torsional oscillations.Actually, if one systematically examines the frequencies, one can find that 0 t ℓ is also well fitted as a function of L, which weakly depends on K 0 , as ℓ with i = 0, 1, 2 are also the adjustable parameters depending on M and R (and N s /N d ) [85].Using this L dependence of 0 t ℓ , we will discuss the QPO frequencies observed in SGR 1806-20 and 1900+14, even though one should also consider the effect of the existence of pasta (see Sec. 5.2 and 5.3, and also [147,148]).In Fig. 6, we show the fundamental frequencies of torsional oscillations with specific values of ℓ for a neutron star model with 1.4M ⊙ and 12 km, adopting the results in [143] for N s /N d .From the left panel one can observe that the 18, 26, 30, and 92.5Hz QPOs in SGR 1806-20 are well identified with ℓ = 3, 4, 5, 15 fundamental oscillations if L ≃ 128.0 MeV, while from the right panel one can observe that the 28, 54, and 84 Hz QPOs in SGR 1900+14 are with ℓ = 4, 8, 13 fundamental oscillations if L ≃ 113.5 MeV.We note that the other low-lying QPO frequencies discovered later can also be identified with the torsional oscillations with different values of ℓ in the same framework [88,131,142] (also see Sec. 5.2 and 5.3).
Since the neutron star mass and radius in SGR 1806-20 and SGR 1900+14 are not constrained, it may be better to consider the stellar models with a typical range of mass and radius.In practice, if the stellar mass and radius are changed in a certain range, the resultant L dependence of 0 t ℓ are shifted up and down in Fig. 6 with the same combination of the ℓ-th oscillations.As a result, the optimal value of L, with which the QPOs observed in SGRs are well identified with the torsional oscillations, is also shifted left and right.In the left panel of Fig. 7, the optimal values of L, with which the QPO frequencies observed in the value of L should be universal, i.e., independent of the astronomical events, one has to simultaneously explain both events, SGR 1806-20 and SGR 1900+14, with the same value of L. Thus, we get a more stringent constraint on L, i.e., 101.1 ≲ L ≲ 131.0 MeV [86], which are shown in the left panel of Fig. 7 with the shaded region.Furthermore, since the symmetry energy at n b = n 0 , S 0 , is approximately estimated as a function of L [120] through one can constrain S 0 to explain the QPOs observed in the SGRs as 35.6 ≲ S 0 ≲ 37.8 MeV.
On the other hand, considering the constraints on L obtained from terrestrial experiments whose fiducial value is L ≃ 60 ± 20 MeV [3,7,10], the resultant constraint on L may be too large, although a larger value of L has also been reported [19,42].So, if any, it may be better to consider an alternative possible correspondence to explain the QPO frequencies observed in the SGRs with the fundamental torsional oscillations.As an alternative possible  [86].This constraint on L additionally gives us the constraint on S 0 as 32.4 ≲ S 0 ≲ 34.4 MeV, using Eq.(11).We note that the fundamental frequencies of ℓ-th torsional oscillations excited inside the phase of spherical nuclei can be expressed as a function of ℓ, L, stellar mass, and radius [149].

Torsional oscillations excited in the spherical and cylindrical nuclei
Since in this review we simply consider that the shear modulus inside the phase of slablike nuclei is zero, as mentioned in Sec. 4, the torsional oscillations in the region composed of spherical and cylindrical nuclei can be excited independently from those in the region composed of cylindrical-hole and spherical-hole nuclei.That is, torsional oscillations are excited independently in two layers across the phase of slablike (lasagna) nuclei as "a lasagna sandwich" [88].In this subsection, we discuss the torsional oscillations excited in the outer layer, i.e., the phase of spherical and cylindrical nuclei, while those excited in the inner layer, i.e., the phase of cylindrical-hole and spherical-hole nuclei, will be discussed in the next subsection (Sec.5.3).
Unlike the phase of spherical nuclei, the understanding of the ratio of N s /N d in the phase of cylindrical nuclei is poor.So, we simply consider the extreme case, i.e., N s /N d = 0 and 1 in the phase of cylindrical nuclei.On the other hand, we adopt the result obtained in [143] for N s /N d in the phase of spherical nuclei, as in the previous subsection.Then, one can determine the frequencies of torsional oscillations by solving an eigenvalue problem.
To see how the frequency of the torsional oscillation changes due to the existence of cylindrical nuclei, in Fig. 9, we show 0 t 2 as a function of L for a neutron star model with 1.4M ⊙ and 12 km, using the EOSs with some parameter sets listed in Table 1 (see Ref. [142] for the concrete parameter sets adopted here), where the left and right panels correspond to the results with N s /N d = 0 and 1, respectively, and the thick-solid line is the fitting with Eq. ( 9).For reference, we also show the L dependence of 0 t 2 excited in the phase of only spherical nuclei discussed in Sec.5.1 with the dashed line.From this figure, one can observe that the modification in 0 t 2 due to the existence of cylindrical nuclei appears only for a neutron star model with a lower value of L. This is because the phase of cylindrical nuclei as well as the other pasta phases becomes narrower with L [120-122], i.e., one can neglect the existence of cylindrical nuclei in the stellar model with larger L. In any case, the fundamental torsional oscillations hardly depend on K 0 and one can discuss the L dependence through the fitting given by Eq. (10).
Using the L dependence of 0 t ℓ with N s /N d = 1 in the phase of cylindrical nuclei, in Fig. 10, we compare the fundamental torsional oscillations with various ℓ for a neutron star model with 1.4M ⊙ and 12 km to the low-lying QPO frequencies observed in SGR 1806-20 (left panel) and SGR 1900+14 (right panel), where the horizontal dashed and dotted lines denote the observed QPO frequencies and the solid lines denote the L dependence of 0 t ℓ .We note that we focus on only the correspondence of the observed low-lying QPOs except for the 26 Hz QPO in SGR 1806-20 here, because the optimal value of L becomes larger than 100 MeV to identify all the observed low-frequency QPOs in terms of the crustal torsional oscillations, as discussed in Sec.5.1, which may be inconsistent with the constraint from existing nuclear experiments.We also note that the 57 Hz QPO additionally discovered in [76] is taken into account this time.From this figure, one sees that the observed QPOs (except for the 26 Hz QPO) in SGR 1806-20 can be identified if L ≃ 73.4 MeV, and those in SGR 1900+14 can be identified if L ≃ 76.1 MeV.
Similarly, one can determine the optimal value of L to identify the observed QPOs with the same set of the fundamental torsional oscillations shown in Fig 10 for various neutron star models with N s /N d = 0 and 1 in the phase of cylindrical nuclei.Assuming that 1.4 ≤ M/M ⊙ ≤ 1.8 and 10 km ≤ R ≤ 14 km as a typical neutron star model, the optimal value of L for identifying the QPO frequencies observed in the SGRs with the same correspondence as shown in Fig. 10, are plotted in Fig. 11 for N s /N d = 0 in the left panel and N s /N d = 1 in the right panel.In this figure, the filled marks with solid lines correspond to the resultant values of L for SGR 1806-20, while the open marks with dashed lines are for SGR 1900+14.Again, since the value of L should be independent of the astronomical events, one has to simultaneously explain both events, SGR 1806-20 and SGR 1900+14, with a specific value of L. Thus, we can get the constraint on L as 53.9 ≲ L ≲ 83.6 MeV for N s /N d = 0 and 58.1 ≲ L ≲ 85.1 MeV for N s /N d = 1, assuming that the central object of SGR 1806-20 and SGR 1900+14 is a neutron star with 1.4 ≤ M/M ⊙ ≤ 1.8 and 10 km ≤ R ≤ 14 km.These constraints on L are shown with the shaded region in Fig. 11.Namely, the uncertainty in N s /N d in the phase of cylindrical nuclei makes only a little difference in the constraint on L, where the allowed L lies in the range of L = 53.9− 85.1 MeV, even if the uncertainty in N s /N d in the phase of cylindrical nuclei is taken into account.This constraint on L gives us the constraint of S 0 as S 0 ≃ 32.0 − 34.4 MeV, using the correlation between L and S 0 given by Eq. (11).
On the other hand, since the torsional oscillations are confined inside the phase of spherical and cylindrical nuclei, we can also discuss the overtone(s) of torsional oscillations.The n-th overtone frequencies of the ℓ-th torsional oscillations, n t ℓ , is theoretically estimated with the crust thickness (or the thickness of elastic region), ∆R, as n t ℓ ∼ v s /∆R [139,150].Meanwhile, ∆R depends on the EOS parameters L and K 0 mainly through the neutron chemical potential at the crust-core boundary [123], when the neutron star mass and radius are fixed.Thus, via the identification of the relatively high-frequency QPO observed in SGR 1806-20, i.e., 626.5 Hz, as the 1st overtone of crustal torsional oscillations, one may obtain information about the EOS parameters [84,142].
Since n t ℓ depends on K 0 and L through ∆R, it is of great use to find a parameter constructed by a combination of K 0 and L, which can characterize n t ℓ .To this end, assuming the combination of (K i 0 L j ) 1/(i+j) with integer numbers i and j, we finally find an appropriate combination of K 0 and L as ς ≡ (K 4 0 L 5 ) 1/9 .
In Fig. 12, we show the 1st overtones of ℓ = 2 torsional oscillations for a neutron star model with 1.4M ⊙ and 12 km constructed with various EOS parameters, where the thick-solid line denotes the fitting formula for the 1st overtones of the torsional oscillations given by [103,142] 12), where the left and right panels respectively denote the results with N s /N d = 0 and 1 in the phase of cylindrical nuclei.The thick-solid line denotes the fitting formula given by Eq. ( 13).Taken from [142].In the left panel, we show the 1st overtones of the ℓ = 2 (solid lines) and 10 (dashed lines) torsional oscillations as a function of ς for a neutron star models with (M, R) = (1.4M⊙ , 10 km), (1.6M ⊙ , 12 km), and (1.8M ⊙ , 14 km) with N s /N d = 1 in the phase of cylindrical nuclei.In the right panel, the ς dependence of the 1st overtone of the ℓ = 2 torsional oscillations for 1.4M ⊙ neutron star models with R = 10 km (solid line), 12 km (dotted line), and 14 km (dashed line) with N s /N d = 1 in the phase of cylindrical nuclei is compared to the 626.5 Hz QPO observed in SGR 1806-20 (dot-dashed line).Taken from [142].
where d (i) ℓn for i = 0, 1, 2 are the adjustable parameters depending on M, R, and N s /N d in the phase of cylindrical nuclei.We will, hereafter, discuss the ς dependence of n t ℓ .
It is well known that the overtone frequencies of torsional oscillations weakly depend on ℓ, unlike the fundamental oscillations [139].As an example, in the left panel, we show the 1st overtones of the ℓ = 2 (solid lines) and 10 (dashed lines) torsional oscillations as a function of ς for a neutron star models with (M, R) = (1.4M⊙ , 10 km), (1.6M ⊙ , 12 km), and (1.8M ⊙ , 14 km) with N s /N d = 1 in the phase of cylindrical nuclei.From this figure, one confirms that 1 t ℓ weakly depends on ℓ, while one also finds that 1 t ℓ strongly depends on the stellar models.So, hereafter, in this review, we focus only on n t 2 to discuss the ς dependence of the overtones.In addition, in the right panel, we show the ς dependence of 1 t 2 for the 1.4M ⊙ neutron star models with different radii with N s /N d = 1 in the phase of cylindrical nuclei together with the 626.5 Hz QPO observed in SGR 1806-20.From this figure, one can observe that the optimal values of ς are 178.5, 149.7, and 107.1 MeV for 1.4M ⊙ neutron stars of R = 10, 12, and 14 km, respectively, if one identifies the 626.5 Hz QPO with the 1st overtone.
In a similar way, one can obtain the optimal values of ς for various neutron star models.In the left panel of Fig. 14, the resultant optimal values of ς are plotted for neutron star models with 1.4 ≤ M/M ⊙ ≤ 1.8 and 10 km ≤ R ≤ 14 km for N s /N d = 1 in the phase of cylindrical nuclei.One can observe that for each R, the optimal ς increases with M. This behavior comes from the fact that ∆R/R (or ∆R with fixed R) decreases with the compactness, M/R [123], which leads to the increases of 1 t 2 (through 1 t 2 ∼ v s /∆R) and that a optimal value of ς becomes larger.Furthermore, we also derive constraint on K 0 through the definision of ς given by Eq. ( 12), i.e., K 0 = (ς 9 /L 5 ) 1/4 , using the optimal values of ς for identifying the 626.5 Hz QPO with the 1st overtone frequency shown in the left panel of Fig. 14 and the optimal values of L for identifying the low-lying QPO frequencies shown in Fig. 11.In the right panel of Fig. 14, such constraints on K 0 are plotted for various neutron star models, where the filled marks with solid lines denote the results for N s /N d = 1, while the open marks with dashed lines denote the results for N s /N d = 0.In the same panel, we also show the experimental constraint on K 0 , i.e., K 0 = 230 ± 40 MeV [4].Adopting this constraint on K 0 as a typical one, although it may still be model dependent, e.g., [5], from the right panel of Fig. 14 one can observe that the neutron star models with M ≃ 1.4 − 1.5M ⊙ for R = 14 km and presumably M ≃ 1.2 − 1.4M ⊙ for R = 13 km are favored by the QPOs observed in SGR 1806-20 up to 626.5 Hz (except for the 26 Hz QPO).That is, a central object in SGR 1806-20 would have a relatively low mass and large radius.We note in passing that a neutron star model with still lower mass and smaller radius than that mentioned above might be acceptable from the right panel of Fig. 14, but one has to assume a larger value of L to construct such a stellar model, which would be presumably inconsistent with the systematic analysis of the mass-radius relation for low-mass neutron stars [44].
Finally, adopting the resultant constraint on the neutron star model of SGR 1806-20, i.e., M ≃ 1.4 − 1.5M ⊙ for R = 14 km and M ≃ 1.2 − 1.4M ⊙ for R = 13 km, the constraints on L shown in Fig. 11 are dramatically improved.Namely, the optimal value of L should be L ≃ 62 − 73 MeV for N s /N d = 1 and L ≃ 58 − 70 MeV for N s /N d = 0 in the phase of cylindrical nuclei.Therefore, we obtain the constraint on L as L ≃ 58 − 73 MeV independently of the uncertainty in N s /N d in the phase of cylindrical nuclei, which is consistent with the existing constraint on L, e.g., [3,7,10].Using the correlation between L and S 0 (Eq.( 11)), one can estimate the corresponding value of S 0 as S 0 ≃ 32.4 − 33.5 MeV.

Torsional oscillations excited in the cylindrical-hole and spherical-hole nuclei
Through Sec.5.1 and 5.2, we discussed that the QPO frequencies observed in SGR 1806-20 and SGR 1900+14 are well identified with the torsional oscillations excited in the phase of spherical and cylindrical nuclei.In particular, by identifying the 626.5 Hz QPO The fundamental frequencies of the ℓ = 2 torsional oscillations excited in the phase of cylindrical-hole and spherical-hole nuclei for a neutron star model with 1.4M ⊙ and 12 km.The thicksolid line corresponds to the fitting formula given by Eq. ( 14).The left and right panels respectively correspond to the results with R = 1 (maximum enthalpy) and R = 0 (minimum enthalpy).Taken from [88].
with the 1st overtone of the torsional oscillations, we showed the possibility that a central object in SGR 1806-20 would be a neutron star with a relatively low mass and large radius, and derived the constraint on L more severe.However, still, we have a missing piece in the identification of the observed QPO frequencies, i.e., the 26 Hz QPO observed in SGR 1806-20.In this subsection, we discuss the possibility of identifying the 26 Hz QPO (and the QPOs additionally discovered in [77]) in SGR 1806-20 with the torsional oscillations excited in the phase of cylindrical-hole and spherical-hole nuclei, keeping the consistency with the identification discussed in Sec.5.2.
To determine the frequencies of torsional oscillations excited in the phase of cylindricalhole and spherical-hole nuclei, one has to know the value of R in Eq. ( 8) in the phases of cylindrical-hole and spherical-hole nuclei, but it is still quite uncertain.So, here we simply consider the extreme cases, i.e., R = 1 for maximum enthalpy and 0 for minimum enthalpy.In Fig. 15, we show the fundamental frequencies of the ℓ = 2 torsional oscillations excited in the phase of cylindrical-hole and spherical-hole nuclei for a neutron star model with 1.4M ⊙ and 12 km, using the EOSs with some parameter sets listed in Table 1 (see Ref. [88] for the concrete parameter sets adopted here), where the left and right panels correspond to the results with R = 1 and 0 in the phase of cylindrical-hole and spherical-hole nuclei, respectively.From this future, we find that the fundamental frequencies excited in the phase of cylindrical-hole and spherical-hole nuclei weakly depend on K 0 and the dependence on L is well fitted with the functional form given by where c(i) ℓ for i = 0, 1, 2 are the adjustable parameters depending on M, R, and R in the phase of cylindrical-hole and spherical-hole nuclei [88].Using this L dependence of 0 t ℓ excited in the phase of cylindrical-hole and spherical-hole nuclei, we will see the correspondence of the observed QPOs.
In the left panel of Fig. 16, we show the identification of the low-lying QPOs observed in SGR 1806-20 with the fundamental frequencies of torsional oscillations excited in the phase of spherical and cylindrical nuclei and the phase of cylindrical-hole and sphericalhole nuclei for a neutron star model with 1.3M ⊙ and 13 km.As shown in the left panel of Fig. 10, the 18, 29, 57, 92.5 Hz QPOs can be identified with the fundamental frequencies of the ℓ = 2, 3, 6, 10 torsional oscillations excited in the phase of spherical and cylindrical nuclei.With this correspondence, we find the optimal values of L as L ≃ 70.8 MeV for N s /N d = 1 and L ≃ 67.5 MeV for N s /N d = 0 in the phase of cylindrical nuclei.We note that in Fig. 16 the solid lines denote the fundamental frequencies of torsional oscillations excited in the phase of spherical and cylindrical nuclei with N s /N d = 1 in the phase of cylindrical nuclei.On the other hand, we also show the fundamental frequencies of ℓ = 2, 3, 4 torsional oscillations excited in the phase of cylindrical-hole and spherical-hole nuclei by the shaded regions, assuming that 0 ≤ R ≤ 1.From this figure, one can observe that the 26 Hz QPO, which cannot be identified with the torsional oscillations in the phase of spherical and cylindrical nuclei, can be identified with the ℓ = 4 fundamental torsional oscillations excited in the phase of cylindrical-hole and spherical-hole nuclei consistently with the optimal value of L given by the identification of the other QPOs with the torsional oscillations in the phase of spherical and cylindrical nuclei.In addition, using this double layer model (lasagna sandwich model), we find that one can identify the QPOs originally discovered in SGR 1806-20 together with the QPOs additionally discovered by a Bayesian procedure, e.g., 51.4,97.3, and 157 Hz QPOs [77], as shown in the right panel of Fig. 16.That is, as shown in the left panel of Fig. 10, the 18, 29, 57, 92.5, and 150 Hz QPOs can be identified with the fundamental frequencies of the ℓ = 2, 3, 6, 10, and 16 torsional oscillations excited in the phase of spherical and cylindrical nuclei.In a similar way, the 157 Hz QPO can be identified with the ℓ = 17 fundamental torsional oscillations in the phase of spherical and cylindrical nuclei, while the 26, 51.4, and 97.3 Hz QPOs are the ℓ = 4, 8, and 15 fundamental torsional oscillations in the phase of cylindrical-hole and spherical-hole nuclei.

Constraint on a neutron star model for GRB 200415A
So far, we have considered the correspondence between the crustal torsional oscillations and the QPO frequencies observed in giant flares.In addition to these observations, another magnetar flare, GRB 200415A, was also detected in the direction of the NGC253 galaxy, where several high-frequency QPOs with varying significance were found at 836, 1444, 2132, and 4250 Hz [78].We note that only high-frequency QPOs have been detected from this event due to the short duration of the observation interval.So, considering the dynamical time of neutron stars, it may be possible to identify these observed QPOs with neutron star oscillations other than the torsional oscillations, such as the fundamental ( f -), gravity (g i -), or shear (s i -) modes (e.g., see Fig. 25).Nevertheless, since these QPOs come from the magnetar flare, here we consider the identification with the crustal torsional oscillations, which is the same framework as discussed in the case of the QPOs observed in giant flares (Sec.5.1, 5.2, and 5.3).
Since the observed QPOs are too high to identify the fundamental torsional oscillations, we consider identifying them with the overtones of torsional oscillations, as in Sec.5.2.Using the fitting of the overtone frequencies given by Eq. ( 13), we plot the ς dependence of n t 2 in the left panel of Fig. 17 for a neutron star model with 1.6M ⊙ and 12 km with N s /N d = 0 in the phase of cylindrical nuclei, where the horizontal shaded regions denote the QPOs observed in GRB 200415A.We note that the overtone frequencies weakly depend on the ratio of N s /N d in the phase of cylindrical nuclei [103], in this subsection, we only consider the case with N s /N d = 0 in the phase of cylindrical nuclei.From this figure, one can observe that the four QPOs can be identified well with the 1st, 2nd, 4th, and 10th overtones of torsional oscillations, which tells us the optimal value of ς is 121.7 MeV.
In a similar way, one can determine the optimal values of ς for various neutron star models.In practice, the optimal values of ς determined with the same combination of the overtones as shown in the left panel of Fig. 17 are plotted in the right panel of Fig. 17.In particular, the neutron star model considered in the left panel is indicated with the arrow.Meanwhile, the value of ς is also estimated with the fiducial value of L and K 0 , i.e., ς = 83.5 − 135.1 MeV with L = 60 ± 20 MeV and K 0 = 240 ± 20 MeV [2], or with the optimal value of L obtained in Sec.5.2 to identify the QPO frequencies observed in the giant flares with the torsional oscillations, such as L = 58 − 73 MeV, with the fiducial value of K 0 , i.e., ς QPO = 104.9− 128.4 MeV.These estimations of ς are also shown in the right panel of Fig. 17 with the shaded region and the enclosed region with the dashed lines.That is, considering these estimations, for example, the stellar model with 1.5M ⊙ and 11 km can be excluded.In consequence, we can get the constraint on the neutron star mass and radius as shown in the left panel of Fig. 18, where the shaded region corresponds to the results with ς = 83.5 − 135.1 MeV, while the bound region by solid lines is those with ς QPO = 104.9− 128.4 MeV.
One may make the constraint on the mass and radius of the neutron star model for GRB 200415A more severe, if one knows the EOSs (or mass-radius relations) corresponding to the estimations of ς.For this purpose, the mass formula for a low-mass neutron star [44] must be crucial.The low-mass neutron star structures may strongly depend on the MeV with the optimal value of L to identify the magnetar QPOs (bound by solid lines).In the right panel, we show the constraint on the neutron star mass and radius obtained from the QPOs in GRB 200415A (the same as the left panel) together with those estimated with the mass formula for a low-mass neutron star as a function of η = (K 0 L 2 ) 1/3 derived in [44] (the right-bottom region), assuming that the central density of neutron stars should be less than twice the nuclear saturation density.Taken from [103].
nuclear saturation parameters, because the central density for such an object is very low.In fact, the mass and gravitational redshift, z, defined by z = 1/ √ 1 − 2M/R − 1 can be well expressed as a function of the stellar central density and a new parameter given by where u c is the central density of the neutron star normalized by the saturation density and η 100 is the value of η normalized by 100 MeV [44].We note that these empirical formulae are valid in the range of 0.9 ≲ u c ≤ 2.0.Since z is a combination of M and R, one can plot the mass and radius relation once one selects the value of η.
The values of ς and ς QPO , which are adopted to make a constraint on the mass and radius of the neutron star shown in the left panel of Fig. 18, are respectively given by L = 60 ± 20 MeV and L = 58 − 73 with K 0 = 240 ± 20 MeV, which correspond to the value of η as η = 70.6 − 118.5 MeV and η QPO = 90.5 − 111.5 MeV, respectively.Using this range of η, one can estimate the stellar mass and radius, whose central density is less than twice the saturation density, as the right-bottom regions in the right panel of Fig. 18.Furthermore, assuming that the radius of the neutron star whose central density is larger than twice the saturation density would be almost constant as shown with the dashed lines corresponding to R = 11.73,12.41, 13.03, and 13.23 km for η = 70.6,90.5, 111.5, and 118.5 MeV, we can get the overlap region with the constraint on the neutron star mass and radius to identify the QPOs with the overtones of torsional oscillations shown in the left panel of Fig. 18.In such a way, we can get the neutron star mass and radius constraint for GRB 200415A.The resultant mass and radius constraint is shown in Fig. 19 with a double-parallelogram, together with the other constraints obtained from the astronomical observations and experimental constraints.
Finally, we comment on the alternative identification of the QPOs observed in GRB 200415A with the overtones of torsional oscillations.To obtain the neutron star mass and radius constraint shown in Fig. 19, we identify the lowest QPO observed in GRB 200415A, QPOs with the overtones of crustal torsional oscillations (double-parallelogram), adopting the fiducial values of L and K 0 (or the optimal value of L to identify the magnetar QPOs).For reference, the other constraints obtained from the astronomical observations and experimental restrictions are also shown.For astronomical observations, as well as those shown in Fig. 1, we plot the constraints of the 1.6M ⊙ neutron star radius, i.e., R 1.6 ≳ 10.7 km [151], the 1.4M ⊙ neutron star radius, i.e., R 1.4 = 11.0 +0.9 −0.6 km [152] or R 1.4 = 11.75 +0.86  −0.81 km [153] obtained from the GW170817; and the massive neutron star mass of PSR J0952-0607, i.e., M = 2.35 ± 0.17M ⊙ [154].Additionally, the mass and radius relations for neutron stars constructed with several EOSs as in Fig. 1 are also plotted.Taken from [103].
i.e., 835.9 Hz QPO, with the 1st overtone of the torsional oscillations, but it may be possible to identify it with the 2nd overtone.In the left panel of Fig. 20, we compare the lowest QPO observed in GRB 200415A with the 1st (dotted line) and 2nd (solid line) overtones for a neutron star model with 1.4M ⊙ and 14 km.In this figure, we also show the range of ς with the fiducial value of L and K 0 (shaded region) and that with the optimal value of L to identify the QPOs observed in the giant flares (the region between the vertical dashed lines).From this figure, one can find that the lowest QPO observed in GRB 200415A cannot be identified with the 2nd overtone at least with this stellar model, consistently with the region of ς estimated with the fiducial value of saturation parameters.In practice, one can identify the four QPOs observed in GRB 200415A with the 2nd, 5th, 8th, and 16th overtones, if ς ≃ 142 MeV, as shown in the right panel of Fig. 20.Still, anyway, it is inconsistent with the range of ς with the fiducial value of saturation parameters.We note that the stellar model considered here is a little larger than the constraint from GW170817, i.e., the 1.4M ⊙ neutron star radius should be less than 13.6 km.But, as shown in the right panel of Fig. 17, the optimal value of ς to identify the observed QPO frequencies becomes larger as the stellar radius decreases with fixed mass.That is, if one considers the 1.4M ⊙ neutron star model with a radius smaller than 13.6 km, the optimal value of ς is more distant from the estimation of ς using the fiducial value of K 0 and L.

Shear and interface oscillations
Owing to the presence of the crust elasticity, the shear (s-) and interface (i-) modes are also excited as well as the torsional oscillations.To see the behavior of the sand i-mode frequencies, again we consider the stellar model with the pasta phases, where the elasticity in the phase composed of slablike nuclei is assumed to be zero.In the same as in the torsional oscillations, one has to consider the effect of the neutron superfluidity, i.e., N s /N d to estimate the value of N s in Eq. (7) and R in Eq. ( 8), but for simplicity we consider only the situation of the maximum enthalpy here, i.e., all the dripped neutron comove with the protons.
In this review, we simply adopt the Cowling approximation to examine the sand i-mode oscillations.Therefore, the perturbation equations for the sand i-modes are MeV.In the right panel, the four observed QPOs are identified with the 2nd, 5th, 8th, and 16th overtones for a neutron star model with 1.4M ⊙ and 14 km.The vertical dashed line denotes the optimal value of ς with this correspondence.Taken from [103].
derived from the linearized energy-momentum conservation laws, which are completely the same as those for the f -, p i -, and g i -modes if one neglects the elasticity.One can find the concrete perturbation equations in [124,155], even though they are not shown explicitly here.We note that the perturbation equations without the Cowling approximation are also found in [61], which can be obtained from the linearized Einstein equations.We note that one can qualitatively discuss the behavior of the frequencies determined with the Cowling approximation, but the quantitative discussion may be problematic at least the f -and p-modes, which deviate ∼ 20% from those determined without the Cowling approximation [156].On the other hand, since the damping rate for the sand i-modes is too small to calculate numerically, the sand i-modes might be determined relatively well even with the Cowling approximation.In fact, it has been shown that the g-modes excited with the density discontinuity, whose damping rate is too small, are accurately determined with the Cowling approximation [56].Anyway, one has to discuss the accuracy of the sand i-mode frequencies with the Cowling approximation somewhere.
As the boundary conditions, one should impose the regularity condition at the center and the condition that the Lagrangian perturbation of pressure should be zero at the stellar surface.In addition, at the interface where the elasticity discontinuously becomes zero, one has to impose the junction conditions, such as the continuity of the radial displacement, the radial and transverse tractions, and Lagrangian perturbation of pressure [155,157].In practice, with our stellar models, one has to impose the junction conditions at the four interfaces, i.e., between the fluid core and the phase of spherical-hole nuclei, between the phase of cylindrical-hole nuclei and the phase of slablike nuclei, between the phase of slablike nuclei and the phase of cylindrical nuclei, and between the surface of the outer crust and envelope.
The transition density at these interfaces becomes crucial for the sand i-mode oscillations.The transition density where the nuclear shape would be changed is determined from the EOS, depending on the saturation parameters, e.g., [120][121][122].On the other hand, the transition density at the outer crust surface is more or less complicated.This is because the density to determine such properties is too low to neglect the thermal effect, even though the thermal effect on the neutron star structure is negligible when the density is very high [158].In this review, we simply consider that the transition density at the surface of the outer crust is 10 10 g/cm 3 .In addition, we set the surface density, ρ s , being 10 6 g/cm 3   The radial profile of the amplitude of the Lagrangian displacement for the i i -modes in the radial (W) and angular direction (V) for a stellar model with 1.44M ⊙ and 10.2 km constructed with the EOS with K 0 = 230 MeV and L = 42.6 MeV.The vertical dotted lines from left to right denote the interface between the fluid core and the phase of spherical-hole nuclei, the interface between the phase of cylindrical-hole and the phase of slablike nuclei, the interface between the phase of slablike nuclei and the phase of cylindrical nuclei, and the interface between the surface of crust and envelope.Taken from [159].
In Fig. 21, as an example, we show the radial profile of the amplitude of the Lagrangian displacement for the i i -modes (i = 1 − 4) in the radial (W) and angular directions (V) for a stellar model with 1.44M ⊙ and 10.2 km constructed with the EOS with K 0 = 230 MeV and L = 42.6 MeV.The classification of i-mode is not uniquely fixed, but here we identify that the i 2 -mode is mainly excited at the interface between the crust surface and envelope; the i 1 -mode is excited at the interface between the crust surface and envelope and at the interface between the phase of slablike nuclei and the phase of cylindrical nuclei; the i 3and i 4 -modes are mainly excited at the interface associated with the cylindrical-hole and spherical-hole nuclei.Since the i-modes are eigenmodes excited due to the existence of the interface where the elasticity discontinuously becomes zero [160], the number of the excited i-modes is generally equivalent to the number of the interface [124,159].With our classification, the order of the i-modes generally becomes i 4 -, i 2 -, i 1 -, and i 3 -mode from the lowest to highest, and typically the frequencies are less than ∼ 100 Hz.Moreover, in Fig. 22, we show how i-mode frequencies depend on the surface density, ρ s , by fixing the other transition density.From this figure, one can observe the i-mode frequencies hardly depend on ρ s .
On the other hand, in Fig. 23, we show the radial profile of the amplitude of the Lagrangian displacement for the s 1 -and s 2 -modes in the radial (W) and angular directions (V) for a neutron star model with 1.4M ⊙ and 12.4 km constructed with the EOS with K 0 = 230 MeV and L = 73.4MeV.Unlike the i-modes, there are an infinite number of the s-modes.As shown in Fig. 23, the s-modes are excited almost only inside the elastic region.Thus, in principle, one may find the s-modes excited inside the phase of cylindrical-hole and spherical-hole nuclei.Even so, we could not find them in the frequency range of less than a few kHz.This situation could be understood as follows.Focusing on the angular  displacement, V, in Fig. 23, one may see that the wavelength of the s i -mode, λ i , is roughly estimated as where ∆R denotes the thickness of the elastic region where the s-mode is confined.Meanwhile, the s i -mode frequency may be estimated as using the shear velocity, v s , [61].With this simple estimation, one can expect that, even if the s-modes are excited inside the phase of cylindrical-hole and spherical-hole nuclei, their frequencies must become much higher than those excited in the phase of spherical and cylindrical nuclei.Namely, ∆R for the phase of cylindrical-hole and spherical-hole nuclei, ∆R CHSH , is much thinner than ∆R for the phase of spherical and cylindrical nuclei, ∆R SpCy , i.e., ∆R CHSH /∆R SpCy ∼ 1/100, which leads to the estimation that the s-mode frequencies excited in the phase of cylindrical-hole and spherical-hole nuclei must be ∼ 100 times higher than those excited in the phase of spherical and cylindrical nuclei.Since the crust thickness, ∆R, also depends on the stellar compactness [119,123], to see the dependence of the sand i-mode frequencies on the neuron star properties, we consider not only the neutron star models constructed with the original OI-EOSs but also those constructed with the EOS, which is constructed in such a way that the OI-EOS for a lower density region (ε ≤ ε t ) is connected to a one-parameter EOS given by p = α(ε For various neutron star models, f i 1 M and f s 1 R are respectively plotted as a function of M/R in the left and right panels, where the circles, diamonds, squares, and triangles denote the results for stellar models constructed with the original OI-EOS, OI-EOS connected to the EOS with α = 1/3, 0.6, and 1, respectively.The solid lines denote the fitting by Eqs.(21) in the left panel and (22) in the right panel, while the values of the EOS parameters (K 0 , L) are also shown on each line.Taken from [159].
for a higher density region (ε ≥ ε t ), where we especially adopt that ε t is equivalent to twice the saturation density [161].In the one parameter EOS, p t is given from the OI-EOS at ε = ε t and α is related to the sound velocity, c s , as c 2 s = α.Thus, α is a kind of indicator for expressing the stiffness of the core (or higher density) region.In this review, we focus on 1/3 ≤ α ≤ 1.
Through a systematical study, we find empirical relations for the i i -mode frequency, f i i , and the s i -mode frequency, f s i .That is, f i i M and f s i R can be expressed as a function of the stellar compactness, x, given by x ≡ M 1.4 /R 12 with M 1.4 = M/1.4M⊙ and R 12 = R/12km, almost independently of the stiffness of the core region, i.e., the value of α, such as where a 0i , a 1i , a 2i , b 0i , and b 1i are the adjusted coefficients depending on the crust properties [124].As an example, in Fig. 24, we show the results for the i 1 -mode (left panel) and the s 1 -mode (right panel), where the solid lines denote the fitting given by Eqs. ( 21) and (22).Thus, once one simultaneously observes the shear and interface mode frequencies from a neutron star, one might extract the neutron star mass and radius via the empirical relations given by Eqs. ( 21) and (22) with the help of the constraint on the crust stiffness obtained from terrestrial experiments.Furthermore, the coefficients in Eqs. ( 21) and ( 22) can be well fitted as a function of a suitable combination of K 0 and L [124].Namely, f i i M and f s i R can be expressed as a function of M/R and a combination of K 0 and L. In practice, the s-mode frequency for a canonical neutron star can be estimated within ∼ 1% accuracy.
In the end, we consider the possibility of applying the s-mode frequency for identifying the QPO observations.As discussed in Sec.5.4, one can identify the high-frequency QPOs observed in GRB 200415A with the overtones of crustal torsional oscillations.Since only high-frequency QPOs have been observed in this event, one may identify them with the other neutron star oscillations.In practice, as shown in Fig. 25 where the f -and s i -mode frequencies for the stellar models with (K 0 , L) = (230, , where the left and right panels correspond to the results for α = 1/3 and 0.6, respectively.The horizontal shaded regions denote the QPO frequencies observed in GRB 200415A [78].The vertical lines denote the optimal value of M/R to identify the four observed QPOs with the s i -modes, i.e., M/R ≃ 0.183 for α = 1/3 and M/R ≃ 0.195 for α = 0.6.Taken from [159]. stellar compactness, together with the QPO frequencies observed in GRB 200415A, one can identify them with the s-modes.Namely, the observed QPOs are identified with the s 1 -, s 2 -, s 4 -, and s 8 -mode frequencies if M/R ≃ 0.183 for α = 1/3 (left panel) or M/R ≃ 0.195 for α = 0.6 (right panel).One could discuss the constraint on the saturation parameters through the identification of the QPOs with the torsional oscillations in Sec.5, because the torsional oscillations are confined inside the elastic region.On the other hand, since the s-mode oscillations depend on not only the crust region (elastic region) but also the (fluid) core region, it may be more difficult to discuss the crust properties via the identification of QPOs with the s-modes.Anyway, to get a severe constraint on the neutron star properties or to validate a theoretical model, we have to wait for the next signals from a neutron star.

Discussion and Conclusions
In this review, we discussed the neutron star oscillations excited by the crust elasticity, such as the torsional, interface, and shear oscillations, adopting the Cowling approximation.Since the torsional oscillations are axial-type oscillations, one can fully discuss them even with the Cowling approximation.On the other hand, since the interface and shear oscillations are polar-type oscillations, one can qualitatively discuss them with the Cowling approximation but one may have to check the accuracy of the interface and shear mode frequencies by comparing them to those taken into account metric perturbations.In any case, through the identification of the QPOs observed in the magnetar flares with the torsional oscillations, we successfully obtain the constraint on the nuclear saturation parameters, especially the density dependence of the nuclear symmetry energy, L, as L ≃ 58 − 73 MeV.On the other hand, via systematical study of the interface and shear oscillations, we find empirical relations almost independent of the stiffness of the core region.
To determine these frequencies, the effect of the superfluidity is important input physics as well as the shear modulus.To deal with the superfluidity, in this review we introduce the parameter, such as N s /N d or R, which is still quite uncertain.In addition, the understanding of the shear modulus inside the pasta phases may still be poor, where one has room for improvement.One may have to improve these properties to consider a more realistic situation.
As shown in this review, the observed QPO frequencies can be identified with the crustal torsional oscillations, but still one has to solve the time-dependence of the QPO excitations.For this purpose, one may have to examine the crustal oscillations beyond linear analysis, e.g., taking into account the nonlinear coupling, and discuss the energy transfer from specific modes to other modes.Moreover, in this review, we consider only the oscillations excited due to the crust elasticity.Meanwhile, the existence of quark matter inside the neutron star is theoretically suggested.In such a situation, the hadronic matter would transition into quark matter, and the mixed phase may appear at the transition region, where the structure becomes similar to the pasta structure in the crust [162,163].So, such a mixed phase behaves as a liquid crystal, i.e., the elasticity becomes non-zero.That is, if such a phase exists, one can expect that the torsional oscillations and/or shear and interface oscillations are excited in the hadron-quark mixed phase [164,165].Through the observations of these oscillation modes, one may probe the properties of quark matter.

Figure 1 .
Figure1.Several constraints on the neutron star mass and radius obtained from the astronomical observations and terrestrial experiments.MSP J0740+6620 is one of the massive neutron stars discovered up to now, whose mass is M/M ⊙ = 2.08 ± 0.07[24], while the radius (and mass) of this object are also constrained via NICER as well as PSR J0030-0451.The tidal deformability restricted from the GW170817 event tells us the upper limit of the radius of 1.4M ⊙ neutron star, i.e., R 1.4 ≤ 13.6 km[26].The observation of the X-ray bursts from neutron stars also gives us the constraint on the neutron star mass and radius, although they depend on the theoretical model[40].In addition, one can theoretically exclude the top-left shaded region due to the causality[41].On the other hand, the terrestrial experiments, e.g., by PREX-II[19], by SπRIT[42], and at RCNP[43], give us the constraint on the stellar models constructed with a lower central density (bottom-right region), adopting the mass formula for low-mass neutron stars[44], where the constraints expected with the fiducial value of the nuclear saturation parameters, i.e., L = 60 ± 20 MeV[3,7,10] and K 0 = 240 ± 20 MeV[2], is also shown (see Eq. (1) for the definition of L and K 0 ).Furthermore, for reference, the mass and radius relations constructed with several EOSs, i.e., SLy4[45,46], SKa[47], SkI3[48], SkMp[49], DD2[50], Shen[51], and Togashi[52], are plotted.Taken from[37].

Figure 2 .
Figure 2. Example of the energy density profile inside the neutron star with 1.4M ⊙ and 12.4 km constructed with a specific EOS with L = 73.4MeV and K 0 = 230 MeV, where L and K 0 are the nuclear saturation parameters given by Eq. (1).The bottom panel is an enlarged view of the above panel.Taken from[124].

Figure 3 .
Figure 3. Shear modulus inside the neutron star model given by Fig.2, where Sp, Cy, CH, and SH respectively denote the phases composed of spherical, cylindrical, cylindrical-hole, and spherical-hole nuclei.The right panel is an enlarged view of the shaded region in the left panel.Taken from[124].

Figure 5 .
Figure 5.For a neutron star model with 1.8M ⊙ and 14 km, the fundamental frequency of ℓ = 2 torsional oscillations is plotted as a function of L, as changing the ratio of N s /N d with each 0.2 from 0 to 1, where N s /N d is assumed to be constant inside the inner crust.Meanwhile, the marks denote the frequencies, using the density-dependent N s /N d derived in[143], and the dashed line is their fitting with Eq. (9).Taken from[85].

Figure 6 .Figure 7 .
Figure 6.The fundamental frequencies of the torsional oscillations with several values of ℓ for a neutron star model with 1.4M ⊙ and 12 km are shown as a function of L. The QPO frequencies observed in SGR 1806-20 and SGR 1900+14 are compared with them in the left and right panels.Taken from[86].

9 Figure 8 .
Figure 8. Alternative possible correspondence of the QPO frequencies observed in SGR 1806-20 (left panel) and SGR 1900+14 (right panel) with the fundamental torsional oscillations for a neutron star model with 1.4M ⊙ and 12 km.Taken from [86].

1 Figure 9 . 17 Figure 10 .
Figure 9.The fundamental frequencies of the ℓ = 2 torsional oscillations excited in the phase of spherical and cylindrical nuclei for a neutron star model with 1.4M ⊙ and 12 km are shown with various marks.The left and right panels are results with N s /N d = 0 and 1, respectively.The thicksolid line denotes the fitting with Eq. (9), while the dashed line denotes the L dependence of 0 t 2 exited in the phase of only spherical nuclei discussed in Sec.5.1.Taken from[142].

Figure 11 .
Figure11.The optimal value of L for identifying the low-lying QPO frequencies (except for the 26 Hz QPO) observed in SGR 1806-20 (filled marks with solid lines) and SGR 1900+14 (open marks with dashed lines) for various neutron star models, where the left and right panels are the results with N s /N d = 0 and 1, respectively.Taken from[142].

1 Figure 12 .
Figure 12.The 1st overtones of ℓ = 2 torsional oscillations for a neutron star model with 1.4M ⊙ and 12 km constructed with various EOS parameters are shown as a function of ς defined by Eq. (12), where the left and right panels respectively denote the results with N s /N d = 0 and 1 in the phase of cylindrical nuclei.The thick-solid line denotes the fitting formula given by Eq. (13).Taken from[142].

1 ( 1 .Figure 13 .
Figure 13.In the left panel, we show the 1st overtones of the ℓ = 2 (solid lines) and 10 (dashed lines) torsional oscillations as a function of ς for a neutron star models with (M, R) = (1.4M⊙ , 10 km), (1.6M ⊙ , 12 km), and (1.8M ⊙ , 14 km) with N s /N d = 1 in the phase of cylindrical nuclei.In the right panel, the ς dependence of the 1st overtone of the ℓ = 2 torsional oscillations for 1.4M ⊙ neutron star models with R = 10 km (solid line), 12 km (dotted line), and 14 km (dashed line) with N s /N d = 1 in the phase of cylindrical nuclei is compared to the 626.5 Hz QPO observed in SGR 1806-20 (dot-dashed line).Taken from[142].

Figure 14 .
Figure 14.In the left panel, the optimal value of ς to identify the 626.5 Hz QPO observed in SGR 1806-20 with the 1st overtone of the ℓ = 2 torsional oscillations are plotted for various neutron star models with 1.4 ≤ M/M ⊙ ≤ 1.8 and 10 km ≤ R ≤ 14 km for N s /N d = 1 in the phase of cylindrical nuclei.In the right panel, the constraint on K 0 obtained by combining the optimal values of ς shown in the left panel (and also those obtained with N s /N d = 0) and the optimal values of L shown in Fig.11is plotted for each stellar model, where filled marks with solid lines (open marks with dashed lines) denote the results for N s /N d = 1 (0) in the phase of cylindrical nuclei.Meanwhile, the shaded region denotes the experimental constraint on K 0 obtained in[4], i.e., K 0 = 230 ± 40 MeV.Taken from[142].
Figure15.The fundamental frequencies of the ℓ = 2 torsional oscillations excited in the phase of cylindrical-hole and spherical-hole nuclei for a neutron star model with 1.4M ⊙ and 12 km.The thicksolid line corresponds to the fitting formula given by Eq. (14).The left and right panels respectively correspond to the results with R = 1 (maximum enthalpy) and R = 0 (minimum enthalpy).Taken from[88].

15 Figure 16 .
Figure16.The correspondence of the QPO frequencies observed in SGR 1806-20 and the fundamental frequencies of torsional oscillations with various values of ℓ excited in the phase of spherical and cylindrical nuclei with N s /N d = 1 in the cylindrical nuclei (solid lines) and the phase of cylindricalhole and spherical-hole nuclei for 0 ≤ R ≤ 1 (shaded regions) for a neutron star model with 1.3M ⊙ and 13 km.We focus on only the low-lying QPOs in the left panel, while in the right panel we also consider the QPOs additionally discovered in[77] together with the QPOs considered in the left panel.The vertical lines denote the optimal values of L to identify the low-lying QPO frequencies observed in SGR 1806-20 except for the 26 Hz QPO with the torsional oscillations excited in the phase of spherical and cylindrical nuclei discussed in Sec.5.2, where the thick line (L ≃ 70.8 MeV) and thin line (L ≃ 67.5 MeV) respectively correspond to the optimal values of L with N s /N d = 1 and 0 in the phase of cylindrical nuclei.Taken from[88].

Figure 17 .
Figure 17.In the left panel, we show the identification of four QPO frequencies observed in GRB 200415A [78] (shaded horizontal regions) with several overtones of torsional oscillations excited in the phase of spherical and cylindrical nuclei for a neutron star model with 1.6M ⊙ and 12 km with N s /N d = 0 in the phase of cylindrical nuclei, where ς = 121.7 is the optimal value with this correspondence.In the right panel, we show the optimal values of ς to identify the four QPOs in GRB 200415A with the same combinations of overtones as in the left panel for various neutron star models.The shaded region is the range of ς = 85.3 − 135.1 MeV estimated with the fiducial value of L = 60 ± 20 MeV and K 0 = 240 ± 20 MeV, while the region in the dashed lines is the range of ς QPO = 104.9− 128.4 MeV with the optimal values of L discussed in Sec.5.2, i.e., L = 58 − 73 MeV, together with K 0 = 240 ± 20 MeV.The arrow in the right panel indicates the neutron star model considered in the left panel.Taken from [103].

5 RFigure 18 .
Figure 18.In the left panel, the constraint on the neutron star mass and radius for GRB 200415A obtained from the optimal values of ς to identify the high-frequency QPOs with 1st, 2nd, 4th, and 10th overtones of torsional oscillations, adopting the range of ς = 85.3 − 135.1 MeV estimated with the fiducial value of L and K 0 (shaded region) or ς QPO = 104.9− 128.4 MeV with the optimal value of L to identify the magnetar QPOs (bound by solid lines).In the right panel, we show the constraint on the neutron star mass and radius obtained from the QPOs in GRB 200415A (the same as the left panel) together with those estimated with the mass formula for a low-mass neutron star as a function of such that M/M ⊙ = 0.371 − 0.820u c + 0.279u 2 c − (0.593 − 1.25u c + 0.235u 2 c )η 100 ,(16)z = 0.00859 − 0.0619u c + 0.0255u 2 c − (0.0429 − 0.108u c + 0.0120u 2 c )η 100 ,

Figure 19 .
Figure19.Mass and radius of GRB 200415A constrained by identifying the observed high-frequency QPOs with the overtones of crustal torsional oscillations (double-parallelogram), adopting the fiducial values of L and K 0 (or the optimal value of L to identify the magnetar QPOs).For reference, the other constraints obtained from the astronomical observations and experimental restrictions are also shown.For astronomical observations, as well as those shown in Fig.1, we plot the constraints of the 1.6M ⊙ neutron star radius, i.e., R 1.6 ≳ 10.7 km[151], the 1.4M ⊙ neutron star radius, i.e., R 1.4 = 11.0 +0.9

Figure 20 .
Figure 20.An alternative possibility to identify the observed QPOs in GRB 200415A with the overtones of torsional oscillations, i.e., the lowest QPO frequency is identified with the 2nd overtone instead of the 1st overtone.In the left panel, the lowest QPO frequency (shaded horizontal region) is compared with the 1st (dotted line) and 2nd (solid line) overtones as a function of ς for a neutron star model with 1.4M ⊙ and 12 km.The vertical shaded region corresponds to ς = 85.3 − 135.1 MeV with fiducial values of L and K 0 , while the region between the vertical dashed lines is ς QPO = 104.9− 128.4 MeV.In the right panel, the four observed QPOs are identified with the 2nd, 5th, 8th, and 16th overtones for a neutron star model with 1.4M ⊙ and 14 km.The vertical dashed line denotes the optimal value of ς with this correspondence.Taken from[103].

Figure 21 .
Figure21.The radial profile of the amplitude of the Lagrangian displacement for the i i -modes in the radial (W) and angular direction (V) for a stellar model with 1.44M ⊙ and 10.2 km constructed with the EOS with K 0 = 230 MeV and L = 42.6 MeV.The vertical dotted lines from left to right denote the interface between the fluid core and the phase of spherical-hole nuclei, the interface between the phase of cylindrical-hole and the phase of slablike nuclei, the interface between the phase of slablike nuclei and the phase of cylindrical nuclei, and the interface between the surface of crust and envelope.Taken from[159].

Figure 22 .Figure 23 .
Figure 22.Dependence of the i-mode frequencies on the surface density, ρ s , for a neutron star model with M = 1.44M ⊙ constructed with K 0 = 230 MeV and L = 42.6 MeV.The other transition density between the phases with and without elasticity is fixed.Taken from [159].

Figure 24 .
Figure 24.For various neutron star models, f i 1 M and f s 1 R are respectively plotted as a function of M/R in the left and right panels, where the circles, diamonds, squares, and triangles denote the results for stellar models constructed with the original OI-EOS, OI-EOS connected to the EOS with

Figure 25 .
Figure25.The f -and s i -mode (i = 1 − 10) frequencies are shown as a function of M/R for neutron star models constructed with the EOS parameter of (K 0 , L) = (230,42.6),where the left and right panels correspond to the results for α = 1/3 and 0.6, respectively.The horizontal shaded regions denote the QPO frequencies observed in GRB 200415A[78].The vertical lines denote the optimal value of M/R to identify the four observed QPOs with the s i -modes, i.e., M/R ≃ 0.183 for α = 1/3 and M/R ≃ 0.195 for α = 0.6.Taken from[159]. .