{'fraction_non_alphanumeric': 0.049766503633700845, 'fraction_numerical': 0.02116665285031363, 'mean_word_length': 3.748097612175282, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 17, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 6, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Compressible perturbation of Poiseuille type flow

The paper examines the issue of stability of Poiseuille type flows in regime of compressible Navier-Stokes equations in a three dimensional finite pipe-like domain. We prove the existence of stationary solutions with inhomogeneous Navier slip boundary conditions admitting nontrivial inflow condition in the vicinity of constructed generic flows. Our techniques are based on an application of a modification of the Lagrangian coordinates. Thanks to such approach we are able to overcome difficulties coming from hyperbolicity of the continuity equation, constructing a maximal regularity estimate for a linearized system and applying the Banach fixed point theorem.


Introduction
The mathematical description of compressible flows is important from the point of view of applications, domains such as aerodynamics and geophysics are the most natural to be mentioned here. On the other hand, complexity of the equations describing the flow delivers very interesting mathematical challenges. In spite of active research in the field, we are still far from the complete mathematical understanding of compressible flows. The only general existence results are available for weak solutions with homogeneous boundary conditions [7], [12]. As far as regular solutions are concerned, we have so far only partial results assuming either some smallness of the data, or its special structure. The problems have been investigated mainly with homogeneous boundary conditions ( [4], [17]). For the overview of the state of art in the theory one can consult the monograph [18].
From the point of view of the aforementioned applications it seems very important to investigate the problems with large velocity vectors, which lead in a natural way to inhomogeneous boundary conditions. Due to the hyperbolic character of the continuity equation the density must be then prescribed on the inflow part of the boundary. Existence issues for such inflow problems are investigated in [9], [10], [19], [20], [21], [22] and [26]. The mentioned group of problems can be regarded as questions of stability of particular constant flows.
In the present article we would like to examine the issue of stability of Poiseuille type flow in pipe-like domain in compressible regime. The unperturbed flow is a solution to the compressible Navier-Stokes system with homogeneous slip boundary conditions for given constant friction, constant density and constant external force (gravitation-like term). The Poiseuille flow is a special symmetric solution to the incompressible Navier-Stokes equations in cylindrical domains. Here it is viewed as a solution to the compressible Navier-Stokes system with constant density and constant external force, parallel to axis of the cylinder, given by the pressure. Hence the pressure, unknown in the incompressible model, is recognized as a given force. Such change of 'observer' looks acceptable from the mechanical point of view. Thanks to that interpretation we obtain a natural physically reasonable flow in compressible regime. The mathematical objective of this article is to establish stability of such flow under some structural assumptions limiting the magnitude of admissible perturbations.
Let us define the system. We consider steady flow of a viscous, barotropic fluid in a bounded, cylindrical domain in R 3 , described by the Navier -Stokes system supplied with inhomogeneous Navier slip boundary conditions. The complete system reads ρv · ∇v − µ∆v − (µ + ν)∇div v + ∇π(ρ) = ρ F in Ω, div (ρv) = 0 in where Ω = [0, L] × Ω 0 with Ω 0 ⊂ R 2 of class C 2 , Γ denotes the boundary of Ω (see Fig.1), v is the velocity field of the fluid, ρ its density, µ and ν are viscosity constants satisfying µ > 0 and (ν + 2µ) > 0, f ≥ 0 is the friction coefficient which may be different on different components of the boundary Γ, π = π(ρ) is the pressure given as a function, at least C 1 , of the density and F is an external force. T denotes the Cauchy stress tensor of the form T(v, π) = 2µDv + νdivvId − πId where D = 1 2 (∇v + ∇v T ) is the symmetric gradient. Next, n and τ k are outer normal and tangent vectors to ∂Ω. Boundary data ρ in , b, d will be discussed later. The boundary Γ is naturally split into three parts: (1.2) Thanks to the chosen geometry of the domain, the above decomposition is easily illustrated by the figure 1.
We shall say few words about the physical interpretation of the system (1.1), in particular about the choice of boundary conditions (1.1) 3,4 . We would like to model a flow through a pipe. We assume that the fluid obeys Navier slip conditions on the walls of the pipe (Γ 0 component of the boundary), hence natural conditions on Γ 0 are d ≡ 0 and b k ≡ 0. However, the mathematical requirements impose a need to prescribe the boundary conditions on Γ in and Γ out . From the physical viewpoint these parts are artificial, this is the area where the parameters of the velocity and density are measured. This gives us a freedom of choice of the type of boundary conditions on the inflow and outflow part, which can be fit to the mathematical approach, hence we choose inhomogeneous slip condition. Note that as the friction coefficient goes to infinity, then the relations (1.1) 3,4 , at least formally, become the standard Dirichlet conditions describing the whole velocity vector at the boundary. Since the velocity does not vanish on the boundary, the hyperbolicity of the continuity equation impose a need to prescribe the density on the inflow part, which lead to the condition (1.1) 5 . The velocity field determines the characteristics of the continuity equation and in particular the total mass Ω ρ dx is determined implicitly by (1.1) 5 .
Our goal here is to analyze a perturbation of the Poiseuille type flow where x 1 points the direction of the axis of the cylinder. It is one of the classical examples of laminar flows satisfying the incompressible Navier-Stokes equations in cylindrical domains. In the classical literature the flow is considered with homogeneous Dirichlet condition on the boundary. V P is then found as a solution to corresponding elliptic problem with Dirichlet boundary condition on each x 1 -cut of Ω. Some explicit formulas on V P in certain domains are well-known ( [8], [11]). In the case of slip boundary conditions that are subject of our analysis in this paper the flow (1.3) can be also found on each cut of the cylinder as a solution to elliptic problem with corresponding boundary condition (see Lemma 1 below). In certain domains it can also be expressed with explicit formulas (see [14] and the example below). Since we are interested in a general cylindrical domain, we will not have such formula but we show that the solution of the form (1.3) exists provided that Ω 0 is sufficiently regular and, under the slip boundary conditions (1.1) 3,4 , V P does not vanish on the boundary.
As we already said, the pressure in the Poiseuille flow can be regarded as an external force parallel to the axis of the cylinder. Like in the classical Poiseuille flow, we assume the pressure to be a linear function of x 1 . Hence it is natural to assume the formπ = ω(f )x 1 , where ω(f ) is a negative constant (the sign describes the directon of the flow). Then we see that V P , the first coordinate ofV = (V P , 0, 0), is a solution to the elliptic problem Now to prove Lemma 1 it is enough to apply the maximum principle to the system (1.7). We show the proof in the Appendix, at this stage we should have a closer look at the dependence ω(f ). This dependence can be justified considering the compatibility condition for the system (1.7), that reads Ω 0 ω dx = f ∂Ω 0 V P dσ. As we will see, the dependence ω(f ) determines the assumptions we will have to make on the viscosity. Note that for ω = 0 the only solution is V P = 0, provided f = 0. For f = 0 (perfect slip) we should obtain a constant flow, hence we put ω = 0 for f = 0. The linear structure of (1. 7) makes ω(f ) continuous. Thus, we conclude ω(f ) → 0 for f → 0, and hence by (1.6) This is an important conclusion since in order to show the energy estimate we have to control ∇V with the viscosity, and so we allow the viscosity to be low provided that the friction on Γ 0 is small, what is a realistic assumption (see also the remarks after the formulation of Theorem 1). Let us illustrate the dependence ω = ω(f ) with the following example.
Example. Take Ω 0 = B(0, 1) ⊂ R 2 and µ = 1. Then, due to axial symmetry of the domain, it is natural to look for We require that ∆v = v rr + 1 r v r depend only on f . Moreover, we expect to obtain a constant flow for f = 0 (perfect slip) and classical Poiseuille profile for f = ∞. The above considerations lead to the family of solutions where k f = (π−2)f +4π 2 and T F is the flux of the flow v f through Ω 0 .
For a perfect slip case f = 0 we obtain a constant flow v 0 = 2 T F k f and for a no-slip case f = ∞ we get a classical Poiseuille profile v ∞ = T F (1 − r 2 ). On the boundary we have f +k f what is a strictly positive constant. Finally, In particular ω(f ) < 0 for f > 0 and ω(0) = 0. Before we formulate our main result, we need one observation concerning the boundary conditions. Note that, since V P is found on every x 1 -cut of Ω, we can impose the boundary conditions (1.4) 3,4 only on Γ 0 . On the other hand, in order to define small perturbations as a solution to (1.1) we have to measure the distance (in appropriate norms) between the solution to (1.1) and the Poiseuille flowV . Hence we need to consider the traces of the quantities from the boundary conditions of (1.1) with the functionV instead of v. Since our analysis acts on a finite cylinder we define these traces at the bottoms Γ in ∪ Γ out : where tr denotes the trace operator. By (1.4), n · T(V ,Π) · τ k + fV · τ k = 0 and n ·V = 0 at Γ 0 .
We are now in a position to formulate our main result. To this end it is convenient to define the quantity which measures the distance of the data from the Poiseuille flow: The main result of the paper reads Theorem 1. Assume that the boundary data of (1.1) is close to the Poiseuille flow (V ,Π), more precisely, let D 0 defined above be small enough. Assume that the friction f is large enough on Γ in and p > 3. Assume further that there exists a constant κ > 0 such that This solution is unique in the class of small perturbations of (V ,ρ).
Let us make some remarks concerning our main result. The condition on the viscosity (1.12) seems to be a serious constraint, but as we will see from the proofs we just need it to control the gradient of the Poiseuille flow, what yields this assumption natural (see also the remark in the proof of Lemma 6). We recall that ∇V depends on the friction f on Γ 0 , and in particular (1.8) holds. It follows that for small values of friction on Γ 0 it is enough to assume that the viscosity is large enough, but only compared to the friction. This assumption is reflected in the condition (1.12). Theorem 1 admits the case of perfect slip f | Γ 0 = 0, and in such case (1.12) reduces to µ > 0, so no lower bound on the viscosity is required. In this case there is no bound on the size ofV . However in this caseV would be a constant flow. We shall recall that the friction at Γ in is chosen independently to f at Γ 0 . From the point of view of modelling, the data at Γ in is given, hence it is important to focus the attention at Γ 0 . The assumption p > 3 is required for , it is required to control the pointwise boundedness of ∇v and the density.
Let us explain the main idea of the proof. We will follow an idea of Lagrangian type coordinates [5], [13], [16], [24]. The continuity equation is of hyperbolic type and contains a term u·∇w (where u and w are perturbations to the velocity and density introduced in the next section), which makes serious troubles for the issues of existence in case of inhomogeneous boundary conditions, see [9], [21], [18], [19], [20]. Here we overcome this obstacle by changing the system of coordinates in such a way that this term disappears (2.7). We obtain a more complex system but with structure suitable for an application of the Banach fixed point theorem. On the other hand our solutions are regular enough, thus we are able to go back to the original system keeping the well posedness of the original model. Our approach works since we are equipped with the maximal regularity estimate for a linearization of the equations in the Lagrangian coordinates -Theorem 2. This tool gives a complete control of the regularity of solutions.
The rest of the paper is organized as follows. In Section 2 we introduce the perturbations as unknown variables obtaining the system (2.5). Next we introduce the Lagrangian-type coordinates that lead to the system (2.11) and we derive the necessary estimates for the Lagrangian transformation. In Section 3 we deal with the linearization of (2.11). For the linear system we show the estimate in W 2 p (Ω) × W 1 p (Ω). It is given by Theorem 2. The first step is the energy estimate (3.3). Then we consider the vorticity of the velocity and the Helmholtz decomposition to reduce the continuity equation to a sort of transport equation (3.22) that enables us to find the bound on w W 1 p (Ω) . This result together with the properties of the Lamé system lets us conclude Theorem 2. In the second part of this section we apply the estimates to solve the linear system and hence show that T given by (2.27) is well defined. In Section 4 we show the contraction principle for T . To this end we consider the system for the difference of two solutions and write it in a form (4.1) which has a structure of (3.1). The contraction results from the estimate (3.30) and bounds on the norms on of the r.h.s. of the system for the difference. At the end of Section 4 we apply the Banach fixed point theorem to solve the system (2.11) and conclude the proof of Theorem 1.
Let us finish this introductory part with some remarks concerning notation. By C we shall denote a constant that is controlled, but not necessarily small. E shall denote a constant that can be arbitrarily small provided the data is small enough. Sometimes we will write E(·) to underline that we need the smallness of certain quantity. The functional spaces on Ω will be denoted without the symbol of the set, for example we will write W k p instead of W k p (Ω) for standard Sobolev spaces of functions intergable with the p-th power with derivatives up to order k, W 1−1/p p (∂Ω) denotes the Slobodeckij spaces, defining regularity of traces from W 1 p (Ω), [1]. Finally, we will need to consider the density in the space L ∞ (0, L; L 2 (Ω 0 )). For simplicity we denote it as L ∞ (L 2 ). We do not use different notation for scalar and vector valued functions, while matrix valued functions are written in bolded font. The coordinates of a vector are denoted by (·) , i.e. u = (u (1) , u (2) , u (3) ).

Preliminaries
In this section we introduce perturbations of the Poiseuille flow (V ,ρ) as unknown variables, what leads to the system (2.5). Then we introduce a change of variables that straightens the characteristics of the continuity equation. We obtain the system (2.11). The simplified form of the continuity equation in this Lagrangian framework makes it possible to apply the Banach fixed point theorem to the system (2.11).

Reformulation of the problem
We come back to the main system (1.1). Since we are interested in solutions that are small perturbations of (V ,ρ), it is convenient to consider the perturbations as unknown functions. For technical reasons it is better to have u · n = 0 on the boundary. Hence we start introducing . It can be found as u 0 = ∇φ where φ solves a Neumann problem. We assume that d −V · n W 1−1/p p (Γ) is small enough for to hold for some θ 1 > 0. In fact this is not really a restriction as we consider small perturbations of (V ,ρ) and in Lemma 1 we have shown that V P is separated from zero. Now we take In particular we want our perturbed flow v to have the first component also separated from zero. This is quite natural constraint if we consider small perturbations ofV . With the above definition of u this constraint reads Next we introduce the perturbation of the density (recall thatρ ≡ 1): where γ = π (1) and From now on we focus on the system (2.5). Notice thatV · ∇V = 0, hence the term −w(u + u 0 +V ) · ∇u + u 0 +V ) is a higher order term and so the form of F and G implies immediately the following lemma.
Lemma 2. Let F (u, w) and G(u, w) be given as above, then where E denotes a small, compared to µ, positive constant.

