3-D axisymmetric transonic shock solutions of the full Euler system in divergent nozzles

We establish the stability of 3-D axisymmetric transonic shock solutions of the steady full Euler system in divergent nozzles under small perturbations of an incoming radial supersonic flow and a constant pressure at the exit of the nozzles. To study 3-D axisymmetric transonic shock solutions of the full Euler system, we use a stream function formulation of the full Euler system for a 3-D axisymmetric flow. We resolve the singularity issue arising in stream function formulations of the full Euler system for a 3-D axisymmetric flow. We develop a new scheme to determine a shock location of a transonic shock solution of the steady full Euler system based on the stream function formulation.


Introduction
In [14, Chapter 147], authors, using an approximate model, describe a transonic shock phenomenon for a compressible invicid flow of an ideal polytropic gas in a convergent-divergent type nozzle called de Laval nozzle: If a subsonic flow accelerating as it passes through the convergent part of the nozzle reaches the sonic speed at the throat of the nozzle, then it becomes a supersonic flow right after the throat of the nozzle. It further accelerates as it passes through the divergent part of the nozzle. If an appropriately large exit pressure p e is imposed at the exit of the nozzle, then at a certain place of the divergent part of the nozzle, a shock front intervenes, the flow is compressed and slowed down to subsonic speed. The position and strength of the shock front are automatically adjusted so that the end pressure at the exit becomes p e . This phenomenon was rigorously studied using radial solutions of the full Euler system in [31] (it was shown that in a divergent nozzle, for given a constant supersonic data on the entrance of the nozzle and an appropriately large constant pressure on the exit of the nozzle, there exists a unique radial transonic shock solution satisfying these conditions). Motivated by this phenomenon, there were many studies on the stability of transonic shock solutions in divergent nozzles (structural stability of radial transonic shock solutions in divergent nozzles under multi-dimensional perturbations of an entrance supersonic data and exit pressure) and related problems.
The stability of one-dimensional transonic shock solutions in flat nozzles was first studied. This subject was studied using the potential flow model in [7,8,9,26,27] and further studied using the full Euler system in [11,29,6,5,25,30,12,13]. These results showed that onedimensional transonic shock solutions in flat nozzles are not stable under a perturbation of a physical boundary condition (supersonic data on the entrance or density, pressure or normal velocity on the exit) and, even if one-dimensional transonic shock solutions in flat nozzles are stable, their shock locations are not uniquely determined unless there exists the assumption that a shock location passes through some point on the wall of the nozzle, as it can be expected from the behavior of one-dimensional transonic shock solutions in flat nozzles (as a shock location changes, the value of the subsonic part of an one-dimensional transonic shock solution in a flat nozzle does not change). After that, the stability of radial transonic shock solutions in divergent nozzles was studied. This subject was first studied using the full Euler system in [19,28]. In these results, authors, by considering a perturbation of radial transonic shock solutions in divergent nozzles, could show that a shock location is uniquely determined for given an exit pressure without the assumption that a shock location passes through some point of the nozzle but they only had the result under the assumption that the tip angle of the nozzle is sufficiently small. After that, this subject without restriction on the tip angle of the nozzle was studied. In [3], the authors studied this subject using the non-isentropic potential model introduced in [3]. And they obtained the stability result for radial transonic shock solutions in divergent nozzles. This subject was also studied using the full Euler system. This for the 2-D case was done in [24,18,21]. In these papers, the authors had the stability result for radial transonic shock solutions in divergent nozzles. Especially, the authors in [21] had the result for flows having C 1,α interior and C α up to boundary regularity, so that they could consider a general perturbation of a nozzle. This for the 3-D case for axisymmetric flows with zero angular momentum components was done in [20]. The authors in this paper also had the same result. This for the general 3-D case was done in [10,23]. The authors in [10,23] also had the same result but under S-condition introduced in [10,23]. Recently, this subject for the general 3-D case for flows having some friction term was studied in [32].
In this paper, we study the stability of radial transonic shock solutions in divergent nozzles under small perturbations of an incomming radial supersonic solution and a constant exit pressure using the full Euler system for the 3-D case for axisymmetric flows. We consider axisymmetric flows with non-zero angular momentum components. (This is a difference from [20].) We consider a divergent nozzle having no restriction on the tip angle of the nozzle and do not have any assumption on an incomming supersonic solution.
The main new feature in this paper is to develop a new iteration scheme to determine a shock location for a transonic shock solution of the steady full Euler system in a divergent nozzle and resolve the singularity issue arising in the stream function formulations of the full Euler system using an elliptic system approach.
To deal with the stability of 3-D axisymmetric transonic shock solutions of the full Euler system, we use a stream function formulation for the full Euler system for an axisymmetric flow. This formulation shows the fact that an initial shock position and a shape of a shock location (see the definitions below the proof of Theorem 2.17) are determined in different mechanisms clearly. Based on this formulation and using the fact that the entropy of the downstream subsonic solution of a radial transonic shock solution on a shock location monotonically increases as a shock location moves toward the exit of the nozzle (see Lemma 2.10), we develop a new scheme to determine a shock location of a transonic shock solution of the full Euler system in a divergent nozzle: 1. Pseudo Free Boundary Problem 2. Determination of a shape of a shock location (see below the proof of Theorem 2.17).
In technical part, we resolve the singularity issue arising in stream function formulations of the full Euler system. A stream function formulation when it is formulated by using the Stokes' stream function (see (2.4.1)) has a singularity issue at the axis of symmetry. We resolve this singularity issue by formulating a stream function formulation using the vector potential form of the stream function (see §2.4) and solving a singular elliptic equation appearing in this stream function formulation as an elliptic system (see §3.2). The stream function formulation formulated by using the vector potential form of the stream function contains a singular elliptic equation. We transform this singular elliptic equation into a form of an elliptic system and, by solving the elliptic system form as an elliptic system, solve the singular elliptic equation. (We also use this approach to prove the orthogonal completeness of eigenfunction of an associated Legendre problem of type m = 1 with a general domain (see Lemma 4.3)). Using the stream function formulation formulated by using the vector potential form of the stream function, we obtain the stability result for flows having C 1,α interior and C α up to boundary regularity. This paper is organized as follows. In Section 2, we present definitions and a basic lemma used throughout this paper and introduce our problem and result. In this section, we introduce the stream function formulation used in this paper. In Section 3, we solve the Pseudo Free Boundary Problem. In this section, we study a linear boundary value problem for a singular elliptic equation and an initial value problem of a transport equation appearing in the Pseudo Free Boundary Problem, and prove the unique existence of solutions of the Pseudo Free Boundary Problem. In Section 4, we show the existence and uniqueness of transonic shock solutions. In Section 5, we present some computations done by using the tensor notation given in §3.2.

