A complete characterization of sharp thresholds to spherically symmetric multidimensional pressureless Euler-Poisson systems

A general system of pressureless Euler-Poisson (EP) equations has the following form,

with smooth initial data $(\rho_{0}\geq 0,\mathbf{u_{0}})$. The constant parameters $k,c\geq 0$ are the forcing coefficient and background state respectively. The sign of the forcing coefficient $k$ signifies the type of particles being modeled and its magnitude gives a measure of the strength between them. When the force in-between particles is repulsive, for example in the case of charge flow, then $k>0$. Within the pressureless setup (1.1), $k<0$ is relevant in the case of interstellar clouds where the pressure gradient becomes negligible compared to the gravitation forces, see [11]. In the pressureless setup with same charge particles ($k>0$), the background state is, in practicality, a profile, that is, a function of the spatial variable, $c=c(x)$, see [13]. The background models the doping profile for charge flow in semiconductors. However, the sheer complicacy of the system has restricted researchers to consider the background as a constant. To our knowledge, [3] is the only work that considers background as a profile in obtaining critical thresholds.

A locally well-posed PDE system exhibiting critical threshold phenomena is the one wherein existence of global-in-time solutions is dependent on whether the initial data crosses a certain threshold manifold. This critical threshold manifold divides the phase space of initial data into two mutually exclusive regions or sets. If the initial data lies completely in one of the sets (also called subcritical region), there is global solution. However, if some part of initial data lies outside this set, or in other words, in the supercritical region, a breakdown occurs and solution loses its smoothness in finite time.

As one can expect, the question of proving the existence and subsequently finding the critical threshold manifold is simpler in one dimension ($N=1$). In this case, the threshold is a curve on the $(u_{0x},\rho_{0})$ plane and the subcritical region is given by,

A vast amount of literature exists for critical thresholds to (1.1) and other similar systems. The existence of such curve was first identified and analyzed in [9] for EP systems. The authors analyzed the one dimensional and multidimensional (with spherical symmetry) cases. A series of works then followed for EP as well as other systems, [2, 3, 4, 5, 7, 12, 16, 17, 18, 19, 20, 21, 22, 24, 25, 26, 27, 28, 29] among many others. As is known for general hyperbolic conservation laws that for any initial data, singularities develop as one moves forward in time, [15]. This is due to the convective forces present. Addition of source terms however, can result in a set of initial data that leads to global solutions. To extract a set of initial data (subcritical region) for which there is global well-posedness, is an interesting territory to be explored. The ‘goodness’ of the external forces can balance or even outweigh the ‘bad’ convective forces and result in a ‘large’ subcritical region. An example where the external forces completely obliterate the convective forces is the strongly singular Euler-Poisson-alignment (EPA) system. In [14], the authors show that the system (1.1) with $N=1$ and an additional nonlocal alignment force in the momentum equation results in global-in-time solutions for any initial data. This result is remarkable and in stark opposition to the conventional results for EP systems.

For threshold results in one dimensional EP and EPA systems, one can refer to [2, 3, 4, 5, 8, 14, 26, 28]. In [26], the authors include pressure and exclude background. It is quite difficult to obtain thresholds for the full EP system (with pressure) due to the strict hyperbolicity resulting in two characteristic flow paths. Critical threshold for EP systems with pressure and $c>0$ is largely an open problem.

For pressureless EP systems, the key is to obtain bounds on the gradient of velocity. Owing to the local existence result, global-in-time solutions to (1.1) can be ‘extended’ upon a priori smooth local solutions. The local existence for such systems is well-known, see [23] for a proof. We state the result here.

As a result, in 1D, one needs to bound $|u_{x}|$ to guarantee global well-posedness. However, things are different in higher dimensions where the gradient of velocity is in fact an $N\times N$ matrix. One finds out that an ODE along the characteristic path can be obtained for the divergence of the gradient matrix. However, analyzing the divergence is not enough. One needs to control the spectral gap (sum of all the absolute differences of eigenvalues), see for example [12], in order to guarantee smooth solutions for all time. This is the tricky part. Given the inherent property of hyperbolic balance laws, it is in general easier to obtain sufficient conditions for breakdown of solutions as compared to conditions for global existence. In [6, 7], the authors obtain bounds on the supercritical region for (1.1) for $k<0$.