Change of variables
With our smallness assumptions it is quite natural to solve (2.5) with a fixed point argument. However, a direct application of this method fails because of the nonlinear term u · ∇w in the hyperbolic continuity equation. The idea to overcome this problem is to introduce a change of variables such that this awkward term vanishes. We look for the appropriate transformation as x = ψ u+u 0 (z) satisfying the identity In the following lemma we construct the mapping ψū for arbitrary functionū small in W 2 p with vanishing normal component on the boundary Γ 0 . Lemma 3. Let ū W 2 p be small enough andū · n| Γ 0 = 0. Then there exists a diffeomorphism x = ψū(z) defined on Ω such that Ω = ψū(Ω) and (2.7) holds with u + u 0 =ū.

Proof.
A key point in the proof is the fact that V P ≥ c > 0. In particular we are able to divide (2.7) by V P obtaining whereũ =ū V P . Sinceū, V P ∈ W 2 p we haveũ ∈ W 2 p and, since we are interested in small perturbations we can assume that Now we can follow the proof from [20] and look for ψ(z 1 , z 2 , (2.9) Due to (2.8) we solve (2.9) for (z 2 , z 3 ) ∈ Γ in following [20] and show that there exists a set Ωū such that ψ(Ωū) → Ω is a diffeomorphism. It remains to show that Ωū = Ω. To this end we examine the derivatives of ψ. We have Now take x n → x 0 ∈ Γ 0 and z n = φ(x n ). Then we have The latter is parallel to Γ 0 sinceũ · n| Γ 0 = 0, hence we have φ(Γ 0 ) = Γ 0 (precisely we say it in a sense of tangent spaces since φ is defined only on Ω). To examine the behavior of tangent vectors on Γ out notice that whereτ i are given by appropriate entries of Dψ, the important fact is that the first coordinate vanishes. Hence if we take x n → x 0 and z n = φ(x n ), this time with x 0 ∈ Γ out , then what is parallel to Γ out . It shows that φ(Γ out ) = Γ out . Since φ(Γ 0 ) = Γ 0 by the definition of φ, we conclude that Ωū = Ω and complete the proof. Now we proceed with transformation of the system (2.5). As a vector fieldū satisfying the assumptions of the Lemma we take u + u 0 where u is the solution of (2.5). So far we don't know if this solution exists, our goal is to show its existence. Hence our approach can be regarded as working in a kind of Lagrangian coordinates [16]; assuming that the solution exists we rewrite the system in the new variables ψ u+u 0 induced by u through (2.7). Then in the new coordinates we hope to be able to apply a fixed point method to show the existence of a solution. Since (u, w) are perturbations that are assumed to be small, we can assume that u W 2 p is small enough that the assumptions of Lemma 3 are satisfied. Hence the solution gives a well defined transformation x = ψ u+u 0 (z), that we denote for simplicity by ψ. If we show in addition its uniqueness in a class of small perturbations, then denoting φ = ψ −1 we have ψ(Ω u ) = Ω and we come back to the original coordinates where our solution solves (2.5). Rewriting the system (2.5) in coordinates z yields (2.11) The functionsF andG involve commutators from the change of variables. More precisely, Here the first variable in the commutator R(·, ·) denotes a function and the second is a differential operator. For example, R(w, ∇) := ∇ x w − ∇ z w and its i-th coordinate reads We shall not give here precise formulas for the other commutators. Instead, we are now ready to give some heuristic arguments that will show what regularity we expect from the change of variables φ. To this end note that the commutators of the operators of order k depends on the derivatives of φ up to order k. More precisely, commutators of order one contain only components of the form ∇ z u · ∇ x φ, while second-order commutators contain the terms and the terms of lower order. Hence in order to find the estimates on F (u, w) Lp and G (u, w) W 1 , what will be satisfied provided that φ ∈ W 2 p due to the imbedding W 1 p ∈ L ∞ (recall that p > 3). We should also note that, for simplicity of notation, the functions F (u, w) and G(u, w) in (2.12) and (2.13) denote exactly the same quantities as before. Hence we should keep in mind that they contain differential operators and so now they also contain some commutators that we will have to control to repeat the estimate (2.6).
In order to show such estimate we should have a closer look at the derivatives of φ. With our construction of ψ = φ −1 , it is easier to consider first the derivatives of ψ that are given by (2.10) withũ = u+u 0 V P where u is the solution to (2.5), not to be confused withũ =ū V P from the proof of Lemma 3 whereū was an arbitrary function. In order to find the bounds on the commutators we will need the following smallness results for our change of variables:

17)
where E is sufficiently small number comparing to µ, depending on norms of perturbations measured by D 0 defined in (1.11).
Proof. The core of the proof is in the imbedding W 1 p ⊂ L ∞ . We start with (2.14). We estimate L p norm of the entries of E (2.10). These results quite directly from the form of E, but needs certain attention as E depends on ψ implicitly. The entries of E without integrals will be small provided that ψ is bounded what obviously holds true. For the entries involving integrals we change the order of integration and derivative obtaining By Jensen inequality we have Integrating the last inequality over Ω we get The smallness ofũ in W 2 p and the imbedding W 1 p ⊂ L ∞ gives (2.14). To show (2.15) we differentiate the entries of E, let us focus on entries with integrals. We have (we omit the sum over k): Again by Jensen inequality, and Integrating the above inequalities over Ω we arrive at The bound on F (u, w) Lp + G(u, w) W 1 p results from (2.6) in Lemma 2 (there are also some commutators since F and G involve differential operators, but these can be estimated as follows). We briefly justify the bounds on the commutators inF . To start with, the first order commutators contain the given functionV , the functions u and w and the derivatives of φ, but only of the form Hence applying Lemma 4 we get The second order commutators contain the second order derivatives of φ and the first order derivatives only of the form (2.22). Hence the application of Lemma (4) yields and we conclude the bound on F Lp . In order to estimate G (u, w) W 1 p we differentiate the commutator Applying again Lemma 4 we get and the proof is complete. From now on we focus on the system (2.11) instead of (2.5). It is of crucial importance for us that we can solve (2.11) in the domain Ω, what results from our choice of the transformation ψ. Now we are in a position to define the operator to which we want to apply the Banach fixed point theorem in order to solve the system (2.11). Namely, we set (u, The point is that the term ∂ x 1 w + u · ∇w is replaced by ((V P • ψū) + u (1) )∂ z 1 w, and for this term we find a bound in W 1 p , what is necessary to show the contraction property of T . From now on for simplicity we will write ψū instead of ψū +u 0 .

A priori bounds and solution of the linear system
In this section we deal with the linear system: with given functions F, G, B k ,ū. Note that we take the superposition of V P and ψū to obtain V P •ψ u for the original system and hence V P when we pass to the original system of coordinates. The same remark concerns the function F , but it does not change anything in the computations. We need to solve the system (3.1) to show that T is well defined by (2.27). To this end we need the appropriate estimates that we show in the first part of this section. In the second part the linear system is solved.

A priori bounds
In this section we show the estimate in W 2 p × W 1 p for the solution of the linear system (3.1) in the maximal regularity regime. The first step is the energy estimate. It is given by the following Lemma 6. Let (u, w) be a solution to the system (3.1) with given (F, G, B,ū) for some θ 3 > 0. Assume that the friction f is large enough on Γ in and the viscosity µ and the friction on Γ 0 satisfy (1.12). Then Proof. We start with two basic observations. First, since V P x 1 = 0, by (2.7) we have . Now recall that the constant θ in (1.5) is independent of the smallness of perturbation. Hence we can assume the ū W 2 p is small compared to θ and by the imbedding W 1 p ⊂ L ∞ we have We keep in mind that (1.12) and (1.6) hold, what implies that (3.5) will be controlled by the viscosity, more precisely even by a constant which decreases with increasing viscosity due to (1.6). In the remaining of this section we will write V P instead of V P • ψū. The fact that we consider the superposition does not influence the computations as we have (3.5). We apply the identities and (3.7) Now we multiply (3.1) 1 by u and integrate. Using the above identities, with application of the boundary conditions (3.1) 3,4 we arrive at (3.8) Note that by (3.5) we have The other terms on the r.h.s. are all 'good' terms. The Γ out term on the l.h.s. is nonnegative and the Γ in term will be positive for f large enough. To deal with the term Ω w div u dx we apply the continuity equation to express div u obtaining: By (3.2) the integral over Γ out will be nonnegative, hence we have To derive the W 1 2 -norm of u we apply the Korn inequality [23,27]: where C K = C K (µ, ν, f, Ω) and C K is increasing with µ. A sketch of the proof of (3.10) one can find in the Appendix, note that in (3.10) we use only information at Γ in , a part of the boundary, but still it is sufficient to control the whole norm of W 1 2 . Combining (3.8), (3.9) and (3.10) we get We have to deal with the last term on the r.h.s. It is impossible to show it has a good sign, hence the only way is to estimate it directly with where C P is the constant from the Poincaré inequality in V . Inserting the above to (3.11) we get where we recall that E is a small constant and C is a data-dependent constant, not necessarily small. Now, C K is increasing with µ, while C P does not depend on µ. Moreover, (1.6) implies that E(ū, ∇V P ) will be decreasing when µ increases. Finally, (1.8) implies that for small values of the friction f on Γ 0 it is enough to assume that the viscosity is large only compared to f | Γ 0 to control E(ū, ∇V P ). We conclude that the constant on the l.h.s. of (3.13) will be positive provided that µ satisfies (1.12).
Here it is a good point to emphasize the necessity of sufficient magnitude of the viscosity coefficient for large f at Γ 0 . This assumption is somehow natural, although in the case V P = const it is not required [19]. We have to control (3.12), and largeness of dissipation may only give us this chance. Note that for the Dirichlet boundary condition [9], although the constant flow is considered, such assumption is required, too.
To complete the proof of (3.3) we find a bound on w L∞(L 2 ) . To this end we refer to the next subsection, where we solve the linear system. Notice that w = S(G − div u), where S is defined in (3.35), and so (3.14) Combining (3.13) and (3.14) we conclude (3.3).
In the next step we show higher bound on the vorticity of the velocity. To this end we take the vorticity of (3.1) 1 . Denoting α = rot u we get The boundary conditions (3.15) 2,3 are derived from differentiation of (3.1) 4 in tangential directions and application of (3.1) 3 , see [15], [19]. The above system gives the estimate ( [27], Theorem 10.3 with µ = 0): Applying the interpolation inequality (5.7) to u W 1 p and the energy estimate we get for any > 0. Now consider the Helmholtz decomposition of the velocity where ∂φ ∂n | Γ = 0 and div A = 0. We see that the field A satisfies the system rot A = α in Ω, div A = 0 in Ω, A · n = 0 on Γ.

(3.18)
This is the classical rot-div system and from [23] we have A W 2 p ≤ C α W 1 p , what by (3.16) can be rewritten as for any > 0. Now we substitute the Helmholtz decomposition to (3.1) 1 . We get Since ∆φ = div u, we can write We underline that we are now at the level of a priori estimates and (3.21) should be treated as the definition ofH. In fact we can think ofH as a kind of effective viscous flux like in the theory of weak solutions to compressible Navier-Stokes equations ( [7], [12]). Now (3.20) can be rewritten as ∇H =F . Combining the last equation with (3.1) 2 we arrive at The equation (3.22) makes it possible to estimate w W 1 p and ∂ z 1 w W 1 p in terms of H. Next we can find the bound on H using interpolation and the energy estimate. The first step is in the following lemma: Lemma 7. Let w solve (3.22) with H ∈ W 1 p and w| Γ in = w in . Then Proof. To estimate w Lp we multiply (3.22) by |w| p−2 w and integrate. Using the boundary conditions we getγ The boundary term on the l.h.s. is positive and the constant in the last term on the r.h.s. is small (note that we take only the z 1 derivative of V P ). Hence the above implies In order to find a bound on ∂ z i w we differentiate (3.22) with respect to z i . If we assume that w ∈ W 1 p then (3.22) implies ∂ z 1 w ∈ W 1 p , since W 1 p is an algebra. Thus we differentiate (3.22) with respect to z i , multiply by |∂ z i w| p−2 ∂ z i w and integrate. We have Proof. To close the estimate (3.30) it remains to find the bound on u W 2 p . To this end notice that in particular u satisfies the Lamé system: (3.31) Lemma 11 applied to the above system yields Applying the interpolation inequality to the term u W 1 p and then the energy estimate (3.3) we get Combining this estimate with (3.24) and (3.29) with appropriate δ we conclude (3.30).
Note that Theorem 2 gives more than we need to solve (3.1), namely the bound in W 1 p of ∂ z 1 w. We will use this result to show the contraction property for the operator T in Section 4.

Solution of the linear system
With the estimates that we obtained we are ready to solve the system (3.1). First we define the weak solution and show its existence. Next, applying the estimate (3.30) we show its regularity under the appropriate regularity of the data.

Weak solution
By the weak solution to (3.1) we mean a couple (u, w) ∈ V × L ∞ (L 2 ) such that (3.34) To find the weak solution we apply the Galerkin method. Hence we introduce an orthonormal basis of ω k ⊂ V and finite dimensional spaces We look for the approximations of the velocity of the form u N = N i=1 c N i ω i . Taking into account the continuity equation we have to define the approximations of the density in an appropriate way. Namely, we set w N = S(G N − div u N ), where S : L 2 (Ω) → L ∞ (L 2 ) is defined as and satisfies the estimate The construction of S is quite straightforward. For a continuous v we set (3.37) Next we show directly the estimate (3.36), which enables us to extend S on L 2 (Ω) using a standard density argument. Now we proceed with the Galerkin scheme.
where F N and G N are orthogonal projections of F and G on V N , we arrive at a system of N equations (3.39) Now, if u N satisfies (3.38) for k = 1 . . . N and w N = S(G N − div u N ), then a pair (u N , w N ) satisfies (3.33) -(3.34) for (v, φ) ∈ (V N ×C ∞ (Ω)), φ| Γout = 0. We will call such a pair an approximate solution to (3.33) -(3.34). To solve the system (3.38) we apply the following well-known result (the proof can be found in [25]): Lemma 8. Let X be a finite dimensional Hilbert space and let P : X → X be a continuous operator satisfying ∃M > 0 : (P (ξ), ξ) > 0 for ξ = M. Then ∃ξ * : ξ * ≤ M and P (ξ * ) = 0.
We define P N : V N → V N as In order to apply Lemma 8 we show that P (ξ N ), ξ N > 0 on some sphere in V N with radius dependent on the norms of the data. To this end we follow the proof of the energy estimate for (3.1). This is in fact standard approach in the Galerkin method: the energy estimate combined for some (u, w) ∈ W 1 2 ×L ∞ (L 2 ). We easily to verify that (u, w) is a weak solution. First, passing to the limit in (3.33) for (u N , w N ) we see that u satisfies (3.33) with w. On the other hand, taking the limit in (3.34) we verify that w = S(G − div u). We conclude that (u, w) satisfies (3.33) -(3.34), thus we have the weak solution. To show the boundary condition on the density we can rewrite the r.h.s of (3.35) as and, treating x 1 as a 'time' variable, adapt Di Perna -Lions theory of transport equation ( [6]) that implies the uniqueness of solution to (3.42) in the class L ∞ (L 2 ) (note that this is the reason we work with weak solutions with the density in L ∞ (L 2 ) instead of usual L 2 ). This completes the proof of existence of weak solution.

Strong solution
The following result gives strong solution to the linear system (3.1) for the data of appropriate regularity.
Assume further that f is large enough on Γ in and µ fulfills (1.12). Then there exist a unique solution (u, w) ∈ W 2 p × W 1 p to the system (3.1) and the estimate (3.30) holds.
Proof. To show appropriate regularity of the weak solution for the regular data it is enough to apply the estimate (3.30) provided that we handle the singularities of the boundary at the junctions of the wall Γ 0 with inlet Γ in and outlet Γ out . To this end we apply the result on the elliptic regularity of the Lamé system with slip boundary conditions, Lemma 11 in the Appendix. Notice that we can apply this method since we work in the fixed domain Ω due to the appropriate choice of the change of variables. Otherwise we would end up in a free boundary problem and the solution of the linear system would be highly nontrivial.

Contraction
andR denotes all the differences between corresponding commutators inF . We can afford using such abbreviation and estimateR without any additional computation if we just notice that the commutators are linear with respect to the functions and hence if certain estimate in terms of the function holds for a commutator, then the same estimate in terms of the difference hold for the difference of the commutators, for example can be dealt with using symmetry arguments. But first we have to reformulate slightly the system (5.3). Having the weak solution, and the bound in W 1 2 for the velocity, we are able to consider only the case f ≡ 0 and B ≡ 0. The general case will be a consequence of this particular one. Assume that Γ in ⊂ {x 1 = 0}, then we define the operator E v as extending a vector field defined for {x : x 1 ≥ 0} on the whole space as and on the plane x 1 = 0 the extension E v as preserves the slip boundary conditions for f ≡ 0 and 3. Now we are allowed to localize the equations in a vicinity of Γ in obtaining a system in a smooth domain. Then the theory [2,3] gives the full regularity of E v as (v) in neighborhood of the junctions. Application of analogous antisymmetric extension on Γ out completes the proof. Proof. The interpolation inequality in L p ([1], Theorem 2.11) and the imbedding W 1 p ⊂ L ∞ for p > 3 yields what entails (5.7) after application of the Cauchy inequality.
Proof. The known result based on properties of the kernel of D(·) [11,28] yields In order to obtain (5.8) we shall prove that u L 2 ≤ 1 2C 2 ∇u L 2 + M ( D(u) L 2 + u| Γpart L 2 (Γpart) ). (5.10) Compactness argument and features of D implies (5.10). We use that fact that the only solution to the system D(u * ) = 0, u| Γpart = 0 is u * ≡ 0. To show it we use the fact that if D(u) = 0 then Hence u is an affine map and it is easy to verify verify that if at least one of the coefficients A, B, C is nonzero then the rank of matrix of u is 2, hence Ker u is a line. On the other hand, Γ part is a two dimensional submanifold as it is an open, regular subset of ∂Ω. But Γ part ⊂ Keru, hence we conclude that A = B = C = 0, hence d i = 0.