Problem and Theorem
2.1. Preliminary. In this paper, we consider a 3-D steady compressible invicid flow of an ideal polytropic gas. The motion of this flow is governed by the following full Euler system:      div(ρu) = 0, div(ρu ⊗ u + pI) = 0, div(ρuB) = 0 for an ideal polytropic gas: if M > 1, then a flow is called supersonic, if M = 1, then it is called sonic and if M < 1, then it is called subsonic. It is generally known that types of the system varies depending on the value of M . If M > 1, then the system is a hyperbolic system and if M < 1, then the system is an elliptic-hyperbolic coupled system. When the flow passes through a domain having a certain geometric structure or satisfies a certain boundary condition, it may have a discontinuity across a surface in the domain in the direction of the flow. Such a discontinuity is called a shock.
A shock solution of (2.1.1) is defined as follows.
In this paper, we deal with a 3-D axisymmetric transonic shock solution of (2.1.1). For precise statement, we define an axisymmetric domain and axisymmetric functions used in this paper. For later use, we present a lemma used to deal with regularities of axisymmetric functions.
Using this spherical coordinate system, an axisymmetric domain and axisymmetric functions are defined as follows.
Definition 2.7. In this paper, when a velocity field u is represented as u = u r e r + u θ e θ + u ϕ e ϕ , u ϕ e ϕ is called the angular momentum component of u.
For later use, we present the following lemma that shows when an axisymmetric function in C k as a function of the spherical coordinate system is in C k as a function of the cartesian coordinate system. This lemma is obtained from [22,Corollary 1]. (ii) f e θ and f e ϕ are in C k (Ω) for k ∈ 0, 1, 2, . . . if and only if f is in C k as a function of spherical coordinate system in Ω and ∂ 2m θ f = 0 for all 0 ≤ m ≤ ⌊ k 2 ⌋. In this paper, we use the same function notation when we represent an axisymmetric function as a function on the cartesian coordinate system or spherical coordinate system.
To introduce our problem and for our later analysis, we study a radial transonic shock solution of (2.1.1) in N .
Fix positive constants (ρ in , u in , p in ) satisfying M in (:= u in / γp in ρ in ) > 1. Let (ρ,ūe r ,p) be a radial shock solution of (2.1.1) in N with a shock Γ t := {r = t} ∩ N for some t ∈ [r 0 , r 1 ] satisfying (ρ,ūe r ,p) = (ρ in , u in e r , p in ) on Γ en := ∂N ∩ {r = r 0 , 0 ≤ θ < θ 1 }. with (ρ,ū,p)(r 0 ) = (ρ in , u in , p in ) (2.2.4) in D − t , where D − t := {r 0 < r < t}, ′ is the derivative with respect to r andB :=ū 2 2 + γp (γ−1)ρ , and is a solution of (2. 3) with (2.2.5) has a unique solution (ρ,ū,p). By these two facts and the local unique existence theorem for ODE, (2.2.3) with (2.2.5) has a solution (ρ,ū,p) satisfyingM < 1 andM ′ < 0 in D + t . Therefore, a radial shock solution (ρ,ūe r ,p) uniquely exists in N for each t ∈ [r 0 , r 1 ] and it is a radial transonic shock solution. From this fact andM | D + t < 1 in D + t , we obtainρ| ′ One can see that the values of (ρ,ū,p)| D + t at a fixed location r in D + t are determined by the three conserved quantities in the right-hand sides of the equations in (2.2.10). This combined with the fact that the conserved quantity forS in D + t given in (2.2.10) varies depending on t (obtained from (2.2.10) by using (M | t for any t ∈ [r 0 , r 1 ]) implies that the values of (ρ,ū,p)| D + t at a fixed location r in D + t vary depending on t. To represent this dependence, we write (ρ,ū,p)| D + t (r) andS| D + t (r) as (ρ,ū,p)| D + t (r; t) andS| D + t (r; t), respectively.
The conserved quantity forS in D + t satisfies the following monotonicity. Lemma 2.10. Let r 0 , r 1 , t be positive constants such that r 0 ≤ t ≤ r 1 . Suppose that (ρ,ū,p) is as above. Then there holds dS| > 0 for any t ∈ [r 0 , r 1 ].
Proof. DifferentiateS| D + t (t; t) with respect to t. Then we have One can easily check that g(1) = 1 and g ′ (x) > 0 for all x > 1. By this fact, From Lemma 2.10, we obtain the following result.
Proof. By the definitions ofB andS and the first and third equation of (2.2.10), Differentiate this with respect to t. Then we get t , Lemma 2.10 and the second equation of (2.2.10), we obtain from (2.2.12) the desired result.
The above proposition implies that for any given p c ∈ [p 1 , p 2 ] where p 1 := p 0 | D + r 1 (r 1 ; r 1 ) and p 2 := p 0 | D + r 0 (r 1 ; r 0 ), there is a unique shock location Γ t in N such that (ρ, ue r , p) satisfies p| D + t (r 1 ; t) = p c . Hereafter, we fix a constant p c ∈ (p 1 , p 2 ) and denote t ∈ (r 0 , r 1 ) such that a radial transonic shock solution of (2.1.1) satisfying (2.2.2) and having a shock location Γ t satisfies p(r 1 ) = p c by r s . Also, we denote a solution (ρ, u, p) of (2.

2.3.
Problem. Using the radial transonic shock solution given in the previous subsection, we present our problem.
In this paper, we use the following weighted Hölder norm. For a bounded connected open set Ω ⊂ R n , let Γ be a closed portion of ∂Ω. For x, y ∈ Ω, set δ x := dist(x, Γ) and δ x,y := min(δ x , δ y ).
Problem 1 (Transonic shock problem). Given an axisymmetric supersonic solution (ρ − , u − , p − ) of (2.1.1) in N satisfying the slip boundary condition where n w is the unit normal vector on Γ w , and an axisymmetric exit pressure p ex on Γ ex : for a positive constant σ, find a shock location Γ f := N ∩ {r = f (θ)} and a corresponding subsonic solution (iv) and the exit pressure condition Remark 2.12. It is generally known that a supersonic solution of (2.1.1) is governed by a hyperbolic system. We assume that (ρ − , u − , p − ) in Problem 1 exists. Remark 2.13. To simplify our argument, we assumed in Problem 1 that (ρ − , u − , p − ) satisfies (2.3.2). This assumption will be used to reduce (2.1.1) and (2.1.3) (see §2.5). The result for Problem 1 (Theorem 2.17) does not change if we consider a general perturbation of (ρ − 0 , u − 0 e r , p − 0 ) in Problem 1. We study Problem 1 using a stream function formulation of the full Euler system for an axisymmetric flow. We introduce a stream function formulation used in this paper in the next subsection.

Stream function formulation.
Let Ω be an open simply connected axisymmetric set in R 3 . Let (ρ, u) be axisymmetric C 1 functions in Ω satisfying the first equation of (2.1.1) and |ρu| > 0. For such (ρ, u), the Stokes' stream function for an axisymmetric flow of the full Euler system is defined by where S x is a simply connected C 1 surface in Ω whose boundary is a circle centered at z-axis, parallel to xy-plane and passing through x, and ν is the unit normal vector on S x pointing outward direction with respect to the cone-like domain made by connecting ∂S x and the origin by straight lines. By the first equation of (2.1.1), the value of this function at x is independent of the choice of S x . Since ∂S x is axisymmetric, V is axisymmetric in Ω.
By the first equation of (2.1.1), V is a constant on each stream surface of ρu in Ω. Here, the stream surfaces of a vector field ρu in Ω by a set of surfaces made by collecting all the streamlines of ρu initiating from a point on a circle in Ω centered at z-axis and parallel to xy-plane. By |ρu| > 0 in Ω, V is a constant on each stream surface of ρu in Ω and V on each different stream surface of vector field ρu in Ω is different from each other. From these facts, we have that if we apply ∇ ⊥ , where ∇ ⊥ = 1 2πr sin θ e r ∂ θ r − e θ ∂ r which satisfies ∇ ⊥ h · ∇h = 0 and ∇ ⊥ h · e ϕ = 0 for a scalar function h, to V , then we have a vector field in Ω tangent to the stream surfaces in Ω and having no e ϕ component. Apply ∇ ⊥ to V . Then we have where u r = u · e r and u θ = u · e θ . Using (2.4.2), we can reformulate the full Euler system for an axisymmetric flow. But if we do this, then there is a singularity issue that can be seen in the relation ||ρu r e r + ρu θ e θ || α,N = ||∇ ⊥ V || α,N ≤ C||V || 1,α,N for any constant C. To avoid this issue, we use the following form of the stream function.
Let Φe ϕ be an axisymmetric vector field in Ω satisfying where r is a parametrization of ∂S x in a counter clockwise direction. Then by the definitions of Φe ϕ and V , It is easily checked that We call Φe ϕ the vector potential form of the stream function.
We reformulate the full Euler system for an axisymmetric flow using (2.4.6). For our later analysis, when we reformulate the full Euler system using (2.4.6), we use the following form of the full Euler system representing the relation between ∇S and ∇ × u clearly which is obtained under the assumption that (ρ, u, p) ∈ C 1 and ρ > 0.
Assume that Γ in (2.1.3) is an axisymmetric C 1 surface. Let τ 2 and τ 1 in (2.1.3) be e ϕ and the unit tangent vector field on Γ perpendicular to e ϕ and satisfying ν · (τ 1 × e ϕ ) > 0, respectively, where ν is the unit normal vector field on Γ pointing toward Ω + . By the definition of V given in (2.4 By the definition of L, the third equation of (2.1.3) can be written as From the second, third and fifth equation of (2. is a function defined in (2.2.9) and variables with lower indices ± denote variables in Ω ± , respectively. Combining these reformulated equations of (2.1.3) and the fifth equation of (2.1.3), we have the following stream function formulation of (2.1.3): Then we present the stream function formulations of (2.1.1) and (2.1.3) satisfied by (ρ + , u + , p + ) in Problem 1. By (2.4.15) and (2.4.20), the stream function formulations of (2.1.1) and (2.1.3) satisfied by (ρ + , u + , p + ) in Problem 1 and reduced by using the assumption that (ρ + , u + , p + ) where ν f is the unit normal vector on Γ f pointing toward N + f and τ f is the unit tangential vector on Γ f perpendicular to e ϕ and satisfying ν f · (τ f × e ϕ ) > 0. Here, (Φ − e ϕ , L − , S − ) and (Φ + e ϕ , L + , S + ) are (Φe ϕ , L, S) given by the definitions of Φe ϕ , L and S for (ρ, u, p) = (ρ − , u − , p − ) and (ρ + , u + , p + ), respectively.
Using the equations and boundary conditions obtained above, Problem 1 is restated as follows: Let S 2,θ 1 := {(x, y, z) ∈ R 3 | r = 1, 0 ≤ θ < θ 1 }. A function f representing a shock location Γ f can be considered as a function on S 2,θ 1 . Using this fact and the stereographic projection from (0, 0, −1) onto the plane z = 1 passing through S 2,θ 1 , we see that f can be regarded as a function on Λ where Λ := {(x, y) ∈ R 2 | x 2 + y 2 < 2 tan θ 1 2 }. Thus, f can be regarded as a function on Λ or (0, θ 1 ). In this paper, we regard f in both ways. To simplify our notation, we use the same function notation when we represent f as a function on Λ or (0, θ 1 ).
Our result of Problem 2, the main result in this paper, is given as follows.
Hereafter, we say that a constant depends on the data if a constant depends on (ρ in , u in , p in ), p c , γ, r 0 , r 1 , θ 1 and α.
The following result of Problem 1 is obtained from Theorem 2.16.
To describe our process of proving Theorem 2.16, we define some terminologies. Let f : . Then by f s (0) = 0, f s is uniquely determined by f ′ s . We call f (0) and f ′ s the initial shock position and the shape of a shock location, respectively. Using these terminologies, our process of proving Theorem 2.16 is described as follows.
1. For given an incomming supersonic solution, an exit pressure and a shape of a shock location (ρ − , u − , p − , p ex , f ′ s ) in a small perturbation of (ρ − 0 , u − 0 e r , p − 0 , p c , 0), show that there exists a pair of an initial shock position f (0) and a subsonic solution (Φe ϕ , L, S) of (2.5.1)-(2.5.3) satisfying all the conditions in Problem 2 except (2.5.5), and that this solution is unique in the class of functions in a small perturbation of (r s , Φ + 0 e ϕ , L + 0 , S + 0 ). 2. For given an incomming supersonic solution and an exit pressure as in Step 1 or in a much small perturbation of (ρ − 0 , u − 0 e r , p − 0 , p c ) if necessary, show that there exists f ′ s in a small perturbation of 0 as in Step 1 such that (f (0), Φe ϕ , L, S) determined by (ρ − , u − , p − , p ex , f ′ s ) in Step 1 satisfies (2.5.5), and that for given (ρ − , u − , p − , p ex ) in a small perturbation of (ρ − 0 , u − 0 e r , p − 0 , p c ), a solution (f, Φe ϕ , L, S) of Problem 2 is unique in the class of functions in a small perturbation of (r s , Φ + 0 e ϕ , L + 0 , S + 0 ). Once Step 1 and Step 2 are done, then f = f (0) + f s and (Φe ϕ , L, S) obtained through Step 1 and Step 2 satisfies all the conditions in Problem 2. Thus, Theorem 2.16 is proved if Step 1 and Step 2 are done. We will deal with Step 1 and Step 2 in Section 3 and Section 4, respectively.
Note that the fact that for transonic shock solutions of the full Euler system, an initial shock position and a shape of a shock location are determined in different mechanisms (an initial shock position is determined by the solvability condition related to the mass conservation law and a shape of a shock location is determined by the R-H conditions) was pointed out in [24] and the authors in [24,24,18,21,20,10,23] prove the stability of transonic shock solutions of the full Euler system using iteration schemes based on this fact. In this paper, we also prove the stability of transonic shock solutions of the full Euler system using a scheme based on this fact. But we do this using a different scheme. In our scheme, a non-local elliptic equation appearing in [24,18,21,20,10,23] does not appear.

Pseudo Free Boundary Problem
As a first step to prove Theorem 2.16, we will solve the Pseudo Free Boundary Problem below. This problem naturally arises from the requirement that a subsonic solution in Problem 2 must satisfy (2.5.12). From the linearized equation of (2.5.11), it is seen that an iteration scheme for a fixed boundary problem does not give a subsonic solution satisfying (2.5.12) in general. Thus, an iteration scheme for a fixed boundary problem is not a proper scheme to find a subsonic solution in Problem 2. To find a subsonic solution satisfying (2.5.12), we let f (0) be an unknown to be determined simultaneously with a subsonic solution. Using this variable, we adjust the value of a subsonic solution so that this solution can satisfy (2.5.12). The main ingredient for this argument to hold is the monotonicity of the entropy of the downstream subsonic solution of a radial transonic shock solution on a shock location with respect to the shock location. This will be seen in the proof of Proposition 3.1.
There exists a positive constant σ 3 depending on the data such that if σ ∈ (0, σ 3 ], then Problem 3 has a solution (f (0), Φe ϕ , L, S) satisfying where C is a positive constant depending on the data. Furthermore, this solution is unique in the class of functions satisfying (3.0.2).

Linearization and reformulation of (A) and (B).
We linearize (A) with respect to (A) satisfied by (Φ + 0 e ϕ , L + 0 , S + 0 ). Since ρ in the first and fourth equation of (A) is given using an implicit relation, to obtain the linearized equations of (A), we first linearize ρ with respect to ρ + 0 .
Proof. With (3.1.5), (3.1.6), (3.1.15) and the fact that if an axisymmetric vector field a on an axisymmetric connected open set Ω is in C k (Ω), then ∇ × a ∈ C k−1 (Ω) and a · e r ∈ C k (Ω) (the second one is obtained from Lemma 2.8), we estimate (3.1.8) and (3.
, respectively. Then we obtain the desired result.
By the fact that (ρ − , u − , p − ) and f = f (0) + f s are axisymmetric, T en,f defined in (3.1.19) can be regarded as a function of θ. As a function of θ, T en,f can be written as To estimate ||T en,f || ). Then we obtain , ν f − e r can be written as Substitute this expression of ν f − e r into ν f − e r in (3.1.23) and then estimate (3.1.23) in .
From (A ′ ), (B ′ ), the Pseudo Free Boundary Problem is naturally derived. We explain this below.
For a given (ρ − , u − , p − , p ex ), find (Ψe ϕ , A, T ) satisfying (A ′ ), (B ′ ) using an iteration scheme for a fixed boundary problem (for example, in a fixed domain N + f (0)+fs , for a given Ψe ϕ , solve (B ′ ), substitute the resultant A and T and the previously given Ψe ϕ into the right-hand sides of (A ′ ), obtain a new Ψe ϕ by solving this (A ′ ) and show that a new Ψe ϕ is equal to the given Ψe ϕ using a fixed point argument). Then since (Ψe ϕ , A, T ) we find in this way does not satisfy the fifth equation of (A ′ ) in general, this iteration scheme does not give a subsonic solution of (2.5.1)-(2.5.3) satisfying (2.5.12) in general. From the facts that the entropy at a point on a shock location in the subsonic side is conserved along the streamline passing through that point and the entropy of the downstream subsonic solution of a radial transonic shock solution in a divergent nozzle on a shock location monotonically increases as a shock location moves toward the exit (see Lemma 2.10), we see that we can find (Ψe ϕ , A, T ) satisfying the fifth equation of (A ′ ) by varying S on Γ ex by adjusting f (0). From this fact, Problem 3 is derived.
(A ′ ) and (B ′ ) are of the form of one linear boundary value problem for a singular elliptic equation (this will be seen in the next subsection) and two initial value problems of a transport equation whose coefficient is an axisymmetric and divergence-free vector field, respectively. We will study these problems, seperately, in §3.2 and §3.3.

Linear boundary value problem for a singular elliptic equation. Fix the righthand sides of the first and fourth equation in (A
where f , F and h i e ϕ for i = 1, 2 are functions given in Lemma 3.5. Since (3.2.1) is expressed as 2) as a boundary value problem for an elliptic system. The following is the main result in this subsection.
. Suppose that f is as in Lemma 3.2 and satisfy f ′ (θ 1 ) = 0. Also, suppose that F : Finally, suppose that h 1 e ϕ : Γ f → R 3 and h 2 e ϕ : Γ ex → R 3 are axisymmetric functions in Then the boundary value problem (3.2.1), (3.2.2) has a unique axisymmetric C 2,α where C is a positive constant depending only on (ρ + 0 , u + 0 , p + 0 ), γ, r s , r 1 , θ 1 and α, and . This form will be used in the proof of Lemma 3.10.
To avoid the singularity issue in (3.2.1), (3.2.2), we deal with (3.2.1), (3.2.2) as a boundary value problem for a vector equation. From ∇ × (∇ × (Ψe ϕ )) = −∆(Ψe ϕ ), we expected that (3.2.1) can be transformed into a form of an elliptic system. We, motivated by the work in [4], thought that if (3.2.1) can be transformed into a solvable elliptic system form, then the unique existence and regularity of solutions of (3.2.1), (3.2.2) can be obtained by obtaining those of solutions of the elliptic system form of (3.2.1), (3.2.2) as a boundary value problem for an elliptic system.
For this argument to hold, it is needed to find a solvable elliptic system form of (3.2.1). For computational convenience to find such a form and for our later argument (reflection argument in the proof of Lemma 3.9 and Lemma 3.10), we use the following tensor notation.
Tensor notation Let a ⊗ b = ab T for a, b ∈ R 3 . Then a ⊗ b is a linear map from R 3 to R 3 and any linear map from R 3 to R 3 can be represented using this operator. This notation can be extended so that using the extension of this operator, we can represent any linear map from R 3×3 to R 3×3 . For any a, b, c, d ∈ R 3 , let a ⊗ b ⊗ c ⊗ d be an operator satisfying where e, f ∈ R 3 . Then a ⊗ b ⊗ c ⊗ d is a linear map from R 3×3 to R 3×3 and any linear map from R 3×3 to R 3×3 can be represented using this operator.
By direct computation done by using the above tensor notation, we found the following form of (3.2.1) where I is the identity map from R 3×3 to R 3×3 and I ⊗ e r ⊗ e r ⊗ I is a linear map from R 3×3 to R 3×3 satisfing (I ⊗ e r ⊗ e r ⊗ I)(a ⊗ b) = (b · e r )a ⊗ e r for any a, b ∈ R 3 (see the definition of I ⊗ a ⊗ b ⊗ I for any a, b ∈ R 3 in (5.0.4)). By M + 0 < 1 in N + f and the boundedness of (ρ + 0 , u + 0 , p + 0 ) in N + f for N + f ⊂ N + rs−δ 1 , there exist positive constants µ and M such that Hence, (3.2.8) is a form of a solvable elliptic system for a dirichlet boundary condition. We obtain the unique existence and regularity of solutions of (3.2.1), (3.2.2) by obtaining those of solutions of (3.2.8), (3.2.2) as a boundary value problem for an elliptic system. The result of the unique existence and regularity of solutions of (3.2.8), (3.2.2) as a boundary value problem for an elliptic system is given in the following lemma.
One can see that that Lemma 3.5 is obtained from Lemma 3.7. To prove Lemma 3.5, in the remainder of this subsection, we prove Lemma 3.7.
We first prove the unique existence of weak solution of (3.2.13), (3.2.14).
Lemma 3.8. Under the assumptions as in Proposition 3.5, the boundary value problem (3.2.13), (3.2.14) has a unique weak solution for a constant C > 0 and µ||U ♯ || 2 With these facts, we apply the Lax-Milgram Theorem to (3.2.15). Then we obtain that there exists a unique This finishes the proof.
We next prove that this weak solution is in Lemma 3.9. Under the assumptions as in Proposition 3.5, let U ♯ be a weak solution of the boundary value problem (3.2.13), (3.2.14). Then for any β ∈ (0, 1), and N + f , and τ is the modulus of continuity of A in N + f given as Under the assumptions as in Proposition 3.5, let U ♯ be a weak solution of (3.2.13), (3.2.14). Then and F i for i = 1, 2, 3 are constants given in Lemma 3.10.
We will prove Lemma 3.9 and Lemma 3.10 using the method of freezing the coefficients (Korn's device of freezing the coefficients) (see [15,Chapter 3]). Since N + f is a Lipshitz domain, U ♯ ∈ C β (N + f ) and U ♯ ∈ C 1,α (N + f ) can be proved by showing that (i) there are positive constants C and R such that for all x 0 ∈ N + f , and (ii) there are positive constants C and R such that We prove (i) and (ii) by obtaining (3.2.17) and (3.2.18) at each point x 0 in N + f for C and R independent of x 0 using the method of freezing the coefficients. When we do this, there exists some difficulty. For the case of we can obtain the integral estimates for the fixed coefficients equation using the Cacciopolli inequality and the quotient difference method, and obtain (3.2.17) and (3.2.18) at x 0 ∈ N + f or Γ f ∪Γ + w ∪Γ ex using these estimates and the method of freezing the coefficients (see [1,Chapter 6]). But for the case of x 0 ∈ Γ f ∩ Γ + w or Γ + w ∩ Γ ex , we cannot obtain the integral estimates for the fixed coefficients equation using the Cacciopolli inequality and the quotient difference method. Thus, we cannot obtain (3.2.17) and (3.2.18) at x 0 ∈ Γ f ∩ Γ + w or Γ + w ∩ Γ ex using the standard method of freezing the coefficients. We resolve this difficulty by developing some reflection argument that holds for a linear boundary value problem on a Lipschitz domain whose all corners are perpendicular for an elliptic system whose the domain part of principal coefficients is diagonal with respect to the coordinate systems representing the walls near the corners of the domain. This will be seen in the proof of Lemma 3.9 and Lemma 3.10.
Hereafter, we use the following notation: T : a one dimensional torus with period 2π, , ϕ ∈ T} r, t : a radius of a ball in the spherical coordinate system, To prove Lemma 3.9 and Lemma 3.10, we prove the following lemma.
Then for any t such that 0 < t ≤ r, there hold where C is a positive constant depending on µ, M, f (θ 1 ) and θ 1 and The result is obtained by using the reflection argument.
). There exists a unique weak solution of (3.2.23), We denote the weak solution of (3.2.23), (3.2.24) by W .
The following Corollary is obtained from Lemma 3.11 in the same way that Corollary 3.11 is obtained from Lemma 3.10 in [17]. We omit the proof.
We first prove Lemma 3.9.
Proof of Lemma 3.9. We prove Lemma 3.9 by proving that (3.2.17) holds for all x 0 ∈ N + f for some positive constants C and R.
Next, we prove Lemma 3.10. We prove Lemma 3.10 using the method of freezing the coefficients and the reflection arguments in the proof of Lemma 3.9. When we do this, there exists some problem: since F ♯ is not in L p (N + f ) for q = 3 1−α nor has the form divG with G ∈ C α (N + f ), we cannot get the power of t required in (3.2.18) from the integral estimate of F ♯ directly. We obtain this power by delivering θ-derivatives imposed on some functions in F ♯ to the functions multiplied to those functions in the integral form of F ♯ using integration by parts and estimating the resultant integral form of F ♯ . To make our argument clear, we present the detailed proof.
Using integration by parts, we change this equation into Step 1 in the proof of Lemma 3.9, transform this equation. Then we obtain A,d andŨ ♯ are functions given below (3.2.31), x is the cartesian coordinate representing N + f , y = Π(x) and (r, θ, ϕ) and (r,θ,φ) are the spherical coordinate systems for x and y, respectively. As we did in Step 2 in the proof of Lemma 3.9, set ξ = 0 outside of Ξ −1 (D * r (x * 0 )) for 0 < r < min(θ 1 , π, r 1 −f (θ 1 )) where x * 0 = (f (θ 1 ), θ 1 , ϕ 0 ) for some ϕ 0 ∈ T and then transform this equation using Ξ. Then we obtain for all ξ ∈ H 1 0 (D * r (x * 0 )) whereŨ * ,d * , andM are functions defined below (3.2.33) and we used div (r,θ,φ) ( 1 detMM TÃ • Ξ −1M D (r,θ,φ)h * )(x * 0 ) = 0. Fix the principal coefficients of the left-hand side of the above equation. Then we get ). Let W be the weak solution of (3.2.35), (3.2.36). As we did in Step 3 in the proof of Lemma 3.9, subtracting the weak formulation of (3.2.35), (3.2.36) from the above equation. And then take ξ = V where V =Ũ * − W to the resultant equation. Then we have Using the Sobolev and Hölder inequality and the facts that we obtain from the above equation In the proof of Lemma 3.9, we showed that for any ε ∈ (0, 1), there exists R 4 > 0 such that Using this inequality and U ♯ ∈ C 0 (N + f ) obtained from Lemma 3.9, we apply Lemma 2.
When x * 0 is in Γ +. * w ∩ Γ * ex , we obtain (3.2.42) using similar argument without the process of } and far away from the corners Γ * f ∩ Γ +, * w and Γ +, * w ∩ Γ * ex , we obtain (3.2.42) using the standard method of freezing coefficients to the spherical coordinate representation of (3.2.13), (3.2.14) with integration by parts argument above. When 18) with α and C replaced by α − ε for any ε ∈ (0, 1) and C(||U ♯ || 2 using the standard method of freezing coefficients to (3.2.13), (3.2.14) with integration by parts argument above. Note that if we estimate the integral form of F ♯ using integration by parts argument above, then there is no singularity issue. From this result, we obtain DU ♯ ∈ C α−ε (N + f ∩ {θ ≤ 2θ 1 3 }). Combining these two regularity results for DU ♯ , we obtain DU ♯ ∈ C α−ε (N + f ). Using the regularity result for DU ♯ and U ♯ ∈ C 0 (N + f ), we obtain from (3.2.41) for a constant R 6 > 0. Using the arguments used when we obtained DU ∈ C α−ε (N + f ), we obtain ( From this result, we obtain the desired result. This finishes the proof. Using the scailing argument given in the proof of Proposition 3.1 in [2] with the results in Theorem 5.21 in [16] and Lemma 3.10, we can obtain the following result. We omit the proof. Lemma 3.13. Under the assumption as in Lemma 3.5, let U ♯ be a weak solution of (3.2.13), (3.2.14).
Finally, we prove that the C 2,α (−1−α,Γ + w ) (N + f ) solution of (3.2.13), (3.2.14) is of the form Ψ(r, θ)e ϕ . This statement is proved using the argument as in the proof of Proposition 3.3 in [4] (Method II). Although the arguments to prove this statement are almost same with those in the proof of Proposition 3.3 in [4] (Method II), since (3.2.13), (3.2.14) is different from the problem in Proposition 3.3 in [4] and similar arguments will be used later in the proof of Lemma 4.3, we present the detailed proof. Proof. Let U = U r e r +U θ e θ +U ϕ e ϕ be the C 2,α (−1−α,Γ + w ) (N + f ) solution of (3.2.13), (3.2.14). Then, where C and C * are positive constants depending on N + f and α and given in (3.2.12), respectively, and satisfies the following spherical coordinate representation of (3.2.13), (3.2.14)
Since the coefficients of (3.2.46), F and the boundary conditions in (3.2.47) are independent of ϕ, (U n r , U n θ , U n ϕ ) satisfies (3.2.46), (3.2.47) for all n ∈ N ∪ {0}. By this fact and the definition of (U * r , U * θ , U * ϕ ), (U * r , U * θ , U * ϕ ) satisfies (3.2.46), (3.2.47). Since (U * r , U * θ , U * ϕ ) is independent of ϕ, (3.2.46) satisfied by (U * r , U * θ , U * ϕ ) is given as  where Ψe ϕ : N + f → R 3 and Q en : Γ f → R are axisymmetric functions. Here, (3.3.1) is a transport equation whose coefficient is an axisymmetric and divergence-free vector field. Thus, the stream function of the coefficient vector field of (3.3.1) can be defined (see (2.4.1)). We find a solution of (3.3.1), (3.3.2) and obtain the regularity and uniqueness of solutions of (3.3.1), (3.3.2) using the stream function of the coefficient vector field of (3.3.1) and the solution expression given by using the stream function in the following lemma.
where c * is a positive constant depending on ρ + 0 , u + 0 , r s , r 1 , δ 1 and δ 6 , and ν f is the unit normal vector on Γ f pointing toward N + f . Suppose that Ψe ϕ : N + f → R 3 is an axisymmetric function in C 1,α (N + f ) satisfying
. This condition will be used to construct the stream surfaces of ∇ × ((Φ + 0 + Ψ)e ϕ ) in N + f in the proof of Lemma 3.15.
Note that by the facts that ∂ θ (V (f (θ), θ)) > 0 for θ ∈ (0, θ 1 ) and ∂ θ V (r 1 , θ) > 0 for θ ∈ (0, θ 1 ) obtained from (3.3.3), the surface θ = h (r ♯ ,θ ♯ ) (r) intersects with Γ f and Γ ex once, respectively. Collect θ = h (r ♯ ,θ ♯ ) (r) for Using this form of (3.3.1), it can be checked that 1) if Q is in C 1 (N + f ) and Q = constant on any curve on any level surface of V in N + f whose ϕ argument is fixed (in the case when a level surface of V is N + f ∩ {x = y = 0}, Q = constant on N + f ∩ {x = y = 0}), then Q is a solution of (3.3.1) and that 2) if Q is a C 1 solution of (3.3.1), then Q = constant on any curve on any level surface of V in N + f whose ϕ argument is fixed. Denote the θ-argument of the intersection points of Γ f and the level surface of V in N + f passing through x ∈ N + f by L(x). Since each level surface of V in N + f intersects with Γ f where the value of V on Γ f is equal to the value of V on the level surface, L(x) is given by Since Q en is axisymmetric, this can be written as Q = Q en (L).
Let Q = Q en (L). One can see that Q is a constant on any level surface of V in N + f and satisfies (3.3.2). Thus, by 1)

