Shifted shock formation for the 3D compressible Euler equations with damping and variation of the vorticity

For the system(1.1), there are many studies both in $1D$ and multi-dimensions. For the one-dimensional Euler equations with damping, the global existence of smooth solutions with small data was proved by Nishida[21] and Slemrod[27] showed that for small data ($L^{\infty}$ sense), the $1D$ Euler equations with damping admit a global smooth solution while for large data, the equations can develop a shock in finite time. Later, these results were generalized by many authors, see [11, 22, 18, 17, 16, 23] and the references therein.
     In multi-dimensional case, the global existence and $L^{2}$ estimates for the solutions to the $3D$ isentropic Euler system with damping was obtained by Kawashima[13]. Later, Sideris-Thomases-Wang in[26] showed that the size of the smooth initial data plays a key role for the lifespan to the $3D$ isentropic Euler equations with damping $a\rho u$. If the initial data are small, then damping can prevent the development of singularities; while if the initial data are large, the damping is not strong enough to prevent the formation of singularities in finite time. However, the authors in[26] only obtained the finite lifespan of the solutions and did not show the shock formation (including the shock mechanism) in finite time. The long time behavior of the solutions was obtained and generalized by many authors, see [29, 12, 24, 28, 26] and the references therein.
     In [6], Christodoulou achieved a significant breakthrough in understanding the shock mechanism for hyperbolic systems in three dimensions. Specifically, he considered the classical, non-relativistic, compressible Euler’s equations in three spatial dimensions, taking the data to be irrotational and isentropic. By imposing certain conditions on the initial data (i.e., the short pulse assumption), he obtained a complete geometric description of the maximal classical development. Notably, he provided a detailed analysis of the behavior of the solution at the boundary of the domain of the maximal classical solution, including a comprehensive description of its geometry. This work represents a significant advancement in understanding of shock formation in hyperbolic systems in three dimensions.
     In [25], Yu-Miao applied Christodoulou’s framework to the quasilinear wave equation in three spatial dimensions. By constructing a family of short pulse initial data that was first introduced by Christodoulou in [5], they demonstrated how the solution breaks down near the singularity and gave a sufficient condition on the initial data which leads to the shock formation in finite time. Additionally, J.Speck and J.Luk applied Christodoulou’s framework to the 2D Euler system without the assumption of irrotation. In [14], they studied plane-symmetric initial data with short pulse perturbation. For such initial data, they showed the shock mechanism such that the first derivatives of $u$ and $\rho$ blow up while $u$ and $\rho$ remain bounded near the shock. Then, in [15], they generalized these results to the 3D case. Specifically, they considered the 1D Euler equation of a simple small-amplitude solution as a plane-symmetric solution in 3D. They perturbed this solution in the $(x_{2},x_{3})$ directions as a nearly plane-symmetric initial data for 3D isentropic Euler equations. They proved that the shock formation mechanism is stable under small and compactly supported perturbations with non-trivial vorticity and provided a precise description of the first singularity.
     Recently, Buckmaster, Shkoller, and Vicol published several results regarding shock formation in the multi-dimensional Euler system. Specifically, in their paper [3], they studied the 2D isentropic Euler equations under azimuthal symmetry (which differs from the 1D problem) with smooth initial data of finite energy and nontrivial vorticity. By utilizing modulated self-similar variables, they were able to obtain point shock forms in finite time, with explicit computation of the blow-up time and location. Furthermore, the solutions near the shock exhibit a cusp type.
     In a subsequent paper, [2], Buckmaster, Shkoller, and Vicol generalized their results to the 3D isentropic case for the ideal gas without any symmetry assumptions. In addition to the findings reported in [3], they were able to demonstrate the precise direction of blow-up and the geometric structure of the tangent surface of the shock profile. They also provided homogeneous Sobolev bounds for the fluid variables with an initial datum of large energy. Later, in[4], the authors extended their results to the full Euler equations. In this work, they primarily investigated the evolution and creation of vorticity, demonstrating that the vorticity remains bounded up to the shock formation. Notably, one of the primary differences between the isentropic and non-isentropic cases is the baroclinic torque in the equation of vorticity, which arises from the interaction between sound waves and entropy waves. Additionally, they constructed a set of irrotational data that results in the instantaneous creation of vorticity, which remains non-zero up to the shock.
     It is important to note that the point shock in the three aforementioned results is stable. That is, for any small, smooth, and generic perturbation of the given initial data, the corresponding Euler system results in a smooth solution that blows up in a small neighborhood of the original shock time and location. However, in [1], the authors demonstrated the existence of an open set of initial data that leads to the formation of an unstable shock. The primary difference between stable and unstable shocks is the set of background solutions for the self-similar variable $W$, specifically the solutions of the various self-similar Burgers equations (3.43). For further details, one may refer to [9].
     The major tool used in the works mentioned above is the method of self-similar coordinates. This method was first introduced by Y.Giga and R.Kohn in [10] to study the asymptotic behavior of the solution to the nonlinear heat equation near the point of singularity. Y.Giga and R.Kohn considered the following nonlinear heat equation

