Global Solutions of the Compressible Euler-Poisson Equations with Large Initial Data of Spherical Symmetry

We are concerned with a global existence theory for finite-energy solutions of the multidimensional Euler-Poisson equations for both compressible gaseous stars and plasmas with large initial data of spherical symmetry. One of the main challenges is the strengthening of waves as they move radially inward towards the origin, especially under the self-consistent gravitational field for gaseous stars. A fundamental unsolved problem is whether the density of the global solution forms a delta measure ({\it i.e.}, concentration) at the origin. To solve this problem, we develop a new approach for the construction of approximate solutions as the solutions of an appropriately formulated free boundary problem for the compressible Navier-Stokes-Poisson equations with a carefully adapted class of degenerate density-dependent viscosity terms, so that a rigorous convergence proof of the approximate solutions to the corresponding global solution of the compressible Euler-Poisson equations with large initial data of spherical symmetry can be obtained. Even though the density may blow up near the origin at a certain time, it is proved that no delta measure ({\it i.e.}, concentration) in space-time is formed in the vanishing viscosity limit for the finite-energy solutions of the compressible Euler-Poisson equations for both gaseous stars and plasmas in the physical regimes under consideration.


Introduction
We are concerned with the global existence theory for spherically symmetric solutions of the multidimensional (M-D) compressible Euler-Poisson equations (CEPEs) with large initial data. CEPEs govern the motion of compressible gaseous stars or plasmas under a self-consistent gravitational field or an electric field, which take the form: for t > 0, x ∈ R n , n ≥ 3, where ρ is the density, p is the pressure, M ∈ R n represents the momentum, and Φ represents the gravitational potential of gaseous stars if κ > 0 and the plasma electric field potential if κ < 0. When ρ > 0, U = M ρ ∈ R n is the velocity. By scaling, we always fix κ = ±1 throughout this paper; that is, κ = 1 for the gaseous star and κ = −1 for the plasma. The pressure-density relation is p = p(ρ) = a 0 ρ γ , where γ > 1 is the adiabatic exponent. Again, by scaling, constant a 0 > 0 may be chosen to be a 0 = (γ−1) 2 4γ . We consider the Cauchy problem for (1.1) with the Cauchy data: as |x| → ∞, (1.2) subject to the asymptotic condition: In (1.2), the initial far-field velocity has been assumed to be zero in (1.2) without loss of generality, owing to the Galilean invariance of system (1.1). Since a global solution of CEPEs (1.1) normally contains the vacuum states {(t, x) : ρ(t, x) = 0} where the fluid velocity U (t, x) is not well-defined, we use the physical variables such as the momentum M(t, x), or M(t,x) √ ρ(t,x) (which will be shown to be always well-defined globally), instead of U (t, x), when the vacuum states are involved. The global existence for problem (1.1)-(1.3) is challenging, mainly owing to the possible appearance of cavitation and concentration, besides the formation of shock waves, in the solutions, which lead to the lack of higher order regularity of the solutions, so that our main focus has to be finite-energy solutions for CEPEs (1.1). To solve this existence problem, we consider the vanishing viscosity limit of the solutions of the compressible Navier-Stokes-Poisson equations (CNSPEs) with carefully adapted density-dependent viscosity terms in R n : ∂ t M + div M⊗M ρ + ∇p + ρ∇Φ = εdiv µ(ρ)D( M ρ ) + ε∇ λ(ρ)div( M ρ ) , ∆Φ = κρ, (1.4) where D( M ρ ) = 1 2 ∇( M ρ ) + (∇( M ρ )) ⊤ is the stress tensor, the Lamé (shear and bulk) viscosity coefficients µ(ρ) and λ(ρ) depend on the density (that may vanish on the vacuum) and satisfy µ(ρ) ≥ 0, µ(ρ) + nλ(ρ) ≥ 0 for ρ ≥ 0, and parameter ε > 0 is the inverse of the Reynolds number. Formally, as ε → 0+, CNSPEs (1.4) converge to CEPEs (1.1). However, its rigorous mathematical proof has been one of the most challenging open problems in mathematical fluid dynamics; see Chen-Feldman [7], Dafermos [15], and the references cited therein. Many efforts have been made in the analysis of CEPEs (1.1). We focus mainly on some relevant time-dependent problems. Some important progress has been made on the M-D CEPEs with κ = −1 (plasmas) in Guo [25], Guo-Ioescu-Pausader [27], Guo-Pausader [28], and Ionescu-Pausader [36], in which they proved the global existence of smooth solutions around a constant neutral background under irrotational, smooth, and localized perturbations of the background with small amplitude. For the 3-D gaseous stars problem (κ = 1), a compactly supported expanding classical solution was discovered by Goldreich-Weber [23] in 1980; see also [21,57]. Hadzic-Jang [30] proved the nonlinear stability of the Goldreich-Weber solution under small spherically symmetric perturbations for the adiabatic exponent γ = 4 3 , while the problem for γ = 4 3 is still widely open. When the initial density is small and has compact support, Hadzic-Jang [31] constructed a class of global-in-time solutions of the 3-D CEPEs in the Lagrangian coordinates for γ = 1 + 1 k , k ∈ N\{1} or γ ∈ (1, 14 13 ). More recently, Guo-Hadzic-Jang [26] constructed an infinite-dimensional family of collapsing solutions of CEPEs (1.1). We also refer [44,54,56] for the local well-posedness of smooth solutions.
On the other hand, owing to the strong hyperbolicity and nonlinearity, the smooth solutions of (1.1) may break down in finite time, especially when the initial data are large (cf. [11,57]). Therefore, the weak solutions have to be considered for the Cauchy problem with large initial data. For gaseous stars (i.e., κ = 1) surrounding a solid ball, Makino [58] proved the local existence of weak solutions for γ ∈ (1, 5 3 ] with spherical symmetry; also see Xiao [71] for global weak solutions for a class of initial data. For the compressible Euler equations, we refer to [6,8,10,12,13,16,17,32,38,43,50,51] and the references therein. For CNSPEs (1.4), many efforts have also been made regarding the global existence of solutions. For CNSPEs (1.4) with constant viscosity, some global existence results for weak solutions for viscous gaseous stars (i.e., κ = 1) have been obtained; see also [19,37,40,41,55] and the references cited therein. For CNSPEs with density-dependent viscosity terms, Zhang-Fang [73] obtained a unique global weak solution for a spherically symmetric vacuum free boundary problem with γ > 1 for a small perturbation around some steady solution; later, the global existence of spherically symmetric weak solutions was proved by Duan-Li [18] for the 3-D problem for κ = 1 and γ ∈ ( 6 5 , 4 3 ] with stress free boundary condition and nonzero initial density for arbitrarily large initial data. Recently, Luo-Xin-Zeng [52,53] proved the existence and large-time stability of spherically symmetric smooth solutions of the 3-D viscous problem (with κ = 1) for a small perturbation around the Lane-Emden solution for γ ∈ ( 4 3 , 2). For the global existence of solutions of the compressible Navier-Stokes equations, we refer to [20,33,39,49,59] for the case with constant viscosity, [45,69] for the case with density-dependent viscosity, and the references cited therein. In particular, we remark that the BD entropy estimate developed in [4] for the derivative estimate of the density plays a key role in [45,69]. Such an estimate is based on the new mathematical entropy -the BD entropy, first discovered by Bresch-Desjardins [1] for the particular case (µ, λ) = (ρ, 0), and later generalized by Bresch-Desjardins [2] to include any viscosity coefficients (µ, λ) satisfying the BD relation: λ(ρ) = ρµ ′ (ρ) − µ(ρ); also see [3]. The BD-type entropy will also be used in this paper.
The idea of regarding inviscid gases as viscous gases with vanishing physical viscosity can date back to the seminal paper by Stokes [67]; see also the important contributions in [15,35,63,64]. Most of the known results are for the vanishing viscosity limit from the compressible Navier-Stokes to Euler equations. The first rigorous convergence analysis of the vanishing physical viscosity limit from the barotropic Navier-Stokes to Euler equations was made by Gilbarg [22], in which he established the mathematical existence and vanishing viscosity limit of the Navier-Stokes shock layers. For the convergence analysis confined in the framework of piecewise smooth solutions; see [24,34,72] and the references cited therein. For general data, due to the lack of L ∞ uniform estimate, the L ∞ compensated compactness framework [16,17,50,51] does not apply directly for the vanishing viscosity limit of the compressible Navier-Stokes equations. Motivated by this, Chen-Perepelitsa [8] established an L p compensated compactness framework for the whole physical range of adiabatic exponent γ > 1 (also see [43] for γ ∈ (1, 5 3 ) and the references therein), and then applied it to prove rigorously the vanishing viscosity limit of the solutions of the 1-D compressible Navier-Stokes equations to the corresponding finite-energy solutions of the Euler equations for large initial data. Most recently, Chen-Wang [12] established the vanishing viscosity limit of the compressible Navier-Stokes equations with general data of spherical symmetry and obtained the global existence of spherically symmetric solutions of the compressible Euler equations with large data, in which it was proved that no delta measure is formed for the density function at the origin.
For problem (1.1)-(1.3), owing to the additional difficulties arisen from the possible appearance of concentration and cavitation, besides the involvement of shock waves, it has been a longstanding open problem to construct global finite-energy solutions with large initial data of spherical symmetry. The key objective of this paper is to solve this open problem and establish the global existence of spherically symmetric finite-energy solutions of (1.1): ρ(t, x) = ρ(t, r), M(t, x) = m(t, r) x r , Φ(t, x) = Φ(t, r) for r = |x|, (1.5) subject to the initial condition: (ρ, M)(0, x) = (ρ 0 , M 0 )(x) = ρ 0 (r), m 0 (r) x r −→ (0, 0) as r → ∞, (1.6) and the asymptotic boundary condition: Since Φ(0, x) can be determined by the initial density and the boundary condition (1.7), there is no need to impose initial data for Φ. To achieve this, we establish the vanishing viscosity limit of the corresponding spherically symmetric solutions of CNSPEs (1.4) with the adapted class of degenerate density-dependent viscosity terms and approximate initial data of similar form to (1.6). For spherically symmetric solutions of form (1.5), systems (1.1) and (1.4) become and respectively. The study of spherically symmetric solutions can date back to the 1950s and has been motivated by many important physical problems such as stellar dynamics including gaseous stars and supernovae formation [5,65,70]. In fact, the most famous solutions of CEPEs (1.1) are the Lane-Emden steady solutions [5,48], which describe spherically symmetric gaseous stars in equilibrium and minimize the energy among all possible configurations. More precisely, for the 3-D case, there exists a compactly supported steady spherical symmetry solution with finite mass for γ ∈ ( 6 5 , 2). For the time-dependent system, the central feature is the strengthening of waves as they move radially inward near the origin, especially under the self-gravitational force for gaseous stars. The spherically symmetric solutions of the compressible Euler equations may blow up near the origin [14,46,61,70] at certain time in some situations. However, it has not well understood how the spherically symmetric solutions of CEPEs (1.1) with self-gravitational force (which drags the gas particles to the origin) blow up when the initial total-energy is finite. A fundamental unsolved problem is whether a concentration is formed at the origin; that is, the density becomes a delta measure at the origin, especially when a focusing spherical shock is moving inward towards the origin under self-consistent gravitational field.
In this paper, we establish the global existence of finite-energy solutions of problem (1.1)-(1.3) for CEPEs with spherical symmetry as the vanishing viscosity limits of global weak solutions of CNSPEs (1.4) with corresponding initial and asymptotic conditions, which indicates that no delta measure is formed for the density of the solution of problem (1.1)-(1.3) in the limit indeed. To achieve these, the main point is to establish appropriate uniform estimates in L p and the H −1 loccompactness of the entropy dissipation measures for the solutions of CNSPEs (1.9) subject to the corresponding initial and asymptotic conditions. Owing to the possible appearance of cavitation, the singularity of geometric source terms at the origin, as well as the gravitational force for the gaseous star case, the global solutions of CNSPEs (1.9) with large initial data are not smooth in general. Thus, we start with the construction of approximate smooth solutions of the truncated approximate problem (3.1)-(3.6) for CNSPEs (1.4), where the origin and the far-field are cut off, and a stress-free boundary condition is imposed.
In general, we have two basic candidates for the boundary conditions of the approximate problem (3.1)-(3.6): One is to use the Dirichlet boundary conditions: u(t, a) = u(t, b) = 0 as in [12], in which case it is difficult to obtain the higher integrability on the velocity (see Lemma 3.7) due to the far-field vacuum (since the total mass is finite). Another choice is to use the vacuum free boundary condition. To obtain the higher integrability on the velocity, we need to know the trace information for both u and u r ; however, some obstacles arise for deriving the trace estimates for u and u r in the case of the vacuum free boundary problem. One of our main observations is that the stress-free boundary conditions (3.3)-(3.4) we have adapted in §3 serve our purpose to avoid the difficulties mentioned above. Even though, we still have to overcome the following additional difficulties: (i) Owing to the effect of self-gravitational force for κ = 1, we need condition (2.8) in §2 to close the basic energy estimate for γ ∈ ( 2n n+2 , 2(n−1) n ], which implies that the initial total mass can not be too large when the total initial-energy is fixed. The lower bound condition γ > 2n 2+n is essentially used when we deal with the gravitational potential. (ii) To obtain the derivative estimate of the density, we use the BD entropy introduced in [1,4]; also see [2,3]. To close the bound, we need to control the boundary term 1 n p(ρ ε,b 0 (b))b n for the approximate initial data; see (3.58). To solve this problem, we construct the approximate initial data (ρ ε,b 0 , ρ ε,b 0 u ε,b 0 ) so that 1 n p(ρ ε,b 0 (b))b n are uniformly bounded; see (3.12). The details of constructing approximate initial data are given in Appendix A; see Lemmas A.2-A.10. (iii) For the free boundary problem (3.1)-(3.6) below, a follow-up point is whether the free boundary domain Ω T (see (3.2)) will expand to the whole space as b → ∞; otherwise, it would not be a good approximation to the original Cauchy problem. We solve this difficulty by proving that provided b ≫ 1 for any given T . Condition (3.12) is crucial to prove (1.10); see Lemma 3.4 for details. (iv) To utilize the L p compensated compactness framework [8], we still need to have the higher velocity integrability. We use the entropy pairs (η # , q # ) generated by ψ(s) = 1 2 s|s|. Then we have to deal with the boundary term (q # − uη # )(t, b(t))b(t) n−1 . In general, it is impossible to have a uniform bound for both q # (t, b(t))b(t) n−1 and (uη # )(t, b(t))b(t) n−1 . One of our key observations is the cancelation between q # (t, b(t)) and (uη # )(t, b(t)) via observing that (1.11) With the help of the trace estimates in the basic estimates and the BD entropy estimate, it serves perfectly to obtain the uniform trace estimate for the terms on the right-hand side of (1.11). On the other hand, the trace of u r can be handled by using (3.4); see (3.104) for details.
This paper is organized as follows: In §2, we first introduce the notion of finite-energy solutions of problem (1.1)-(1.3) for CEPEs and then state the main theorems of this paper and several remarks. In §3, we first derive some uniform estimates of the solutions of the free boundary problem (3.1)-(3.6) for the approximate CNSPEs. In §4, we establish the global existence of weak finite-energy solutions of (1.4) with large initial data of spherical symmetry and finite-energy. Moreover, some uniform estimates in L p and the H −1 loc -compactness of entropy dissipation measures for the weak solutions of CNSPEs (1.9) are also obtained. In §5, the vanishing viscosity limit of weak solutions of CNSPEs (1.9) is proved by using the compensated compactness framework [8], which leads to a global finite-energy solution of CEPEs (1.1). In Appendix A, we construct the approximate initial data with desired properties.

Mathematical Problem and Main Theorems
In this section, we first introduce the notion of finite-energy solutions of problem (1.1)-(1.3) for CEPEs in R n+1 We assume that the initial data (ρ 0 , M 0 )(x) and corresponding initial potential function Φ 0 (x) have both finite initial total-energy: and finite initial total-mass: where e(ρ) := a 0 γ−1 ρ γ−1 represents the internal energy, and ω n := 2π n 2 Γ( n 2 ) denotes the surface area of the unit sphere in R n .  (ii) For a.e. t > 0, the total energy is finite: • For κ = 1 (gaseous stars), For the case that κ = 1 (gaseous stars), denote M c (γ) as the critical mass given by depending only on (n, γ). We now state the main theorem of this paper. Remark 2.3. By the Poisson equation, the initial condition on ∇ x Φ 0 is indeed a condition on the initial density ρ 0 . In fact, to make the Poisson equation solvable and make sense, we need the additional condition n+2 ) as required in Theorem 2.2. However, for case κ = 1 (gaseous stars), such an additional condition is not required for γ > 2n n+2 . Remark 2.4. To the best of our knowledge, Theorem 2.2 provides the first global-in-time solution of the M-D CEPEs (1.1) with large initial data. For κ = 1 (gaseous stars), condition γ > 2n n+2 (e.g. γ > 6 5 for n = 3) is necessary to ensure the global existence of finite-energy solutions with finite total mass, which corresponds to the one for the Lane-Emden solutions (cf. [5,48]). It is still an open problem for the existence of global weak solutions for gaseous stars when γ ∈ [1, 2n n+2 ] (e.g., γ ∈ [1, 6 5 ] for n = 3). Remark 2.5. For the steady gaseous star problem, it is well-known that there exists no steady white dwarf star with total mass larger than the Chandrasekhar limit M ch when γ ∈ ( 6 5 , 4 3 ]; see [5]. In this paper, for the 3-D time-dependent gaseous star problem with γ ∈ ( 6 5 , 4 3 ], we also need the restriction on the total mass of the gaseous star: M < M c (γ) with M c (γ) defined in (2.8). It is interesting to clarify whether the delta measure could be formed at some time when M < M c (γ) is violated.