Estimate ||Q||
Thus, to estimate ||Q|| To estimate ||DQ|| Here, we used the definition of ∇ ⊥ V and the spherical coordinate expression of ν f . Using (3.1.1) and (3.3.3), it is easily seen that where c * * is a positive constant depending on c * , r s and δ 1 . Using (3.1.1) and (3.3.5), it is also easily seen that where Cs are positive constants depending on ρ + 0 , u + 0 , r s , α, δ 1 and δ 6 . By these three estimates, one can see that to estimate ||DL|| 0,α,N + f , it is enough to estimate || θ L || 0,α,N + f . To estimate ||DL|| 0,α,N + f , we prove the following claim.
In this proof, Cs and C i for i = 1, 2, . . . denote positive constants depending on the whole or a part of the data, δ 1 , δ 2 , δ 4 , δ 5 and δ 6 unless otherwise specified. Each C in different situations differs from each other.
3. Find f (0) ♯ using (3.4.13). By the definition of Ψ ♯, * e ϕ and (2.4.5), the transport equation in (3.4.14) can be written as From this form of the transport equation in (3.4.14), we see that the stream surface of the vector field ∇ × ((Φ + 0 + Ψ ♯, * )e ϕ ) in N + f ♯ is obtained by stretching or contracting the stream surface of the vector field ∇ × ((Φ + 0 + Ψ * )e ϕ ) in N + f * in r-direction. Using this fact, we obtain that the solution of (3.4.14) is given by . We denote this solution by T ♯ . By (3.1.19), T ♯ is expressed as Substituting this expression of T ♯ into the place of T in (3.4.13) using the fact that Π f ♯ f * (r 1 , θ) = (r 1 , θ), we obtain (L) := 1 r 1 sin θ 1 We find f (0) ♯ satisfying (3.4.15). For this, we estimate |(R)|.
Estimate |(R)|: With (2.3.3) and (3.0.1), we estimate (a) 2 , (a) 3 and Substitute this T ♯ into the place of T in (3.4.11) and (3.4.12). Then we obtain   Take M 1 = 2C 8 and σ 3 = min(σ (5) 3 , 3 ). And then define a map J from By the choice of M 1 and σ 3 , J is a map of P(M 1 ) into itself. Using the standard argument, one can easily check that J is continuous. Thus, the Schauder fixed point theorem can be applied to J . We apply the Schauder fixed point theorem to J . Then we obtain that there exists a fixed point (f ♭ (0),Ψ ♭ e ϕ ) ∈ P(M 1 ) of J .
In this proof, Cs denote positive constants depending on the data unless otherwise specified. Each C in different situations differs from each other.
By the choice of σ 3 , Ψe ϕ = Ψ i e ϕ for i = 1, 2 satisfy (3.3. 3) for f = f i (0) + f s . By this fact and ||f 2 − r s || for all t ∈ [0, 1] for some positive constantc * . Using arguments similar to the ones in the proof of Claim in Lemam 3.15, we can obtain and ||k −1 2 • (tV 2 + (1 − t)Ṽ 1 )|| 1,0,(f 2 (θ),r 1 )×(0,θ 1 ) ≤ C, we estimate the right-hand side of (3.4.34) in C β (N + f 2 ). Then we obtain in a similar way. As we do this,k 1 and k 2 play the role of V 2 andṼ 1 in the estimate of (e) in C β (N + f 2 ) and ϑ is regarded as the argument ofk 1 and k 2 . We have ||(d)|| 0,β,N + Note that if we change L 2 −L 1 in the way that we changed L 2 −L 1 in the estimate of 1 0 h ′ (tL 2 + (1 − t)L 1 )dt(L 2 −L 1 ) in the proof of Claim and estimate the resultant terms in With the above Claim, we estimate (a) in C β (N + f 2 ). Write (a) as .