$$u_{t}-\Delta u-|u|^{p-1}u=0,$$

(1.4)

where $p>1$ and $(x,t)\in\mathbb{R}^{n}\times(-1,0)$. They aimed to show the behavior of the solution near the singularity. To this end, they proposed the following self-similar transformation²²2Note that the transformation degenerates at $t=0$.:

$$y=e^{\frac{1}{2}s}x,\hskip 14.22636pts=-\ln(-t),\hskip 14.22636ptw(y,s)=e^{-\frac{1}{2(p-1)}s}u(x,t),$$

(1.5)

which transforms (1.4) into

$$w_{s}-\Delta w+\frac{1}{2}y\cdot\nabla_{y}w+\frac{1}{p-1}w-|w|^{p-1}w=0.$$

(1.6)

This transformation is motivated by the scaling property of (1.4). That is, for any $u$ solves (1.4), $u_{\lambda}:=\lambda^{\frac{2}{p-1}}u(\lambda x,\lambda^{2}t)$ solves (1.4) as well. Based on the analysis of (1.6), they were able to demonstrate the asymptotic behavior of $u$ near the blow-up point $(0,0)$. Later, such method has been applied to various other dynamical systems, including the Sch$\ddot{o}$dinger equation [19], the Prandtl equations [7], the transverse viscous Burgers equation [8], and the semilinear wave equation [20]. The method of self-similar variables can provide precise information about the singularity of a given system by adding modulation variables to enforce certain orthogonality conditions and track the position of the singularity. It can also be used to characterize the type of singularities present, such as stable fixed points or center manifold and etc…, based on the behavior of the solutions near the singularity. For example, the stable shock in [3], [2], and [4] is based on the fact that the solutions approach the background solution near the shock exponentially with respect to the self-similar variables. For further details, one may refer to [9].