To obtain bounds on subcritical region is more involved. Several workarounds have been used by various researchers by simplifying the EP system. In [16, 21], the authors study and obtain thresholds for the restricted EP equations in two and three dimensions respectively.

Another derivative of the EP system that has been extensively studied is (1.1) with spherical symmetry. It can be shown that if the initial data to (1.1) is such that $\rho(0,\mathbf{x})=f(|\mathbf{x}|),\mathbf{u}(0,\mathbf{x})=g(|\mathbf{x}|)\frac{\mathbf{x}}{|\mathbf{x}|}$ for smooth functions $g,f$, then for any fixed time and as long as the solution is smooth, the symmetry extends, that is, the solutions are $\rho(t,|\mathbf{x}|),u(t,|\mathbf{x}|)\frac{\mathbf{x}}{|\mathbf{x}|}$. Applying this simplification to (1.1), we obtain the following system of equations,

with $r=|\mathbf{x}|>0,N\geq 2$ and subject to smooth initial density and velocity

Even though the system is now simplified as compared to the general EP system, it turns out that it is still tricky to analyze for thresholds, especially when $k>0$. (1.2) was first studied in [9] with $c=0$ and expanding flows ($u_{0}>0$). Sufficient conditions on global existence and finite-time-breakdown were obtained for $N=2,3$ and a sharp condition was obtained for $N=4$. Later, a sharp condition was derived in [29] however, once again, only for zero background and expanding flows. A sufficient blow-up condition was derived in [30] for $c=0$, however, no conclusion was made with regards to global existence.

The dynamics of (1.2) is quite different for $c>0$ compared with $c=0$. One major difference is that in the former, the density is a perturbation around the constant background. The mass is infinite and the following holds,

which is different from the zero background state where the mass is finite and conserved. A major step forward to mitigate the restrictiveness of the past results was taken in [28]. The author there reduced (1.2) to a $4\times 4$ ODE system along a characteristic path, and subsequently proved the existence of critical thresholds. This $4\times 4$ ODE system derived by the author is crucial in our analysis as well. An intriguing discovery by the author was that the Poisson forcing in (1.2b) is enough to avoid concentration at the origin, even if there are points where the initial velocity points towards the origin. To our knowledge, it is the first work with regards to thresholds of spherically symmetric EP system wherein the assumption of expanding flows was dropped. It was also noted that $N=2$ is the critical case and the analysis needs to be done fairly differently as compared to the case when $N\geq 3$. The author obtained partial results for thresholds for $c=0$, that is, bounds on the subcritical and supercritical regions were presented. However, no results were presented for $c>0$ citing the highly chaotic dynamics in such a scenario.

In this paper, we present sharp thresholds for both $c>0$ as well as $c=0$. Several steps are needed to arrive at the precise thresholds. Special care has to be taken with regards to finding the subcritical region for $c>0$. We essentially start out by characterizing the supercritical region and narrowing down the possible initial configurations that might allow for global solutions, eventually extracting out the subcritical region.

Our main results can be stated non-technically as follows:

• •

  For the nonzero background case ($c\not\equiv 0$), we show that the EP system admits a global-in-time smooth solution if and only if the initial data is smooth, lying in a subcritical region, $\Theta_{N}$. Theorem 2.1 contains such critical threshold result. The explicit definition of $\Theta_{N}$ is stated after the result.

• •

  For the zero background case ($c=0$), we show that the EP system admits a global-in-time smooth solution if and only if the initial data is smooth, and lying in certain subcritical region, $\Sigma_{N}$. Theorem 2.6 contains the precise thresholds for dimensions greater than or equal to three. Theorem 2.8 contains the threshold results in dimension two. The explicit definition of $\Sigma_{N}$ is stated after each of the results.

We discover a completely novel nonlinear quantity that helps to analyze the system. At $t=0$, this is given by

In fact, such quantity is very crucial in analyzing and simplifying the representation of subcritical/supercritical regions. The full motivation and usage of this quantity will be thoroughly discussed in the Sections to come. It is found that a simpler breakdown condition can be obtained using this expression, the corresponding result for dimensions greater than or equal to three is provided in Theorem 2.4.