Determination of a shape of a shock location
In the previous section, for given (ρ − , u − , p − , p ex , f ′ s ) in a small perturbation of (ρ − 0 , u − 0 e r , p − 0 , p c , 0), we found (f (0), Φe ϕ , L, S) satisfying all the conditions in Problem 2 except (2.5.5). In this section, to finish the proof of Theorem 2.16, for given (ρ − , u − , p − , p ex ) in a small perturbation of (ρ − 0 , u − 0 e r , p − 0 , p c ) as in the previous section or in a much small perturbation of (ρ − 0 , u − 0 e r , p − 0 , p c ) if necessary, we find f ′ s in a small perturbation of 0 as in the previous section such that a solution of Problem 3 for given (  ≤ σ}, For given (ρ − , u − , p − , p ex , f ′ s ) ∈ B (4) σ for σ ≤ σ 3 , let (f (0), Φe ϕ , L, S) be a solution of Problem 3 satisfying (3.0.2) given in Proposition 3.1. We define where τ f (0)+fs is the unit tangent vector on Γ f (0)+fs perpendicular to e ϕ and satisfying (τ f (0)+fs × e ϕ ) · ν f (0)+fs > 0 and ν f (0)+fs is the unit normal vector field on Γ f (0)+fs . Then A satisfies Cσ where C is a positive constant depending on the data. If for given (ρ − , u − , p − , p ex ) ∈ B (1) , then (f (0) + f s , Φe ϕ , L, S) satisfies all the conditions in Problem 2 and thus the existence part of Theorem 2.16 is proved. We will find such f ′ s using the weak implicit function theorem introduced in [3]. To apply the weak implicit function theorem, we need to prove that A is continuous, A is Fréchet differentiable at (ρ − 0 , u − 0 e r , p − 0 , p c , 0)(=: ζ 0 ) and the partial Fréchet derivative of A with respect to f ′ s at ζ 0 is invertible. We will prove these in the following lemmas.
We first prove that A is continuous.
Proof. The result is obtained by using the standard argument.
Next, we prove that A is Fréchet differentiable.  Proof. In this proof, we prove that A is Fréchet differentiable as a function of f ′ s at 0 with the other variables fixed at (ρ − 0 , u − 0 e r , p − 0 , p c ). The Fréchet differentiability of A as a function of (ρ − , u − , p − , p ex , f ′ s ) at ζ 0 can be proved in a similar way. In this proof, Cs denote positive constants depending on the data. Each C in different inequalities differs from each other. 1 rs or on a part of ∂N + rs by using Π rsfε where f ε := f ε (0) + εf s withf s = θ 0f ′ s . And then subtract the resultant equations from the same equations satisfied by (f 0 (0), Ψ 0 e ϕ , A 0 , T 0 ). Then we obtain , y is the cartesian coordinate system representing N + fε ,θ is θ coordinate for y, (r, θ) is (r, θ) coordinates for N + rs , and Π r fεrs is the r-component of Π * fεrs , ν fε is the unit normal vector field on Γ fε pointing toward N + fε andF 1 is F 1 changed by using the transformation Π fεrs . Divide the above equations by ε and formally take ε → 0 + using εf ′ s → 0 as ε → 0 + and (4.1.3). Then we have Substitute this into (4.1.14) again. Then we havẽ SubstitutingT (f ′ s ) given in (4.1.16) into (4.1.8) and (4.1.9), we get 0f s sin ζdζ andT (f ′ s ) are the partial Fréchet derivatives of f (0),Ψe ϕ ,Ã andT with respect to f ′ s at ζ 0 , respectively.
2. Show that A is Fréchet differentiable as a function of f ′ s at 0 with the other variables fixed at (ρ − 0 , u − 0 e r , p − 0 , p c ).
2. Show that there exists a unique q ∈ H 1 0 ((0, θ 1 ), sin θdθ) satisfying (4.1.24) for all ξ ∈ H 1 0 ((0, θ 1 ), sin θdθ). For given f ∈ L 2 ((0, θ 1 ), sin θdθ), we consider From this fact, we see that U satisfies for all ξ ∈ H 1 0 (D). Using this fact and the fact that the coefficients of ∆ S 2 U in the spherical coordinate system are independent of ϕ, we apply arguments similar to the ones in the proof of Lemma 3.14 to U (here, we use the facts that a bounded sequence in H 2 (D) contains a weakly convergent subsequence and a H 1 0 (D) function satisfying (4.1.26) for all ξ ∈ H 1 0 (D) is unique). Then we have that U only has the form u ϕ (θ)e ϕ .
From this, we can conclude thatΨ m l e ϕ and D(Ψ m l e ϕ ) weakly converge toΨ (f ′ s ) e ϕ and D( Sincef ′ s = 0 by the assumption, there exists k ∈ N such that c k > 0 or c k < 0 in the expression off ′ s in (4.1.32). Without loss of generality, assume that c k > 0 for some k ∈ N. Then since λ k > 0 by Lemma 4.3, p k , that is, the solution of (4.1.35), (4.1.36) for j = k, satisfies p k (r 1 ) > 0. Using this fact, p k (r s ) = 0 and the form of (4.1.35), we can deduce that p k ≥ 0 in [r s , r 1 ]. Thus, p ′ k (r s ) ≥ 0. Write (4.1.29) in the form for all ξ ∈ L 2 ((0, θ 1 ), sin θdθ). Rewrite this as (4.1.40) , for a sufficiently large l such that m l ≥ k, if we take ξ = q k , then the left-hand side of (4.1.40) becomes a negative number (here we used the trace theorem). This contradicts to the assumption that D f ′ s A(ζ 0 )f ′ s = 0. This finishes our proof.
Applying the weak implicit function theorem introduced in [ where C is a positive constant depending on the data. Let σ 1 be a positive constant ≤ σ 1 and to be determined later. Suppose that there exist two solutions (f i , Φ i e ϕ , L i , S i ) for i = 1, 2 of Problem 2 for σ ≤ σ 1 satisfing (4.2.1).

Appendix
In this section, we present some computations done by using the tensor notation given in §3.2. We explain how we transformed an elliptic system in the cartesian coordinate system into a system in the spherical coordinate system in the proof of Lemma 3.9 and Lemma 3.10, and show that (3.2.8) is equivalent to (3.2.1).
Let (q 1 , q 2 , q 3 ) be an orthogonal coordinate system in R 3 . The unit vectors in this coordinate system in the direction of q i for i = 1, 2, 3 are given as 1 h i ∂x ∂q i (=: e q i ) for i = 1, 2, 3 where x = xe 1 + ye 2 + ze 3 and h i := | ∂x ∂q i |. By ∇ = 3 i=1 eq i h i ∂ q i , ∇U where U : R 3 → R 3 can be written as Using this relation and the relation div(a ⊗ b) = ∇ab + adivb, we compute the left-hand side of (3.2.8). Then we get a(r)∆(Ψe ϕ ) + ∇(Ψe ϕ )∇a(r) − ∇e ϕ b(r)∂ r Ψe r − e ϕ div(b(r)∂ r Ψe r ) − ∂ r ρ +