Construction and Uniform Estimates of Approximate Solutions
In order to deal with the difficulties for the appearance of cavitation and the singularity at the origin, besides shock waves, as well as uniform estimates of approximate solutions, we construct our approximate solutions via the following approximate free boundary problem for CNSPEs: On the free boundary r = b(t), the stress-free boundary condition is chosen: On the fixed boundary r = a = b −1 , we impose the Dirichlet boundary condition: The initial condition is We always assume that the initial data functions (ρ ε,b 0 , u ε,b 0 )(r) are smooth, and compatible with the boundary conditions (3.4)-(3.5), and 0 When κ = 1, for given total energy E ε,b 0 > 0, similar to (2.8), we define the critical mass: (3.9) For the initial data (ρ ε 0 , m ε 0 ) imposed in (2.9) satisfying (2.13)-(2.18), it follows from Lemma A.10 in Appendix A that there exists a sequence of smooth functions (ρ ε,b (3.13) Property of (3.12) is important for us to close the BD-type entropy estimate in Lemma 3.3 below.
For strong solutions, it is convenient to deal with IBVP (3.1)-(3.6) in the Lagrangian coordinates.
and t ∈ [0, T ], we define the Lagrangian coordinates (τ, x) as Applying the Euler-Lagrange transformation, (3.19) and the fixed boundary x = M ωn corresponds to the free boundary b(τ ) = r(τ, M ωn ) in the Eulerian coordinates.
Lemma 3.1 (Basic energy estimate). Any smooth solution (ρ, u)(t, r) of problem (3.1)-(3.6) satisfies the following energy identity: where ρ(t, r) has been understood to be 0 for r ∈ [0, a] ∪ (b, ∞) in the second term of the right-hand side (RHS) and the second term of the left-hand side (LHS). In particular, the following estimates hold: where the positive constant C γ > 0 is defined as (3.23) where C(M ) > 0 is some positive constant depending only on the total initial-mass M .
Proof. We divide the proof into seven steps.
For the second term of (3.25)-LHS (i.e., the left-hand side of (3.25)), it follows from (3.17) 1 and (3.18)- (3.19) and integration by parts that For the first term of (3.25)-RHS (i.e., the right-hand side of (3.25)), a direct calculation shows that For the last term of (3.25)-RHS, it follows from (3.19 Plugging (3.29) back to the Eulerian coordinates, we obtain To close the estimates, we need to control the terms involving potential Φ. Noting (3.15), a direct calculation shows that (3.32) On the other hand, it follows from (3.15) that which, together with (3.32), yields that where we need to understand ρ to be zero for r ∈ [0, a) ∪ (b, ∞) in the last equality of (3.33).
For the case that κ = 1 with γ = 2(n−1) n , i.e., n−2 n(γ−1) = 1, we use (3.31) and (3.37) to obtain and where we used Thus, due to the continuity of which contradicts (3.45). Therefore, (3.42) always holds under condition (3.13). Now, under condition (3.13), it follows from (3.42) that For later use, we analyze the behavior of density ρ on the free boundary. It follows from (3.17) 1 and (3.18) that Then we obtain In the Eulerian coordinates, it is equivalent to the form: The density behavior on the free boundary (3.49) is important, which will be used frequently later.

Lemma 3.3 (BD-type entropy estimate).
Under the conditions of Lemma 3.1, for any given T > 0, the following holds for any t ∈ [0, T ]: Proof. We divide the proof into four steps.
Proof. Noting the continuity of b(t), we first make the a priori assumption: A direct calculation by using (3.49), (3.68), and Lemma 3.1 yields that Thus, we have closed our a priori assumption (3.68). Then, using (3.72) and the continuity arguments, we conclude (3.67). Multiplying (3.1) 2 by w(r), we have (ρuw) t + (ρu 2 + p(ρ))w r + n − 1 r Multiplying (3.75) by ρw and performing a direct calculation, we have 4. For I 2 , integrating by parts, we have which yields that To close the estimates for I 1 and I 2 , we need to bound the last term of both (3.80)-RHS and (3.81)-RHS. There are three cases: To use the compensated compactness framework in Chen-Perepelitsa [8], we still need to obtain the higher integrability on the velocity. For this, we require to exploit several important properties of some special entropy pairs. First, taking ψ(s) = 1 2 s|s| in (2.11), then the corresponding entropy and entropy flux are represented as (3.89) A direct calculation shows that where and whereafter C γ > 0 is a universal constant depending only on γ > 1. We regard η # as a function of (ρ, m) to obtain It is direct to check that From [8,9], we know that The following lemma is crucial to control the trace estimates for the higher integrability on the velocity. In fact, we have the boundary parts (uη # )(t, b(t)) and q # (t, b(t)), and it is impossible to have the uniform trace bound (independent of ε) for each of them. Our key point is to identify the cancelation between these two boundary parts. (3.93) Proof. It follows from (3.89) that (3.94) A direct calculation shows that For I 2 , we note that I 2 = 0 if u = 0. Thus, it suffices to consider u = 0. We divide the proof into three cases.
For the second term of (3.101)-RHS, we have For the third term of (3.101)-RHS, For the last term of (3.101)-RHS, it follows from (3.90) that

Existence of Global Weak Solutions of CNSPEs
In this section, for fixed ε > 0, we take limit b → ∞ to obtain global weak solutions of CN-SPEs with some uniform bounds, which are essential for applying the compensated compactness framework in §5 below. We often denote the solutions of (3.1)-(3.6) as (ρ ε,b , u ε,b ) for simplicity of presentation in this section, since To take the limit, we have to be careful, since the weak solutions may involve the vacuum. We use the similar compactness arguments as in [29,60]  Then it is direct to check that the corresponding vector function U ε,b = (ρ ε,b , M ε,b , Φ ε,b ) is a classical solution of CNSPEs: Lemma 4.1. For fixed ε > 0, there exists a function ρ ε (t, r) such that, as b → ∞ (up to a subsequence), a.e. and strongly in C(0, T ; L q loc ), (4.4) for any q ∈ [1, ∞), where L q loc denotes L q (K) for K ⋐ (0, ∞). Proof. It follows from Lemmas 3.1 and 3.3 that Using the mass equation (3.1) 1 and Lemma 3.1, we have , which, together with the Aubin-Lions lemma, yields that ρ ε,b is compact in C(0, T ; L q loc ) for any q ∈ [1, ∞). Notice that, for any K ⋐ (0, ∞) and b 1 , b 2 ∈ (1, ∞), where C T,K > 0 is a constant independent of b 1 and b 2 . Then there exists a function ρ ε (t, r) such that, as b → ∞ (up to a subsequence), ( ρ ε,b , ρ ε,b ) → ( √ ρ ε , ρ ε ) a.e. and strongly in C(0, T ; L q loc ) for any q ∈ [1, ∞). Then (4.4) follows.

Moreover, it follows from Lemma 3.7 and Fatou's lemma that
Next, since (ρ ε,b , m ε,b ) converges almost everywhere, it is direct to know that sequence Moreover, for any given positive constant k > 0, we have (4.14) For k ≥ 1, we cut the L 2 -norm as It is direct to know that ρ ε,b u ε,b I |u ε,b |≤k is uniformly bounded in (ε, b) in L ∞ (0, T ; L q loc ) for all q ∈ [1, ∞). Then it follows from (4.14) that Using (4.13), we have Substituting (4.16)-(4.17) into (4.15), we see that Thus, the result follows by taking k → ∞.
We show that (ρ ε , M ε , Φ ε ) is a global weak solution of the Cauchy problem for CNSPEs (1.4) in R n in the sense of Definition 2.6.

H −1 loc -Compactness.
To use the compensated compactness framework in [8], we need the H −1 loc -compactness of entropy pairs. Lemma 4.12 (H −1 loc -compactness). Let (η, q) be a weak entropy pair defined in (2.11) for any smooth compact supported function ψ(s) on R. Then, for ε ∈ (0, ε 0 ], Proof. To obtain (4.73), we have to make the argument in the weak sense, since (ρ ε , M ε , Φ ε ) is a weak solution of CNSPEs (1.4). In fact, we first have to study the equation for ∂ t η(ρ ε , m ε ) + ∂ r q(ρ ε , m ε ) in the distributional sense, which is more complicated than that in [8,9]. We divide the proof into five steps.
Remark 4.13. Since the weak solution (ρ ε , m ε ) of CNSPEs may contain vacuum states, and the estimate of u ε r is not available, it is difficult to use directly the weak form of CNSPEs to study the H −1 loc -compactness. In fact, we first have to study the equation satisfied by ∂ t η ε + ∂ r q ε in the distributional sense, and then use the equation to prove the compactness result.
Finally, we consider Poisson's equation. Let ξ(x) ∈ C 1 0 (R n ) be any smooth function with compact support. For any t 2 > t 1 ≥ 0, we use (4.72), (5.7), and similar arguments as in (5.17), and then pass limit ε → 0+ (up to a subsequence) to obtain Applying the Lebesgue point theorem, we obtain that, for a.e. t ≥ 0, Lemma A.1 (Sobolev's Inequality). For n ≥ 3, let ∇f ∈ L 2 (R n ) and lim |x|→∞ f (x) = 0. Then where A n is the best constant which is given by with ω k = 2π k 2 Γ( k 2 ) as the surface area of unit sphere in R k .
To keep the L p -properties of mollification, it is more convenient to smooth out the initial data in the original coordinates in R n ; so we do not distinguish functions (ρ 0 , m 0 )(r) from (ρ 0 , m 0 )(x) = (ρ 0 , m 0 )(|x|) for simplicity below.
From now on, we denote C > 0 is a universal constant independent of ε, δ, and b.