This section is devoted to stating precise versions of our main results. We start by specifying that for (1.2), the finite-time-breakdown is manifested in the following form,

for some $r_{c}>0$ and finite $t_{c}$. The first blowup (shock formation) can be concluded from Theorem 1.1. The density blowing up at the same time/position is a consequence of weak hyperbolicity of the pressureless Euler-Poisson systems and is a well-known phenomena for such systems. In the next section, our analysis will allow us to conclude the same. Now we state our results and follow it up with an interpretation. As mentioned in the previous section that the $N=2$ case is critical and needs to be handled separately. Consequently, the representation of subcritical regions is also different for $N=2$ than those for $N\geq 3$. In this section, we also provide the subcritical region sets.

We will describe the subcritical region set with the help of the quantity (1.3). Motivated by this, we set, for a point $(\alpha,x,y,z,\omega)$,

• •

  For $a\in\frac{k}{N-2}\left(\left(\frac{y_{m,N}}{\frac{c}{N}-\frac{y}{\alpha}}\right)^{\frac{N-2}{N}}-1,\left(\frac{y_{M,N}}{\frac{c}{N}-\frac{y}{\alpha}}\right)^{\frac{N-2}{N}}-1\right),$

  (2.2)

  $\displaystyle\begin{aligned} &\frac{1}{\max\{\eta_{1}(0),\eta_{2}(0)\}}<\omega<\frac{1}{\min\{\eta_{1}(0),\eta_{2}(0)\}},\quad\text{for }x\neq 0,\\ &z\in\frac{ky}{\alpha a}\left(\frac{d\eta_{1}}{dt}(0),\frac{d\eta_{2}}{dt}(0)\right),\quad\text{for }x=0.\end{aligned}$

• •

  For $a\in\frac{k}{2}\left(\ln\left(\frac{y_{m,2}}{\frac{c}{2}-\frac{y}{\alpha}}\right),\ln\left(\frac{y_{M,2}}{\frac{c}{2}-\frac{y}{\alpha}}\right)\right),$

  (2.3)

  $\displaystyle\begin{aligned} &\frac{1}{\max\{\eta_{1}(0),\eta_{2}(0)\}}<\omega<\frac{1}{\min\{\eta_{1}(0),\eta_{2}(0)\}},\quad\text{for }x\neq 0,\\ &z\in\frac{ky}{\alpha a}\left(\frac{d\eta_{1}}{dt}(0),\frac{d\eta_{2}}{dt}(0)\right),\quad\text{for }x=0.\end{aligned}$

Here, $0<y_{m,N}<\frac{c}{N}-\frac{y}{\alpha}<y_{M,N}$ are the roots of the equation,

and

with

The explicit expression for $R_{N}(\cdot,\cdot)$ is in (1.6), also in (4.7). The functions

are defined via a second order linear IVP and are positive for all except at infinitely many discrete times.

We now give an interpretation of Theorem 2.1. Similar to existing threshold results, Theorem 2.1 enables us to check pointwise whether or not the given initial data lies in the subcritical region. Not only this, we can in fact construct the subcritical region piece by piece. Our result enables us to divide the $\mathbb{R}^{5}$ space of $(r,u_{0},\phi_{0r},u_{0r},\rho_{0})$ into a $3D$ space and a plane, that is, $(r,u_{0},\phi_{0r})$ (say $P_{1}$) and $(u_{0r},\rho_{0})$ (say $P_{2}$). For each point in $P_{1}$, we can construct the subricitical region on $P_{2}$ through $A_{0}$. Let $(\alpha,x,y):=(r,u_{0},\phi_{0r})$ be any point in $P_{1}$. Note that $y_{m,N},y_{M,N}$ can now be calculated using $\alpha,x,y$ and hence, one can now explicitly find the interval given in Definition 2.2. For any $a$ in that interval, (2.1) outputs a line in $P_{2}(z,\omega)$. In other words, by fixing an $a$, we have a linear relation between $(u_{0r},\rho_{0})$. For this value of $a$, we can find the two constants, $\eta_{1}(0),\eta_{2}(0)$. Indeed for fixed $\alpha,x,y,a$, the functions $\eta_{i}$ are known functions. Finally, (2.2) gives the desired portion of the line on $P_{2}$ that forms part of the subcritical region, as and according to whether $x$ is zero or not. We can do this for all $a$ in the above interval and obtain the complete part of the subcritical region for the fixed point $(\alpha,x,y)$ in $P_{1}$. Carrying out this procedure for each point in $P_{1}$ gives the entire subcritical region.