In section 2, the fundamental ideas of [2] and [4] will be illustrated by applying them to the standard Burgers equation and the Burgers equation with damping. Our results will demonstrate that the modulated variables can accurately track the location and time of the shock, and that the damping effect only shifts the blow-up time while leaving the blow-up location unchanged compared with the undamped case.
     Section 3 is devoted to reformulating of the Euler system (1.1) into its modulated self-similar version (3.28) by adapting the framework of [4]. Additionally, we extend the solution of the self-similar Burgers equation, which is introduced in Section 2, to three dimensions and derive higher-order derivatives of the evolution equations for the self-similar variables $(W,Z,A,K)$.
     In Section 4, we present an explicit construction of the initial datum in [4]. Note that the initial data for variable $W$ is dependent on the background solution $\bar{W}$, which is a refined solution of the self-similar Burgers equation. Furthermore, we outline the self-similar Bootstrap assumptions for various variables and their derivatives, including the modulation variables and the self-similar variables $(W,Z,A,K)$. These assumptions are more restrictive than the initial data constructed earlier and will be recovered in the subsequent sections. Finally, we state the main theorem (Theorem 4.1) of the paper at the end of this section.
     In Section 5, we present multiple estimates based on bootstrap assumptions, including the estimation of damping and forcing terms. Additionally, we state the $\dot{H}^{m}$ energy bounds for $(W,Z,A,K)$, whose proof solely relies on the bootstrap assumptions and the Friedrich’s energy estimates for the hyperbolic system. These estimates lead to higher-order derivatives estimates for $(W,Z,A,K)$, which in turn help to further refine the estimates for the forcing terms.
     Section 6 is a key aspect of our work. In this section, we recover the bootstrap assumptions for the modulation variables under 10 constraints regarding $\partial^{\gamma}W(0,s)$, where $|\gamma|\leq 2$. These constraints are initially satisfied and remain valid over time as long as the estimates for the modulation variables hold. In addition, it will be clearly shown that how the damping effect can influence the information of shock within this section.
     In Section 7, we derive lower bounds for the Lagrangian trajectories and introduced the crucial lemma7.2, which provides the key estimate required for recovering the remaining bootstrap assumptions. Moreover, we investigate the evolution of the vorticity, which is concentrated on the non-blow up direction, and demonstrated how it is affected by the damping effect. Inspired the work in[4], it will be shown that the initial region which genuinely affects the shock formation.
     By following the estimates presented in sections 5 to 7, we are able to recover the bootstrap assumptions for $(W,Z,A,K)$ as outlined in section 8. This recovery is based on the significant lemma 7.2 and the utilization of the weighted framework (8.9). Note that to recover the bootstrap assumption for $\partial_{1}A_{\nu}$, one has to establish the relation between $\partial_{1}A_{\nu}$ and the specific vorticity $\Omega_{\nu}$.
     Section 9 serves to establish the $\dot{H}^{m}$ energy bounds for the self-similar variables, thereby concluding the proof for the main theorem 4.1. Instead of studying the system for $(W,Z,A,K)$, we utilize various Sobolev inequalities to derive the energy bounds for the equations of the velocity $U$, the pressure $P$, and the entropy $H$ (as seen in (3.42)). This approach proves to be more convenient to apply the Friedrich’s energy estimates.

Through the whole paper, the following notations will be used unless stated otherwise.

• •

  Latin indices $\{i,j,k,l,\cdots\}$ take the values $1,2,3,$ and Greek indices $\{\alpha,\beta,\gamma,\cdots\}$ take the values $2,3$. Repeated indices are meant to be summed.

• •

  For a three-component vector $v$, denote the last two components of $v$ simply as $\check{v}$. For example, one can rewrite the gradient operator as $\nabla=(\partial_{1},\check{\nabla}).$

• •

  The convention $f\lesssim h$ means that there exists a universal positive constant $C$ such that $f\leq Ch$. The convention $A\sim B$ means that there exists a universal positive constant such that $A=O(B)$.

• •

  For any function $A(y,s)$, denote $A^{0}$ to be $A^{0}=A^{0}(s)=A(y,s)|_{y=0}$.

Consider the Surface $\Gamma:(x,t,u(x,t))$ in $R^{3}$. Then, the normal vector $N$ of $\Gamma$ and the Gauss curvature $K$ of $\Gamma$ at each point are given as follows

$$N=J^{-1}(-u_{x},-u_{t},1),\hskip 14.22636ptK=\dfrac{-u_{x}^{4}}{(1+u_{t}^{2}+u_{x}^{2})^{2}},$$

(2.17)

where $J=\sqrt{1+u_{x}^{2}+u_{t}^{2}}$. Initially, $N_{0}=\dfrac{1}{\sqrt{2}}(1,0,1)$, $K_{0}=-\dfrac{1}{4}$. In the self-similar coordinate, consider the evolutions of $N,K$ at $y=0$,

	$\displaystyle N$	$\displaystyle=\dfrac{1}{\sqrt{1+e^{2s}U_{y}^{2}+e^{s}U^{2}U_{y}^{2}}}\left(-e^{s}U_{y},-e^{\frac{s}{2}}UU_{y},1\right)\left|{}_{y=0}\right.,$		(2.18)
	$\displaystyle K$	$\displaystyle=-\dfrac{1}{(1+e^{-2s}U_{y}^{-2}+e^{-s}U^{2})^{2}}\left|{}_{y=0}\right..$		(2.19)

