The static spherically symmetric body in relativistic elasticity
aa r X i v : . [ g r- q c ] J u l The static spherically symmetric body in relativistic elasticity
Jörg Frauendiener ∗ Center of Mathematics for Applications, University of Oslo, P.O. Box 1053 Blindern, NO-0316 Oslo, Norway
Alexander Kabobel † Institut für Astronomie und Astrophysik, Universität Tübingen,Auf der Morgenstelle 10, D-72076 Tübingen, Germany (Dated: November 5, 2018)In this paper is discussed a class of static spherically symmetric solutions of the general relativisticelasticity equations. The main point of discussion is the comparison of two matter models given interms of their stored energy functionals, i.e., the rule which gives the amount of energy stored in thesystem when it is deformed. Both functionals mimic (and for small deformations approximate) theclassical Kirchhoff-St. Venant materials but differ in the strain variable used. We discuss the behaviorof the systems for large deformations.
I. INTRODUCTION
Classical elasticity is a phenomenological theory to describe the properties of solids. It is used heavily in appli-cations such as structural mechanics and other engineering disciplines. This also explains the influence of classicalelasticity in the formulation of the popular finite element method for the numerical solution of partial differentialequations (see e.g. [1]). While it is entirely sufficient to restrict oneself to a classical theory for bulk matter to describeevery-day engineering problems it is nevertheless of conceptual interest to formulate a theory of elasticity which iscompatible with the space-time structure established by Einstein’s theories of relativity.The first attempts to merge elasticity with special relativity go back to the early 20th century and there have beenseveral other formulations including Synge [2] and Rayner [3]. The most influential work, however, has been thepaper by Carter and Quintana [4] who formulated the geometric setting for the theory and derived the basic fieldequations. The theory has also been considered from a field theoretical point of view by Kijowski and Magli [5].Recently, the theory has been analyzed from the point of view of the initial value problem formulation by Beig andSchmidt [6]. They showed that the field equations can be put into a first order symmetric hyperbolic form and theyprove among other things that the Cauchy problem for the system is well-posed under various circumstances. Basedon this formulation it is shown in [7] that there exist solutions of the elasticity equations in Newtonian theory and inspecial relativity describing elastic bodies in rigid rotation. In [8] it is proved that there exist solutions of the staticelastic equations for sufficiently weak gravitational interaction. Losert [9] analyzes the case of a self-gravitatingelastic spherical shell and shows existence of solutions in the Newtonian case.In a series of papers[10, 11, 12] Karlovini, Samuelsson and Zarroug adopt the formulation of Carter and Quintanato discuss spherically symmetric equilibrium configurations and their radial perturbations. They also present anexact static and spherically symmetric solution with constant energy density.Our intention in this paper is to discuss two different equations of state in the static and spherically symmetriccontext with two different materials. In [6] the familiar Kirchhoff-St. Venant stored energy functional for hyper-elastic isotropic materials has been extended to the relativistic case. Recall that this functional is quadratic in thestrain variable and contains the Lamé coefficients as two material constants. Kijowski and Magli in [5] use the samefunctional. However, they adopt a different definition for their strain variable which has a non-linear relationshipto the strain used by Beig-Schmidt. Hence, this results in two different stored energy functionals which have theproperty that, by construction, they agree with the classical Kirchhoff-St.-Venant functional in the non-relativistic,small deformation limit.The plan of the paper is as follows. In sec. II we provide the necessary background on the formulation of thetheory of relativistic elasticity. The exposition follows that of [6]. We present the two formulations by Beig-Schmidtand Kijowski-Magli and point out their differences. It turns out that the only difference is in the definition of thestrain variable which accounts for the above mentioned different energy functionals. ∗ Electronic address: [email protected] † Electronic address: [email protected]
In sec. III we specialize to the static and spherically symmetric case and derive the equations which govern thissituation. We show that this system of equations has a unique smooth solution once the central compression of thebody has been specified.Sec. IV is devoted to a study of various models. In order to compare the two energy functionals we considervarious scenarios. We discuss a solid aluminum sphere and a relativistic highly compact material similar to thenucleonic matter inside a neutron star, described with both theories as well as with the classical theory of elasticity.We also discuss how the choice of different natural states affects the solutions.
II. PRELIMINARIESA. Relativistic Elasticity
The relativistic theory of elasticity in the form that we will use in this work has been described in [6]. The kinematicstructure of the theory can be formulated as follows. As the basic variable one considers a (smooth) map f : M → B (2.1)from space-time ( M , g ) [23] to a 3-dimensional manifold B , the ’material manifold’ or ’body manifold’ or simply ‘thebody’. This is a reference manifold which carries some additional structure which will be described later. The bodymanifold can be interpreted as the collection of all point-like constituents (baryons) of the actual body. Coordinateson B are labels for each individual ‘particle’ of the body. The map f is a map from a 4-dimensional to a 3-dimensionalmanifold so its derivative d f must have a kernel. One requires that F = d f has maximal rank at each point so thatthis kernel has dimension one and is spanned by a unit vector field u a . Using small Latin indices for tensors on M and capital Latin indices for tensors on B , the derivative of f may be written as F Aa = ∂ a f A . In local coordinates ( x a ) on M and ( X A ) on B the map f is given by expressions of the form X A = f A ( x a ) . (2.2)The dynamics of the theory is specified by a Lagrangian density ρ which is regarded as a functional of f and itsfirst derivative d f . In addition, it will depend on the metric g . Thus, the action may be written as A [ f , d f ; g ] = Z ρ [ f , d f ; g ] p − det g d x . (2.3)The Euler-Lagrange equations for this action are G A : = p − det g ∂ a (cid:18)p − det g ∂ρ∂ ( ∂ a f A ) (cid:19) − ∂ρ∂ f A =
0. (2.4)From the properties of f one can already derive several useful consequences. Let Ω ABC be a 3-form on B . This 3-formcan be interpreted as defining a measure on B which gives to each subset of B the number of particles contained init. The pull-back of Ω ABC to M along f is a 3-form ω abc on M which is dual to a vector field ω a . It is clear that thisvector field spans the kernel of f so that it must be proportional to u a . Hence we have the formulas ω abc = Ω ABC F Aa F Bb F Cc = ǫ abcd ω d , ω a = nu a = ǫ abcd ω bcd . (2.5)The proportionality factor n is interpreted as the number density of particles (baryon density) constituting the bodyin the state it acquires when embedded into space-time.The energy-momentum tensor of the theory is defined as usual by the variation of the action with respect to themetric T ab = − ∂ρ∂ g ab + ρ g ab . (2.6)A consequence of the diffeomorphism invariance of the Lagrangian is (see [6]) that ∇ b T ba = G A F Aa (2.7)i.e., that the elastic field equations are satisfied if and only if the energy-momentum tensor is divergence free. This isnot necessarily the case in other field theories, such as e.g., for the Maxwell field.The inverse metric g ab defines a contravariant, symmetric and positive definite 2-tensor H AB on the body by push-forward with the map f (the minus sign is due to our signature) H AB = − g ab F Aa F Bb . (2.8)This characterizes the current state of the body which can vary due to the space-time curvature. In order to describethe variation the conventional way is to compare the actual state with a reference state that is given a priori as afixed structure on the body B . This can be done by postulating the esxistence of a (positive definite) reference metric γ AB on the body manifold which characterizes a ‘natural’ state of the body in which – by definition – there is nostrain[24]. The difference E AB = H AB − γ AB between H AB and the inverse γ AB provides a measure of the ‘size’ ofthe strain on the body. Equivalently, one may use the linear map E = E AB = H AB − δ AB obtained by lowering anindex on E AB with γ AB .Writing ρ = ne where e is the energy per particle then the second Piola-Kirchhoff stress tensor is obtained as thederivative τ AB = ∂ e ∂ H AB . (2.9)Thus, specifying e as a function of the strain provides the stress-strain relation i.e., the equation of state for thematerial under consideration. If the stress tensor τ AB does not vanish in the natural state in which there is no strainthen one talks about a pre-stressed state, otherwise the state is called stress-free or relaxed . We will be concernedonly with a relaxed state. Thus, the energy density has a minimum in the natural state. For most applications it isenough to assume that the energy density is at most quadratic in the strain and we will do so here. Invariance undercoordinate transformations in the body implies that it can depend only on the scalar invariants of E and since we arein three dimensions those invariants which are at most quadratic in E are Tr E and Tr ( E ) . Thus, the energy densitycan be written as ρ = nm (cid:18) + n p Tr ( E ) + q ( Tr E ) o(cid:19) = nm (cid:18) + n p ( H AB H AB − H A A + ) + q ( H A A − ) o(cid:19) where m is the rest mass of a particle and p and q are constants. This is the stored energy functional which isassumed in [6]. It describes the so-called Kirchhoff-St. Venant materials. When we refer below to the Beig-Schmidt(BS) formulation we mean the use of this stored energy functional.The fact that there exists a metric on the body implies that there are now two 3-forms available: the 3-form Ω ABC which gives the number of particles in each sub-domain of the body and the volume form V ABC induced by γ AB which gives the volume of the sub-domain. Since the two forms must be proportional we have Ω ABC = n V ABC , (2.10)thus defining the particle density n in the natural state. This can be used to define the mass density ρ = mn in thenatural state. Using the ‘natural’ particle density n we can obtain the following formula ω abc = Ω ABC F Aa F Bb F Cc = n V ABC F Aa F Bb F Cc = n ǫ abcd u d = n ǫ abc . (2.11)In local coordinates where Ω ABC . = p det γ d X and ǫ abc . = √− det h d x with h ab = g ab − u a u b , we have n p det γ det F = n √− det h . (2.12) B. The Kijowski-Magli strain
The main difference between the Beig-Schmidt [6] and Kijowski-Magli [5] formulations is the choice of the variablewhich measures the deformation. Beig-Schmidt use the difference between the actual and the relaxed metrics on thebody while Kijowski-Magli use a logarithmic variable. They claim that this variable has better behavior when largedeformations are studied.With our choice of conventions and notation this variable is S ab = −
12 log ( u a u b − γ ab ) , (2.13)where γ ab = F Aa F Bb γ AB is the pull-back of the reference metric on B to the space-time. Note, that γ ab is positivedefinite so that the tensor inside the parentheses has only positive eigenvalues and the logarithm is well-defined.Kijowski-Magli write down an action functional in terms of this variable. As before, the scalar character of theaction implies that it can depend only on the scalar invariants of S and Kijowski-Magli assume that it is at mostquadratic in S . They introduce the invariants α = Tr S , β =
12 Tr ˜ S (2.14)where ˜ S is the trace-free part of S . Then, they write the action in the form A = Z n (cid:18) m + A α + B β (cid:19) p − det g d x . (2.15)Here, we have adapted the formula of Kijowski-Magli somewhat because we use the particle density n instead of thematter density and consequently we have to interpret e as the energy per particle. When we refer to the Kijowski-Magli (KM) formulation we mean the use of this stored energy functional.When deriving the equations of motion Kijowski and Magli use familiar techniques from Lagrangian field theory.However, their energy-momentum tensor is the canonical one and not the dynamical (symmetric) one which is ob-tained by varying the action with respect to the metric (see [13] for a thorough discussion of this difference). Sincewe are using the latter tensor we cannot simply take over the expression of Kijowski-Magli. Instead, we need toderive this energy-momentum tensor explicitly as given in appendix A. We obtain T ab = ρ u a u b + n α Ah ab + nB ˜ S ab . (2.16) C. Comparison of the two formulations
In this section we want to compare the two presented formulations of relativistic elasticity. We establish thatthey agree on the linearized level and show how they differ for large deformations. In order to compare these twoformulations we introduce the following variable which measures the deformation from a given state ǫ ab = − ( h ac + γ ac ) h cb = − h ab − γ ab .In terms of ǫ ab we can write the KM deformation tensor in the form S = −
12 log ( + ǫ ) .The BS deformation is E AB = H AB − γ AB . We can relate these two difference deformation variables by the followingcomputation E ab = F Aa F Bb E AB = − h cd F Aa F Bb γ AC γ BD F Dd F Cc − γ AC γ BD γ CD = − γ ac ( γ bd + h bd ) h cd = γ ac ǫ bc . (2.17)It follows that E A A = ( H AB − γ AB ) γ AB = ǫ ab h ab = ǫ aa and also E AB E AB = ǫ ab ǫ ab . Thus, the BS-energy density takesthe form ρ BS = nm (cid:18) + n p ǫ ab ǫ ab + q ( ǫ aa ) o(cid:19) .The KM variables α and β can be expressed in terms of ǫ as well. Thus, e.g., α becomes α = −
12 log det ( + ǫ ) and, similarly, β can be expressed as before in terms of S and hence in terms of ǫ . In order to connect with the BSformulation we expand the energy density up to quadratic terms in ǫ . For the expansion of α and β we find S ab ≈ − ǫ ab , α ≈ − ǫ aa , β ≈ ǫ ab ǫ ab − ( ǫ aa ) ,so that the energy density of Kijowski-Magli up to second order in ǫ is ρ KM = n ( m + A α + B β ) ≈ n (cid:18) m + A ( ǫ aa ) + B ǫ ab ǫ ab − B ( ǫ aa ) (cid:19) .The expressions for the energy density in the two formulations agree in this approximation if we put mp = B , mq = A − B .The coefficients in front of the quadratic terms can be related to the classical elastic constants. Introducing the numberdensity n in the natural state one defines the Lamé coefficients λ = n mq and 2 µ = n mp . Then, n mA becomes thebulk modulus K while µ is the shear modulus of the material.Under these circumstances the energy-momentum tensor is given up to first order terms in ǫ ab by T ab = nn ρ u a u b − nn ( µǫ ab − λǫ cc γ ab ) .We contrast this with the exact energy-momentum tensorsBS: T ab = nn ρ u a u b + nn ( µ γ ac ǫ cb + λ ǫ cc γ ab ) (2.18)KM: T ab = nn ρ u a u b + nn ( µ h ab α + λ S ab ) . (2.19)In the case of no deformation i.e., at a point where the body is in the natural state one has γ ab = − h ab .The BS-energy-momentum tensor reflects this relationship. It almost agrees with the linearized energy-momentumtensor except that γ ab appears instead of − h ab . Both theories are quadratic in their respective deformation variablesand therefore describe in some sense a Hookean theory in which stress and strain are proportional. However, therelationship between the two different strain variables is highly non-linear. While the two energy-momentum tensorsagree for small strain they disagree heavily for large deformations. Similarly, the stored energy functionals whichgive the energy per particle as a function of strain are completely different in the two cases when viewed in terms ofthe strain variable ǫ ab . Thus, the two formulations describe materials with different equations of state. Both materialsbehave like the usual Kirchhoff-St. Venant materials for small strain, but have a completely different behavior forlarge deformations. We want to explore some of the consequences of these differences in the remainder of this article.In the KM formulation the strain variable is defined in terms of the difference tensors in space . This results in aterm proportional to h ab in the energy-momentum tensor, i.e. a term which is isotropic in space. In contrast, in [6] isused the difference between the actual and the reference state on the body as the basic variable. This results in a termproportional to γ ab , i.e., isotropic on the body. This has the consequence that it is much easier to describe a fluid as aspecial case of elastic material within the KM framework than in the BS case. III. SPHERICAL SYMMETRY
Now we specialize to spherical symmetry. We take the space-time metric to be the general spherically symmetricand static metric g = e η dt − e ξ dr − r ( d θ + sin θ d φ ) (3.1)and we assume that the body metric γ is spherically symmetric as well[25], i.e., when expressed in polar coordinates γ AB . = e ξ dR + R ( d Θ + sin Θ d Φ ) . (3.2)Since the geometry at the origin should be regular we need e ξ = ξ ( ) =
0. The function ξ and η depend on r while ξ depends only on R . The map f : M → B is assumed to be equivariant and thus without lossof generality it can be expressed as f ( r , θ , φ ) = ( F ( r ) , θ , φ ) (3.3)for some function F ( r ) with F ( ) =
0. Then, the deformation gradient is given as F Aa . = F ′ d r ⊗ ∂ R + d θ ⊗ ∂ Θ + d φ ⊗ ∂ Φ . (3.4)Clearly, because of staticity we must have u a . = e − η ∂ t and from (2.12) nn = F ′ e ξ F sin θ e ξ r sin θ = F ′ e ξ − ξ F r . (3.5)The pull-back γ ab of the reference metric on B is γ ab . = ( e ξ F ′ ) dr + F ( d θ + sin θ d φ ) (3.6)With these formulas and the abbreviations x = F ′ e − ( ξ − ξ ) and y = F / r we can compute the deformation tensor ǫ ab . = (cid:16) x − (cid:17) dr ⊗ ∂ r + (cid:16) y − (cid:17) ( d θ ⊗ ∂ θ + d φ ⊗ ∂ φ ) (3.7)Using this variable and the formulas (2.18) and (2.19) we can find the energy-momentum tensors in both cases. Theyare given explicitly in the appendix B. A. The equations
The Einstein equations in the spherically symmetric and static case are well known, see e.g. [14]. They are G = e − ξ (cid:26) r (cid:16) − e ξ (cid:17) − r ξ ′ (cid:27) = − π ρ , (3.8) G = e − ξ (cid:26) r (cid:16) − e ξ (cid:17) + r η ′ (cid:27) = π P , (3.9) G = e − ξ (cid:26) η ′′ + ( η ′ ) − η ′ ξ ′ − r (cid:0) η ′ − ξ ′ (cid:1)(cid:27) = π Q , (3.10)where we have put T = ρ , P = − T and Q = − T . These are three equations for the unknown functions ξ , η and f (the function ξ which specifies the reference metric is considered as given). A consequence of the Einsteinequations is that the divergence of the energy-momentum tensor vanishes identically ∇ a T ab =
0. (3.11)Under the current conditions this equation has only one non-trivial component P ′ + η ′ ( ρ + P ) + r ( P − Q ) =
0. (3.12)In order to obtain a useful system one replaces the equation (3.10) by (3.12). Furthermore, one integrates (3.8) byintroducing the mass function M ( r ) = π Z r ρ ¯ r d ¯ r (3.13)or, equivalently, the mean density w = M ( r ) / r to obtain e − ξ = − r w ( r ) . (3.14)Inserting this into (3.9) one can solve for η ′ and insert this into (3.12). Then the following system of equations isobtained rw ′ = − w + πρ , (3.15) r η ′ = r π P + w − r w , (3.16) rP ′ = − r π P + w − r w ( ρ + P ) − ( P − Q ) . (3.17)This system is somewhat deceptive, because ρ , P and Q are functions of F and its derivatives. Since they contain F and F ′ in a non-linear way the third equation gives a complicated non-linear equation for F ′′ . Equivalently, we willregard these functions as depending on x and y defined above. Then P ′ = P x x ′ + P y y ′ . From their definition we geta relationship between x and y ry ′ = xe ξ − ξ − y (3.18)which can be used to substitute for y ′ . With this preparation we now have the following final system of equations rw ′ = − w + πρ , ry ′ = x √ − r w e − ξ − y , rx ′ = P y P x (cid:18) y − x √ − r w e − ξ (cid:19) − r π P + w − r w ρ + PP x − P x ( P − Q ) . (3.19)Once a solution of this system is found we can obtain η by integrating (3.16), e − ξ is given by (3.14) and F is foundfrom the definition of x . The functions ρ , P and Q are specified by the choice of the elastic model as functions of x and y , while e − ξ is any given function of r characterizing the natural state of the body. It is only restricted by havingthe value of unity at the origin. B. Behavior at the origin
The equations are singular at the origin r = r =
0. Then the left hand sides vanish and from the right hand sides we get w ( ) = π ρ ( ) , y ( ) = x ( ) = : a , P ( ) = Q ( ) . (3.20)This shows us that the only free datum is the value a . It characterizes the volume change of the body at its cen-ter. Since the body should be compressed we assume that a = lim r → F ( r ) / r >
1. The initial value for w can becomputed from the expression of ρ in terms of x and y . The third condition states that in the center the radial andthe tangential stresses should be equal and this is a condition on the matter model which cannot be influenced byspecifying initial conditions. The fact that the central compression is enough to characterize a solution uniquely isphysically reasonable and corresponds to the fact that a static fluid configuration is uniquely characterized by thecentral pressure.In order to show that with these initial conditions there exist regular unique solutions we apply the theorem byRendall and Schmidt [15]. The verification of the conditions necessary for that theorem are somewhat lengthy andwe refer the reader to appendix C. It follows from this analysis that for a given value a there exists a unique andsmooth solution of the system of equations (3.19) in a neighborhood of the origin. IV. NUMERICAL MODELING OF SPHERICAL ELASTIC BODIESA. The models
In the rest of this paper we solve the system (3.19) for several specific matter models. We consider two situations,a sphere consisting entirely of an ordinary material such as aluminum and a sphere which consists of material whichresembles the neutron star crust. In both cases we choose the two different energy functionals corresponding tothe BS and KM formulation, respectively. For aluminum we use the values ρ = , λ = µ =
25 GPa.For the neutron star matter we follow the presentation in [16] where the structure of the neutron star crust isdescribed in detail. In the crust of a neutron star the density increases from the outer layer with 10 g/cm to theinner edge where the density is approximately 10 g/cm . While the ground state of the matter is a lattice which hasanisotropic elastic properties it is customary to approximate it by a homogeneous and isotropic elastic material. Thismaterial is under high pressure and hence it is much easier to shear it than to compress it[26]. In fact, the material isoften assumed to be incompressible. The shear modulus of the matter in the neutron star crust has been calculatedin e.g. [17] and we use a value of µ = Pa. The fact that the crust material is almost incompressible means thatthe bulk modulus K is very large compared to the shear modulus and hence that also the Lamé coefficient λ is verylarge. We take it here three orders of magnitude larger than the shear modulus, i.e., λ = Pa.In the relativistic theories there is no canonical choice for the relaxed state of the body. While it seems natural tospecify a flat metric on B this is not necessary. The choice of the relaxed metric has been discussed in the litera-ture [18, 19]. We follow here a suggestion by Carter [20] according to which one can obtain the relaxed state by thefollowing procedure. One assumes the body is heated up until it melts and then one lets it cool down until it solidi-fies again. Assuming that the fluid phase is an ideal fluid then the body settles in a state which can be described by asolution of the perfect fluid equations. Hence, besides a flat metric we also consider the spatial metric correspondingto an incompressible fluid with a constant density ρ , i.e, we put e − ξ = − π ρ r . (4.1)As a third formulation we consider the classical non-relativistic theory of elasticity. The equations for the classicaltheory can be obtained from the relativistic equations as the Newtonian limit, see [6]. The difference to the relativisticequations is that one puts ξ = ξ =
0, so that x = F ′ . Furthermore, eq. (3.12) is replaced by P ′ + η ′ ρ + r ( P − Q ) =
0, (4.2)where η ′ is the gravitational force, determined from the equation η ′ = M ( r ) r = rw . (4.3)The stress components P and Q have the same functional form in terms of x and y as those for the BS-energy-momentum tensor, while ρ is the mass density in the actual state, given by ρ = nn ρ = xy ρ .Hence, the non-relativistic system is rw ′ = − w + πρ , ry ′ = x − y , rx ′ = P y P x ( y − x ) − r P x w ρ − P x ( P − Q ) . (4.4)All the numerical solutions have been obtained using the Runge-Kutta ODE solver suite provided in MATLAB.The calculation is started with an initial value a for x ( ) = y ( ) which is used to calculate the initial value for w ( ) = π /3 ρ ( ) from the energy-momentum tensor. The calculation stops when P vanishes, indicating that theboundary of the body has been found. B. Numerical examples
We first study the aluminum sphere for the three formulations of elasticity. Clearly, for small values of the relativecentral compression δ = a − δ = × kg. The figure shows the relative difference ∆ w = ( w − w N ) / w N betweenthe classical solution w N and the BS and KM solutions w , respectively. The BS solution is indistinguishable fromthe classical solution, the maximum value of the relative difference being 6 × − , while the KM solution alreadyindicates its general property: the system is more tightly bound than in the classical or BS case. Still, in this situationof small relative central compression the maximal difference is only 7.5 × − . −6 r / [km] ∆ w Kijowski−MagliBeig−Schmidt
FIG. 1: Relative difference between BS resp. KM and the classical solution for an aluminum sphere with relative central compres-sion of δ =
1. Aluminum with BS-formulation
Let us now look at the BS model in more detail. The radial pressure is given by (B2) P = ( xy ) x (cid:16) ( λ + µ )( x − ) + λ ( y − ) (cid:17) . (4.5)On the boundary of the body, this expression vanishes. This can happen either when x = y = y ′ that as long as x remains positive we have y ′ > y = y cannot vanish before x vanishes. Thus, on the boundary we have either x = ( x , y ) lies on the ellipse defined by ( λ + µ ) x + λ y = λ + µ . (4.6)In Fig. 2 we show a sequence of such final pairs ( x , y ) obtained from initial values δ in the interval [ ] . Obviously,both cases discussed above can occur. For small relative central compressions the final pair ( x , y ) lies on the ellipseand for increasing compression it moves towards the y-axis until it hits it for an initial value of δ ≈ y . The vanishing of x at the boundary means that the radial distancebetween two adjacent particles there becomes infinite, i.e. the body ruptures. Imagine a large elastic sphere withoutgravitational self-interaction being compressed so that the central compression is above the critical value. Whengravity is switched on, the sphere will be divided into a central piece and a shell at the radius where x vanishes.The equation for x ′ in (3.19) shows that x vanishes with an infinite negative slope because the leading term on theright hand side goes like 1/ x near x =
0. Thus, the solution becomes singular just at the boundary.The two different cases just discussed can also be seen in the behavior of the mass-radius diagram in Fig. 3, wherewe display radius and mass of the aluminum spheres corresponding to relative central compressions δ ∈ [ − , 10 ] .We plot it in double logarithmic and linear axes. The curve shows three different regimes, the classical one where M ∝ R (indicated by the solid line) and an ‘extreme’ regime where M ∝ R , indicated by the dashed line andfinally a ‘linear’ regime with M ∝ R where mass and radius decrease with increasing central compression. Thecross indicates the configuration which is closest to the critical configuration where the radial strain x vanishes.This mass-radius diagram should be compared with Figure 1 from [12]. The similarity of the qualitative behavior isobvious. Karlovini and Samuelsson argue that the branch from the maximal mass towards zero is unstable and wedo find numerical indications of this here as well. Increasing the central compression beyond the value needed for0 −0.2 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.600.511.522.533.544.5 x y FIG. 2: A sequence of pairs ( x , y ) for relative central compressions 0.01 ≤ δ ≤ −2 −1 −4 −2 Radius [R E ] M a ss [ M E ] Radius [R E ] M a ss [ M E ] FIG. 3: Double logarithmic (left) and linear (right) plot of the mass-radius diagram for BS aluminum spheres with relative centralcompression δ , log δ ∈ [ −
3, 3 ] , in units of the earth mass and radius, respectively. The solid line is M ∝ R the maximal mass configuration we observe that we can generate the smaller configurations up to a certain valueof δ depending on the required precision. Beyond this value the solver suddenly settles to a solution which yieldsa configuration in the ‘eye’ inside the mass-radius diagram. This dot in fact contains nine different configurations.The location of the ‘eye’ is roughly at the mass resp. radius for which the radius resp. the mass are maximal on thecurve. The behavior of this system close to the eye should be analyzed in much more detail using more accuratesolution methods.
2. Aluminum with KM formulation
In the formulation of Kijowski-Magli the radial pressure P is given by (B8) P = ( xy ) (( λ + µ ) log x + λ log y ) . (4.7)1As before, at the boundary we have either x = ( λ + µ ) log x + λ log y =
0, (4.8)the case y = y = x + α , α = µλ >
0. (4.9)This curve approaches the y-axis but never intersects it. This indicates that only the case when the final pair lies onthe curve does occur. This is in fact confirmed in Fig. 4 where we show the final pairs ( x , y ) for aluminum spheres −0.2 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6012345678 x y FIG. 4: A sequence of pairs ( x , y ) for relative central compressions 0.001 ≤ δ ≤ with relative central compression δ ∈ [ ] . While in the BS case the value of y can grow arbitrarily, this is notthe case here. In fact, the numerical investigations show that the exhibited value of y ≈ y can achieve. This behavior can be understood when we show the mass-radius diagram for KM aluminum spheresin Fig. 5 which shows a peculiar spiral. The maximal value of y is reached at the same point as the maximal radius.Thus, it is not possible with the KM formulation to create arbitrarily large objects. There exists a maximal mass anda maximal radius for KM aluminum spheres achieved for different objects and there exists a region where a KMaluminum sphere of a given radius can have at least four different masses. It looks like the sequence converges to alimit point. We have not been able to prove this rigorously.
3. The neutron star matter
We have also looked at an exotic material which is somewhat similar to the nucleonic matter that is assumed to bepresent in neutron stars. We show in Fig. 6 and in Fig. 7 the mass-radius diagrams for the neutron star like matterdistributions with the BS and KM stored energy functionals. In both cases the diagrams look qualitatively the sameas those for aluminum except that the size of the configurations are orders of magnitudes different. In the KM casewe find a spiral as before while in the BS case we have the ‘loop’ with a linearly decreasing branch. Again, thisbranch seems to be unstable and the final dot in the diagram is the last for which we could generate a configuration.This shows that there is no qualitative difference in the behavior of aluminum and the exotic matter. This mightchange if one would use the high-pressure formulation developed by Carter and Quintana [4].2 e ] M a ss [ M e ] FIG. 5: Mass-radius diagram for KM aluminum spheres with relative central compression δ , log δ ∈ [ −
3, 3 ] , in units of the earthmass and radius, respectively. The solid line indicates the curve M ∝ R . Radius [km] M a ss [ M S ] FIG. 6: Mass-Radius diagram for BS neutron star like matter with relative central compression δ , log δ ∈ [ −
3, 0.114 ] , in units ofsolar mass and kilometers.
4. The role of the relaxed metric
As discussed above we employ two possible choices for the metric of the relaxed state of an elastic configuration.To compare the two different scenarios we compute configurations with the same relative central compression δ forvalues of δ between 10 − and 1 for the two energy functionals and the two possible materials. In Fig. 8 we show the3 M a ss [ M S ] FIG. 7: Mass-Radius diagram for KM neutron star like matter with relative central compression δ , log δ ∈ [ −
3, 3 ] , in units ofsolar mass and kilometers. −3 −2 −1 −14 −12 −10 −8 −6 −4 −2 δ ∆ w / w BS−ALBS−NSKM−ALKM−NS
FIG. 8: Maximal relative difference in the mean density w between calculation with flat and curved metric for the relaxed statefor given δ with log δ ∈ [ −
3, 0 ] . behavior of the maximal absolute value ∆ ww = max r (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) w c − w f w f (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) of the relative difference of the mean densities w f and w c for the flat and curved cases resp. as a function of δ .Obviously, in the given range of δ the difference between the two configurations is almost negligible. The differenceis larger for the exotic neutron star like material than for aluminum. The maximal difference is reached for the BS-4energy functional with roughly 3%. For increasing δ the differences in three cases reach a maximum and afterwardsdecrease again. With increasing δ the elastic energy in the configuration increases with respect to the gravitationalrest mass energy. Thus, the more the energy of the configuration is dominated by the elastic energy the smaller is theinfluence of the choice of a relaxed state. In any case, what can be learned from Fig. 8 is, that for practical purposesone can safely assume that the metric of the relaxed state is flat. V. CONCLUSION
We have discussed in this work the spherically symmetric body in relativistic elasticity for two different storedenergy functionals. We find that the BS-functional corresponding to the classical Kirchhoff-St. Venant materialsand the KM-functional have entirely different behavior for large deformations even though they agree for smalldeformations. The BS-functional gives rise to a mass-radius diagram which qualitatively is very similar to the onefound by Karlovini and Samuelsson in [12]. They obtain this diagram for a stiff ultra-rigid equation of state inthe Carter-Quintana high-pressure formulation. They find that the decreasing branch is unstable. We can confirmthis numerically and we even see indications of another region of configurations. This is an indication that the BS-functional gives rise to an increasingly stiff equation of state quite in contrast to the KM-functional for which theequation of state becomes increasingly soft. The result of this softness can be seen in qualitatively very differentbehavior of the mass-radius diagram which shows a spiral which approaches a limit point for large deformations.We looked at these functionals for two different materials, the ‘every day’ material aluminum and an artificialexotic material. While the sizes of the individual configurations are very different the qualitative behavior is verysimilar in both cases.In order to analyze in more detail the features in these configurations and in particular the stability properties ofthe different branches of the mass-radius diagrams it might be advantageous to formulate the problem not as aninitial value problem as we have done here. Instead of specifying the central compression and integrating outwardsto (possibly) find the boundary of a configuration one would instead set up a boundary value problem on the bodysubject to the boundary conditions imposed by the symmetry requirements in the center and the vanishing of theradial pressure on the boundary. Steps in this direction have already been made by Losert [9].Of course, our considerations are to a certain point academic because any real material will break at already quitemoderate deformations compared to the ones we have used. But we feel that such questions of principle may shedsome light on the differences between the various possible choices and therefore on the justification of assumptionsmade when relativistic elasticity is used for real problems.
Acknowledgments
The authors are very grateful to Robert Beig and Bernd Schmidt for several very valuable discussions. This workwas supported by a grant from the Deutsche Forschungsgemeinschaft.
APPENDIX A: THE SYMMETRIC ENERGY-MOMENTUM TENSOR FOR THE KIJOWSKI-MAGLI ACTION
The action for the Kijowski-Magli formulation of elasticity is (2.15) A = Z n (cid:18) m + A α + B β (cid:19) p − det g d x .The energy-momentum tensor for this action is obtained (with our conventions) by variation of A with respect tothe inverse metric T ab = − ∂ρ∂ g ab + ρ g ab . (A1)The energy density ρ is specified in terms of the variables α = S aa and β = S ab S ba − ( S aa ) . The tensor S ab wasdefined in (2.13) S ab = −
12 log (cid:16) u a u b − γ ab (cid:17) . (A2)5To find the variation in these variables we first need to compute the variation with g ab in the tensor K ab = u a u b − γ ac g cb .We sketch the calculations here without going too much into the details. We start with the variation of (2.5) to obtain δω d = (cid:16) g ab δ g ab (cid:17) ω d .Since, n = g ab ω a ω b we can now obtain the variation in n δ nn = h ab δ g ab and with ω a = nu a this yields the variation of u b and u a δ u b = (cid:16) u c u d δ g cd (cid:17) u b , δ u a = (cid:16) u c u d δ g cd (cid:17) u a − g ac u d δ g cd . (A3)From this we find δ ( u a u b ) = (cid:16) u c u d δ g cd (cid:17) u a u b − g ac u d u b δ g cd = − h ac u d u b δ g cd (A4)and hence δ K ab = − h ac u d u b δ g cd − γ ac δ g cb . (A5)To compute δα we use the formula Tr log A = log det A (A6)which follows from the corresponding well-known equationdet exp B = exp Tr B ,valid for any matrix B . Thus, we have δα = − δ Tr log K = − δ log det K = −
12 1det K δ det K .Using the formula δ det K = det K Tr ( K − δ K ) we get the final result δα = −
12 Tr ( K − δ K ) = − ( K − ) ba δ K ab . (A7)The variation of β follows from the formula which can easily be derived δ Tr f ( A ) = Tr ( f ′ ( A ) δ A ) .We have β = Tr ( S ) − ( Tr S ) so we first compute δ Tr ( S ) = δ Tr ( log ( K ) ) =
12 Tr ( log ( K ) K − δ K ) = − Tr ( SK − δ K ) = − S ab ( K − ) bc δ K ca .With this result we obtain 2 δβ = − Tr ( SK − δ K ) +
13 Tr S Tr ( K − δ K ) = − Tr ( ˜ SK − δ K ) . (A8)6These formulas can be simplified as follows. Using the fact that S ab u a = K ab u a = u b we get ( K − ) ab δ K ba = − ( K − ) ab ( h bc u d u a δ g cd + γ bc δ g ca ) = − ( K − ) ab γ bc δ g ca = − ( K − ) ab ( u b u c − K bc ) δ g ca = ( − u a u c + g ac ) δ g ca = h ab δ g ab .Therefore, we have δα = − h ab δ g ab and, similarly, δβ = −
12 ˜ S ab δ g ab .Now we can write down the variation of the energy δρ − δ nn ρ = n ( A αδα + B δβ ) = − n (cid:2) A α h ab + B ˜ S ab (cid:3) δ g ab .With this in hand we can now finally write down the energy-momentum T ab = ρ u a u b + n α Ah ab + nB ˜ S ab . APPENDIX B: ENERGY-MOMENTUM TENSORS
Using the abbreviations x = F ′ e − ( ξ − ξ ) and y = F / r the energy-momentum tensor for the Beig-Schmidt theoryis obtained from (2.18) using the expression (3.7) for ǫ ab in spherical symmetry. It has the following non-trivialcomponents T = ρ = xy (cid:18) ρ + µ n ( x − ) + ( y − ) o + λ n ( x − ) + ( y − ) o (cid:19) , (B1) T = − P = − ( xy ) x (cid:16) ( λ + µ )( x − ) + λ ( y − ) (cid:17) , (B2) T = T = − Q = − ( xy ) y (cid:16) ( λ + µ )( y − ) + λ ( x − ) (cid:17) . (B3)The energy-momentum for the Kijowski-Magli formulation in spherical symmetry is obtained in a similar wayfrom (2.19) once the tensor S ab and its invariants α and β are determined. Since S = − log ( + ǫ ) and ǫ is di-agonal in spherical symmetry this is straightforward: S ab . = − ( log x ) dr ⊗ ∂ r − ( log y ) ( d θ ⊗ ∂ θ + d φ ⊗ ∂ φ ) . (B4)Hence, we get α = − ( log x + y ) , (B5)and β = h ( log x ) + ( log y ) i − ( log x + y ) . (B6)With these expressions we obtain the non-trivial components of the energy-momentum tensor T = ρ = xy (cid:18) ρ + µ n ( log x ) + ( log y ) o + λ { log x + y } (cid:19) , (B7) T = − P = − ( xy ) (( λ + µ ) log x + λ log y ) , (B8) T = T = − Q = − ( xy ) (( λ + µ ) log y + λ log x ) . (B9)7 APPENDIX C: EXISTENCE AND REGULARITY OF THE SOLUTIONS
Our goal here is to see whether the theorem by Rendall and Schmidt can be applied to our situation. For easyreference, we cite the theorem here
Theorem C.1 (Rendall and Schmidt)
Let V be a finite dimensional vector space, N : V → V a linear map all of whoseeigenvalues have positive real parts, and G : V × ( − ǫ , ǫ ) → V and g : ( − ǫ , ǫ ) → V smooth maps, where ǫ > . Then, thereexists δ < ǫ and unique bounded C function u : ( − δ , 0 ) ∪ ( δ ) → V which satisfies the equations duds ( s ) + Nu ( s ) = sG ( s , u ( s )) + g ( s ) (C1) Moreover, u extends to a smooth solution of (C1) on ( − δ , δ ) . If N, G and g depend smoothly on a parameter t and the eigenvaluesof N are distinct, then the solution depends smoothly on t. To analyze the behavior of the equations at the center we follow the paper by Park [21] who has found that theequations for the Kijowski-Magli allow smooth regular solutions. We need the same result for the equations comingfrom the formulation of Beig and Schmidt.It is easy to see that a regular spherically symmetric scalar function u ( r ) on a spherically symmetric space-timecan smoothly be extended as an even function of r to negative values of the radius. Thus, they depend smoothly on s = r and their derivative at the center vanishes. Using s as the independent parameter instead of r our system (3.19)is 2 sw ′ + w = πρ , (C2)2 sy ′ + y = x f (C3)2 sP x x ′ = P y ( y − f x ) − s f ( w + π P )( P + ρ ) − ( P − Q ) , (C4)where we have introduced f = ( − sw ) − . The solutions of this system will be regular only if x ( ) = y ( ) = a , P ( ) = Q ( ) and w ( ) = π /3 ρ ( ) . Thus, we can write x = a + sx ( s ) , y = a + sy ( s ) and derive equations for the functions x ( s ) and y ( s ) . This yields the equations2 sw ′ + w = πρ ,2 sy ′ + y = a f − s + x f ,2 sx ′ + x = P y P x (cid:18) − fs a + ( y − f x ) (cid:19) − f P x ( w + π P )( P + ρ ) − s P − QP x .The desired form is 2 su ′ + Nu = sG ( u , s ) + g ( s ) .The first equation is already in this form. In the second equation, the right hand side contains non-linear terms in w (via f ) but these terms have no factor s in front. So we expand f = ( − sw ) − = + sw + O ( s ) and get 2 sy ′ + y = aw + x + O ( s ) .Thus, this equation is in the desired form. Note, that it is not important for this analysis to know the O ( s ) exactly. Itis enough to know that they are well behaved at the origin. The O ( s ) implies that they have the desired factor of s infront.The third equation is more difficult.Here, the right hand side of the equations contains non-linear terms in w , x and y through the dependence of f and the functions P , ρ and Q . Therefore, we have to expand the right hand side8around s = s . So wewrite P = P ( x , y ) = P ( a + sx , a + sy ) = P ( a , a ) + s ( x P x ( a , a ) + y P y ( a , a )) + O ( s ) and similarly for all the other functions on the right hand side. We write P ( a , a ) = P and similarly for the otherfunctions. First, we observe that we need to have P − Q = P x = sx ′ + x + w P + ρ − aP y P x + y P y − Q y P x + x P y − Q x P x = − π P ρ + P P x + O ( s ) .This is the desired form and we can collect the coefficients of the linear terms in the matrix NN = − a − P + ρ − aP y P x P y − Q y P x + P y − Q x P x .For the theorem to apply the eigenvalues of this matrix must be positive. This leads to the following conditions: if P x > P x + P y − Q x − Q y >
0, (C5)7 P x + P y − Q x >
0, (C6)while for P x < [1] T. J. R. Hughes, The finite element method (Dover, 2000).[2] J. L. Synge, Math. Z. , 82 (1959).[3] C. B. Rayner, Proc. Roy. Soc. London A , 44 (1963).[4] B. Carter and H. Quintana, Proc. Roy. Soc. London A , 57 (1972).[5] J. Kijowski and G. Magli, J. Geom. Phys. , 207 (1992).[6] R. Beig and B. G. Schmidt, Class. Quant. Grav. , 889 (2003).[7] R. Beig and B. G. Schmidt, Classical and Quantum Gravity (2006).[8] L. Andersson, R. Beig, and B. G. Schmidt, Static self-gravitating elastic bodies in einstein gravity (2006), URL http://arxiv.org/abs/gr-q /0611108 .[9] C. M. Losert,
Static elastic shells in einsteinian and newtonian gravity (2006), URL http://arxiv.org/abs/gr-q /0603103 .[10] M. Karlovini and L. Samuelsson, Classical and Quantum Gravity , 3613 (2003).[11] M. Karlovini, L. Samuelsson, and M. Zarroug, Classical and Quantum Gravity , 1559 (2004).[12] M. Karlovini and L. Samuelsson, Classical and Quantum Gravity , 4531 (2004).[13] L. B. Szabados, Living Rev. Relativity (2004), URL .[14] R. M. Wald, General Relativity (Chicago University Press, Chicago, 1984).[15] A. Rendall and B. G. Schmidt, Class. Quant. Grav. , 985 (1991).[16] P. Haensel, in Physics of Neutron Star Interiors , edited by D. Blaschke, N. Glendenning, and A. Sedrakian (2001), vol. 578, pp.127–174.[17] T. Strohmayer, H. M. van Horn, S. Ogata, H. Iyetomi, and S. Ichimaru, Astrophysical Journal , 379 (1991).[18] A. Kabobel, Master’s thesis, Universität Tübingen (2001).[19] B. Lukács, Nuovo Cimento
B40 , 169 (1977).[20] B. Carter, personal communication.[21] J. Park, Gen. Rel. Grav. , 235 (2000), URL gr-q /9810010 .[22] R. Penrose and W. Rindler, Spinors and Spacetime , vol. 1 (Cambridge University Press, Cambridge, 1984).[23] In this paper we use geometric units and the conventions of [22].9