The next two results give a complete picture of the zero background case with dimension greater than or equal to three.

In the next theorem, we give a full picture of what happens for initial data which satisfies $A_{0}(r)>-\frac{k}{N-2}$ for all $r>0$.

The quantity $a$ in the Definitions 2.7, 2.9 below is the same as in (2.1).

• •

  For $a\in\left(-\frac{k}{N-2},0\right)$,

  (2.4)

  $\displaystyle\begin{aligned} &z>-\frac{ky}{x}+\frac{\alpha a}{x\eta_{1}(0)},\quad\text{for }x\neq 0,\\ &z>\frac{ky}{\alpha a}\frac{d\eta_{1}}{dt}(0),\quad\text{for }x=0.\end{aligned}$

• •

  For $a=0$,

  (2.5)

  $\displaystyle\begin{aligned} &z=-\frac{ky}{x},\ \omega>\frac{1}{\eta_{2}(0)},\quad\text{for }x<0,\\ &z=-\frac{ky}{x},\ \omega>0,\quad\text{for }x>0.\end{aligned}$

• •

  For $a\in\left(0,\frac{k}{N-2}\left(\left(\frac{-\alpha y_{N}}{y}\right)^{1-\frac{2}{N}}-1\right)\right)$,

  (2.6)

  $\displaystyle\begin{aligned} &-\frac{ky}{x}+\frac{\alpha a}{x\eta_{1}(0)}<z<-\frac{ky}{x}+\frac{\alpha a}{x\eta_{2}(0)},\quad\text{for }x<0,\\ &\omega>0,\quad\text{for }x>0.\end{aligned}$

• •

  For $a\in\left[\frac{k}{N-2}\left(\left(\frac{-\alpha y_{N}}{y}\right)^{1-\frac{2}{N}}-1\right),\infty\right)$,

  (2.7)

  $\displaystyle\begin{aligned} &\omega>0,\quad\text{for }x>0.\end{aligned}$

Here, $0<-\frac{y}{\alpha}<y_{N}=y_{N}\left(\frac{u_{0}(\alpha)}{\alpha},\frac{\phi_{0\alpha}(\alpha)}{\alpha}\right)=\left(\frac{(N-2)}{2k}\right)^{\frac{N}{N-2}}\left(R_{N}\left(\frac{u_{0}(\alpha)}{\alpha},\frac{-\phi_{0\alpha}(\alpha)}{\alpha}\right)\right)^{\frac{N}{N-2}}$. The explicit expression for $R_{N}(\cdot,\cdot)$ is in (4.7) (with $c=0$). The functions

are defined later via linear IVPs.

• •

  For $a\in\left(-\infty,0\right)$,

  (2.8)

  $\displaystyle\begin{aligned} &z>-\frac{ky}{x}+\frac{\alpha a}{x\eta_{1}(0)},\quad\text{for }x\neq 0,\\ &z>\frac{ky}{\alpha a}\frac{d\eta_{1}(0)}{dt},\quad\text{for }x=0.\end{aligned}$

• •

  For $a=0$,

  (2.9)

  $\displaystyle\begin{aligned} &z=-\frac{ky}{x},\ \omega>\frac{1}{\eta_{2}(0)},\quad\text{for }x<0,\\ &z=-\frac{ky}{x},\ \omega>0,\quad\text{for }x>0.\end{aligned}$