Then, as $s\to\infty$ ($t\to 0$), since $U_{y}(0,s)\to-1$, $U(0,s)\to\bar{U}(0)=0$, it holds that

$$N\to(-1,0,0),\hskip 14.22636ptK\to-1.$$

(2.20)

Therefore, shock formation to the Burgers’ equation (2.1) is equivalent to that

• •

  the transformation between the physical variables (the Cartesian coordinates) and the self-similar coordinates becomes degenerate;

• •

  the normal $N(t)$ of the shock front become horizontal at shock point.

We consider the following Cauchy problem for 1D Burgers’ equation with damping:

$$\left\{\begin{split}&\partial_{t}u+u\partial_{x}u=-au,\\ &u(x,t=-1)=f(x),\end{split}\right.$$

(2.21)

where $a\neq 0$ is the damping constant and $f(x)$ is the same as in (2.1). It will be shown that when $a<1$, the damping effect is small and (2.21) behaves like the standard Burgers’ equation such that smooth data can form a shock in finite time. Moreover, the damping will shift the blow-up time according to the sign of $a$ as follows.

• •

  If $0<a<1=-\min\frac{1}{\partial_{x}f}$, the damping will delay the formation of shock and the closer $a$ is to $1$, the larger the blow-up time is;

• •

  if $a<0$, then the blow-up time will be in advanced.

Precisely, the blow-up time can be clarified as $T_{\ast}=-\frac{1}{a}\ln(1-a)-1$⁷⁷7If one assumes generally that $\partial_{x}f(0)=-c$, then the first case becomes that $0<a<\frac{1}{c}$ while $T_{\ast}=-\frac{1}{a}\ln(1-\frac{a}{c})-1$. and in both case, the blow-up location remains the same compared with the undamped case (i.e. $x_{\ast}=0$). When $a\geq 1$, the damping effect is strong enough and (2.21) admits a global solution on $[-1,\infty)$. To this end, define the following self-similar coordinates⁸⁸8The choice of modulation variables $\tau(t)$ and $\xi(t)$ is mainly dependent on the invariance of (2.21) under space translation and time translation. That is, given any $x_{0},t_{0}\in\mathbb{R}$, then for any solution $u$ to (2.21), $\tilde{u}:=u(x-x_{0},t-t_{0})$ solves (2.21) as well.:

$$\begin{split}s&=-\ln(\tau(t)-t),\\ y&=e^{\frac{3}{2}s}(x-\xi(t)),\\ u(x,t)&=e^{-\frac{s}{2}}U(y,s),\end{split}$$

(2.22)

where $\tau(t)$ and $\xi(t)$ represent the blow-up time and location respectively, with the initial condition $\tau(-1)=\xi(-1)=0$. Note that the Jacobian of the transformation is given by

$$\dfrac{\partial(y,s)}{\partial(x,t)}=\left|\begin{array}[]{cc}e^{\frac{3}{2}s}&\ast\\ 0&(1-\dot{\tau})e^{s}\end{array}\right|=(1-\dot{\tau}(t))\dfrac{1}{(\tau(t)-t)^{\frac{5}{2}}}.$$

(2.23)

Hence, the transformation will degenerate at time $T_{\ast}=\tau(T_{\ast})$ (it will be shown that $1-\dot{\tau}(t)>0$) and is a diffeomorphism for $t<\tau(t)$. Then, (2.21) will be transformed into

$$\left(\dfrac{\partial}{\partial s}-\dfrac{1}{2}\right)U+\left(\dfrac{3}{2}y+\dfrac{U}{1-\dot{\tau}}-e^{\frac{s}{2}}\dfrac{\dot{\xi}(t)}{1-\dot{\tau}}\right)\partial_{y}U=-\dfrac{ae^{-s}U}{1-\dot{\tau}},$$

(2.24)

with the initial condition

$$U(y,0)=f(x),U(0,0)=0,\partial_{y}U(0,0)=-1.$$

(2.25)

To obtain a global solution for (2.24) is similar to the former case⁹⁹9In this case, define the characteristics as $$\dfrac{dY}{ds}=\frac{3}{2}Y+\frac{1}{1-\dot{\tau}}U(Y,s),\quad Y(0)=y_{0},$$ (2.26) along which $U\circ Y=e^{\frac{s}{2}}e^{-a(\tau-e^{-s}+1)}U_{0}(y_{0}):=g(s)U_{0}(y_{0})$. Then, $U$ can be solved implicitly as $U(y,s)=g(s)U_{0}(e^{-\frac{3}{2}s}y-\int_{0}^{s}\frac{e^{\frac{3}{2}s^{\prime}}}{1-\dot{\tau}}g(s^{\prime})ds^{\prime}\frac{1}{g(s)}U(y,s))$. To show the global existence of $U$, it suffices to show $$1+\frac{1}{a}(1-e^{-a(t+1)})U_{0}^{\prime}\neq 0,\quad\forall s\in[0,\infty),$$ (2.27) which is guaranteed by the initial data and the definition of $T_{\ast}$ (see(2.30)).. Here, we focus on the evolutions of $\tau(t)$ and $\xi(t)$. Differentiating (2.24) w.r.t $y$ yields

$$\left(\dfrac{\partial}{\partial s}-\dfrac{1}{2}\right)\partial_{y}U+\left(\dfrac{3}{2}+\dfrac{1}{1-\dot{\tau}}\partial_{y}U\right)\partial_{y}U+\left(\dfrac{3}{2}y+\dfrac{U}{1-\dot{\tau}}-e^{\frac{s}{2}}\dfrac{\dot{\xi(t)}}{1-\dot{\tau}}\right)\partial_{y}^{2}U=-\dfrac{ae^{-s}\partial_{y}U}{1-\dot{\tau}}.$$

(2.28)

Now we postulate that $U(0,s)=0$ and $\partial_{y}U(0,s)=-1$ for all $s\in[0,\infty)$. This can be achieved by choosing $\tau(t)$ and $\xi(t)$ suitably since initially, it holds that $U(0,0)=0$ and $\partial_{y}U(0,0)=-1$. Evaluating (2.24) and (2.28) at $y=0$ yields

$$\left\{\begin{split}&\dot{\xi}(t)=0,\\ &\dot{\tau}(t)=ae^{-s}=a(\tau(t)-t),\end{split}\right.$$

(2.29)

which implies

$$\left\{\begin{split}&\xi(t)=\xi(-1)=0,\\ &\tau(t)=t+e^{a(1+t)}+\dfrac{1}{a}(1-e^{a(1+t)}).\end{split}\right.$$

(2.30)

To conclude, it holds that

• •

  if $a\geq 1$, then $\tau(t)\geq t$ for all $t\geq-1$, and one obtains a global classical solution to (2.21);

• •

  if $a<1$¹⁰¹⁰10If one assumes generally that $\partial_{x}f(0)=-c$, then this case becomes ”if $a<\frac{1}{c}$, then a shock forms at $T_{\ast}=-\dfrac{1}{a}\ln(1-\frac{a}{c})-1$.” Therefore, for fixed $a$, small data(i.e.$c\leq\frac{1}{a}$) leads to the global solution while large data leads to the shock formation in finite time., then the transformation between Cartesian coordinates and self-similar coordinates will degenerate at time $T_{\ast}=\tau(T_{\ast})$ which can be computed as $T_{\ast}=-\dfrac{1}{a}\ln(1-a)-1$, and

  $$\lim_{t\to T_{\ast}}\partial_{x}u(\xi(t),t)=\lim_{t\to T_{\ast}}e^{s}\partial_{y}U(0,s)=\lim_{t\to T_{\ast}}\dfrac{-1}{\tau(t)-t}\to-\infty,$$

  (2.31)

  i.e. the solution will form a shock at $(x,t)=(0,T_{\ast})$. In this case, one sees that the blow-up location doesn’t shift.

First, do time-rescaling as

$$\mathrm{t}\rightarrow t=\dfrac{1+\alpha}{2}\mathrm{t},$$

so that all the modulation variables are defined w.r.t the rescaled time (i.e. $t$).
     Given an unit vector $n(t)=(\sqrt{1-n_{2}^{2}(t)-n_{3}^{2}(t)},n_{2}(t),n_{3}(t))=(n_{1},\breve{n}(t))$, one can generate a rotation matrix which transforms $e_{1}$ to $n_{1}(t)$. Precisely, let

$$R(t)=\left(\begin{array}[]{ccc}n_{1}&-n_{2}&-n_{3}\\ n_{2}&1-\frac{n_{2}^{2}}{1+n_{1}}&-\frac{n_{2}n_{3}}{1+n_{1}}\\ n_{3}&-\frac{n_{2}n_{3}}{1+n_{1}}&1-\frac{n_{3}^{2}}{1+n_{1}}\end{array}\right).$$

(3.2)

Then, the rotation matrix $R(t)$ rotates $e_{1}$ to $n(t)$. Define the new coordinates as

$$(\tilde{x},t)=(R^{T}(\mathit{x}-\xi(t)),t)$$

(3.3)

and the corresponding fluid variables as

$$\tilde{\sigma}(\tilde{x},t)=\sigma(\mathit{x},t),\quad\tilde{u}(\tilde{x},t)=R^{T}(t)u(\mathit{x},t),\quad\tilde{k}(\tilde{x},t)=k(x,t).$$

(3.4)

Note that

	$\displaystyle\dfrac{\partial}{\partial t}$	$\displaystyle=\dfrac{1+\alpha}{2}\dfrac{\partial}{\partial\mathit{t}}+\tilde{v}\dfrac{\partial}{\partial\tilde{x}_{i}},$		(3.5)
	$\displaystyle\dfrac{\partial}{\partial\mathit{x}_{i}}$	$\displaystyle=R_{ik}\dfrac{\partial}{\partial\tilde{x}_{k}},$		(3.6)

where

$$\tilde{v}=\dot{Q}\tilde{x}-R^{T}\dot{\xi},\quad\dot{Q}=\dot{R}^{T}R,\quad\dot{Q}_{ij}=-\dot{Q}_{ji}.$$

(3.7)

Hence, the Euler system (1.2) is transformed in $(\tilde{x},t)$ as

$$\begin{split}\dfrac{1+\alpha}{2}\partial_{t}\tilde{\sigma}+(\tilde{u}+\tilde{v})\cdot\nabla_{\tilde{x}}\tilde{\sigma}+\alpha\sigma div_{\tilde{x}}\tilde{u}&=0,\\ \dfrac{1+\alpha}{2}\partial_{t}\tilde{u}-\dot{Q}\tilde{u}+(\tilde{u}+\tilde{v})\cdot\nabla_{\tilde{x}}\tilde{u}+\alpha\tilde{\sigma}\nabla_{\tilde{x}}\tilde{\sigma}&=\dfrac{\alpha}{2\gamma}\tilde{\sigma}^{2}\nabla_{\tilde{x}}\tilde{k}-a\tilde{u},\\ \dfrac{1+\alpha}{2}\partial_{t}\tilde{k}+(\tilde{u}+\tilde{v})\cdot\nabla_{\tilde{x}}\tilde{k}&=0,\end{split}$$

(3.8)

together with the vorticity equation

$$\dfrac{1+\alpha}{2}\partial_{t}\tilde{\zeta}-\dot{Q}\tilde{\zeta}+(\tilde{u}+\tilde{v})\cdot\nabla_{\tilde{x}}\tilde{\zeta}-\tilde{\zeta}\cdot\nabla_{\tilde{x}}\tilde{u}=\dfrac{\alpha}{\gamma}\dfrac{\tilde{\sigma}}{\tilde{\rho}}\nabla_{\tilde{x}}\tilde{\sigma}\times\nabla_{\tilde{x}}\tilde{k}-a\tilde{\zeta}.$$

(3.9)