• •

  For $a\in\left(0,\frac{k}{2}\ln\left(\frac{-\alpha y_{2}}{y}\right)\right)$,

  (2.10)

  $\displaystyle\begin{aligned} &-\frac{ky}{x}+\frac{\alpha a}{x\eta_{1}(0)}<z<-\frac{ky}{x}+\frac{\alpha a}{x\eta_{2}(0)},\quad\text{for }x<0,\\ &\omega>0,\quad\text{for }x>0.\end{aligned}$

• •

  For $a\in\left[\frac{k}{2}\ln\left(\frac{-\alpha y_{2}}{y}\right),\infty\right)$,

  (2.11)

  $\displaystyle\begin{aligned} &\omega>0,\quad\text{for }x>0.\end{aligned}$

Here, $0<-\frac{y}{\alpha}<y_{2}=y_{2}\left(\frac{u_{0}(\alpha)}{\alpha},\frac{\phi_{0\alpha}(\alpha)}{\alpha}\right)=e^{\frac{R_{2}\left(\frac{u_{0}(\alpha)}{\alpha},\frac{-\phi_{0\alpha}(\alpha)}{\alpha}\right)}{k}}$. The explicit expression for $R_{2}(\cdot,\cdot)$ is in (4.7) (with $c=0$).

This paper is arranged as follows. Section 3 entails some preliminary calculations showing how to reduce the full dynamics to a weakly coupled system of four ODEs along characteristics. Section 4 is devoted to the analysis of the system of two of the four ODEs, that are decoupled from the other two, obtaining a closed form of trajectory curves with related properties established for later use. Section 5 entails conditions under which the velocity gradient/density blow up. Section 6 is devoted to constructions of the precise subcritical region leading to global solutions, proving Theorem 2.1. The zero-background case is analysed in Section 7. Finally concluding remarks are given in Section 8.

From Theorem 1.1, we note that to extend the local solutions, we need to obtain bounds on the gradient of velocity. The following result, that was derived by the authors in [25], designates the quantities that need to be bounded for (1.2).

Though the proof has already been given by the authors in [25], we include it here for the sake of completion.

Through this Lemma, we conclude that the Hessians of radial functions are diagonalizable using the same matrix of eigenvectors. From Theorem 1.1, we know that for solutions to persist for all time, the gradient matrix of velocity should be bounded. Given the velocity is radial vector field, we have the following relation for $\mathbf{u}$, $u$ in (1.1b), (1.2b) respectively,

Hence, the two eigenvalues of $\nabla\mathbf{u}$ will be $u_{r}$ and $u/r$. Therefore, to ensure that the gradient of velocity is bounded, we need to control the two quantities: $u_{r},u/r$.

Inspired by this, we will obtain ODE along the characteristic path,

for the desired quantities. Rearranging (1.2a), we obtain,

Note that by (1.2c),

Taking spatial derivative of (1.2b) and using the above equation we obtain,

Next using (1.2b),

Upon integrating (1.2c) from $0$ to $r$,

Local well-posedness requires a boundary condition at the origin. The convention is to assume zero boundary conditions, that is, the above limit is zero. Similarly the density flux ($\rho ur^{N-1}$) approaches zero, which signifies that there is no loss of material at the origin. Taking time derivative of (3.5),

As a consequence, we get $\phi_{rt}=\rho u$. We use this and the previously obtained expression for $\phi_{rr}$ in the below calculation.

Set

Equations (3.2), (3.3), (3.4) and (3.6) can be used to obtain an ODE system,

with initial data $\rho_{0}:=\rho(0,\beta),p_{0}:=p(0,\beta),q_{0}:=q(0,\beta),s_{0}:=s(0,\beta)$ respectively. Here, ^′ indicates differentiation along the path (3.1). Note that we use the notation $\rho$ for the density as in (1.2a), which is a function of time and the spatial variable, as well as for the solution to the ODE (3.8a), which is another function of time only (for a fixed parameter $\beta$). Similarly for the notation $\rho_{0}$. However, it will be clear from context whether we are referring to $\rho$ as the solution of the main PDE system (1.2) or as the unknown in (3.8a).

Note that (3.8c), (3.8d) are decoupled from (3.8a), (3.8b). In fact, a closed form for the trajectory curve for equations (3.8c), (3.8d) can be obtained. We have the following Lemma which describes the behavior of the system (3.8c),(3.8d).