Multicomponent reactive ﬂows: Global-in-time existence for large data

Multicomponent reactive flows arise in many physical applications 
in sciences and engineering. The objective of this work is 
to develop a rigorous mathematical theory based on the principles of 
continuum mechanics.


Introduction
Multicomponent reactive flows arise in many practical applications such as combustion, atmospheric modelling, astrophysics, chemical reactions, mathematical biology etc.The objective of this work is to develop a rigorous mathematical theory based on the principles of continuum mechanics.Accordingly, the physical systems considered in this paper are characterized by the state variables: the total mass density = (t, x), the velocity field u = u(t, x), the absolute temperature ϑ = ϑ(t, x), and the species mass fractions Y k = Y k (t, x), k = 1, ..., n, depending on the time t ∈ (0, T ) and the Eulerian spatial coordinate x ∈ Ω ⊂ R 3 .
The balance equations associated with multicomponent reactive flows express the conservation of mass, momentum, energy, and conservation of species mass (see Giovangigli [22, Chapter 2, Section 2.2]): where p is the pressure, S stands for the viscous stress tensor, E denotes the total energy per unit mass, q is the heat flux, F k denotes the diffusion flux and ω k is the production rate of the k-th species.Motivated by several recent studies devoted to the scale analysis as well as numerical experiments related to the proposed model (see Klein et al. [23]), our analysis is based on the following physically grounded assumptions: [A1] The species formation enthalpies h k are supposed to be constant.
[A2] The viscous stress tensor S is determined through Newton's rheological law where µ > 0, η ≥ 0 are respectively the shear and bulk viscosity coefficients.
[A4] The species diffusion fluxes are determined through Fick's empirical law where D k stands for the diffusion coefficient of the k-th species.
[A5] The model is consistent with the principle of mass conservation, that means, and the second law of thermodynamics requiring the entropy production to be non-negative.
The existence of global-in-time solutions for system (1.1 -1.4), supplemented with physically relevant constitutive relations, was established by Giovangigli [22,Chapter 9,Theorem 9.4.1] on condition that the initial data are sufficiently close to an equilibrium state.The same problem for a simplified model with large initial data was addressed in [4], [5] in the case of one space dimension and one irreversible chemical reaction.Large-data existence theory for the corresponding initial-boundary value problem posed on a spatial domain Ω ⊂ R 3 has been developed in [9], [18], still in the case when system (1.4) reduces to a single equation corresponding to one irreversible chemical reaction.
Yet in most practical applications one has to deal with dozens or rather hundreds of species undergoing many reactions that are, in general, completely reversible (see for instance Bose [3,Chapter 6]).Therefore, there is a fundamental need for developing a relevant existence theory.The main objective of the present paper is to undertake a first step in this direction.Given the complexity of the problem, we consider a simplified model based on hypotheses [A1 -A5] so that the essential lines of thought can be worked out straight-forwardly.
The main ingredients of our approach can be formulated as follows: • A suitable variational formulation of the underlying physical principles based on the second law of thermodynamics, in particular, replacing the energy balance (1.3) by the corresponding equation for the total entropy of the system.
• Physically grounded structural hypotheses imposed on the thermal equation of state for the pressure p.In particular, the effect of radiation, significant in the high temperature regime, is taken into account.
• A priori estimates based solely on boundedness of the initial energy and entropy of the system.As a matter of fact, this step requires the transport coefficients µ, κ, and D k to be effective functions of the absolute temperature.
• The weak stability property of the effective viscous pressure established by Lions [24], its generalization to non-constant viscosity coefficients proved in [14], combined with the approach based on the oscillation defect measures introduced in [19].
The main contribution to the existing theory, and the principal new difficulties to be dealt with can be characterized as follows: • The approximation scheme used to construct the solution is based, on one hand, on the Faedo-Galerkin type approximation to deal with the "fluid" part of the system while the "reaction" part requires uniform estimates based on the invariant regions technique (cf.Chueh et al. [6]).These two approaches being rather incompatible, some extra terms must be introduced in the approximate system.Moreover, in order to ensure strict positivity of the absolute temperature, a singular source term must be added to the approximate thermal energy equation as well as to the corresponding total energy balance, which makes the analysis quite delicate.
• In order to accommodate realistic growth conditions imposed on the transport coefficients, a new technique based on weighted oscillations defect measures must be used.Moreover, the velocity does not (is not known to) belong to the standard "energy space" W 1,2 (Ω; R 3 ) but rather to W 1,p (Ω; R 3 for a certain 1 < p < 2. This fact calls for rather delicate energy estimates obtained in Section 6.In particular, a generalized version of Korn's inequality is shown (see Proposition 6.1) that may be of independent interest.
• The standard entropy with the corresponding balance equation must be considerably modified (cf.Giovangigli [22, Chapter 2, Section 2.6]) in order to handle the reversibility of one or several chemical reactions.
The paper is organized as follows.In Section 2, we collect all the necessary hypotheses imposed on the non-linearities appearing in the constitutive relations.The weak formulation of problem (1.1 -1.4) supplemented with suitable boundary conditions is introduced in Section 3. The main existence result (see Theorem 4.1) is stated in Section 4. The rest of the paper is devoted to the proof of Theorem 4.1.A suitable approximation scheme based on the Faedo-Galerkin method applied to the momentum equation, and a regularization technique used for the remaining field equations, is introduced in Section 5. Section 6 is devoted to the limit passage in the family of approximate solutions.The most delicate part of the analysis concerned with weighted estimates of the oscillations defect measure is carried out in Section 7.

Constitutive relations
Unlike the field equations (1.1 -1.4) reflecting general physical principles, the constitutive relations characterize the specific material and close, at least formally, system (1.1 -1.4).In terms of the mathematical formalism, the constitutive equations can be expressed as typically non-linear functions relating the unspecified quantities in (1.1 -1.4) to the state variables , u, ϑ, and Y k , k = 1, ..., k.

Total energy, Gibb's equation
The total energy E can be written in the form where, in accordance with the basic principles of statistical mechanics (cf.Gallavotti [21]), the specific internal energy e is related to the pressure p and the species mass fractions Y k , k = 1, ..., n, through Gibb's equation where D denotes the differential with respect to the state variables , ϑ, Y k , and (see Giovangigli [22, Chapter 2, Section 2.6.1]).In formula (2.2), the symbol s denotes the specific entropy and, as a matter of fact, (2.2) may be viewed as its definition.Furthermore, in accordance with hypothesis [A1], the formation enthalpies h k as well as entropies s k , k = 1, ..., n, are taken to be constant.

Thermal equation of state
Motivated by Klein et al. [23], we consider the pressure p = p( , ϑ) independent of the species concentrations, more specifically, p satisfies the thermal equation of state where p R denotes the radiation pressure, which is particularly relevant in the high temperature regimes typical for combustion processes (see Bose [3,Chapter 11, Section 11.1]).Similarly, we set where e F = e F ( , ϑ), p F ( , ϑ) are interrelated through a general caloric equation of state characteristic for mixtures of mono-atomic gases.
Furthermore, we shall assume that p F = p M +p E , where p M is the classical molecular pressure obeying Boyle's law, while p E is the pressure of electron gas constituent behaving like a Fermi gas in the degenerate regime of high densities and/or low temperatures (see Chapters 1,15 in Eliezer et al [11]).Thus necessarily p F takes the form p F = ϑ where (2.9) Condition (2.9) reflects the fact that the specific heat at constant volume is strictly positive and uniformly bounded.The reader may consult [16] for more details and further discussion.

Transport coefficients
The well accepted physical stipulation asserts that viscosity of a gas is independent of the density (see Becker [1]).Here we suppose the viscosity coefficients µ and η are continuously differentiable functions of the temperature satisfying In agreement with Proposition 7.5.7 in Giovangigli [22,Chapter 7], the species diffusion coefficients coincide for all k, specifically, where, in accordance with the general axioms postulated in [22,Chapter 7], we assume that D is a continuously differentiable function of ϑ such that Finally, we take κ = κ F (ϑ) + κ R ( , ϑ), (2.16) where κ F , κ R are continuously differentiable functions satisfying Similarly to the above, the presence of the extra heat conductivity coefficient κ R is related to the effect of radiation (see Oxenius [26]).

Species production rate
The species production rates ω k are continuous functions of the absolute temperature ϑ and the species mass fractions Y 1 ,...,Y n .For the sake of simplicity, we shall assume that Furthermore, in agreement with postulate [A5], we suppose where the latter condition is enforced by the second law of thermodynamics (see Section 3 below).Finally, we suppose 3 Weak formulation

Entropy equation
For technical reasons, it seems more convenient to replace (1.3) by the total entropy balance.A straightforward manipulation yields By virtue of the second law of thermodynamics, the entropy production rate must be non-negative for any admissible process.In particular, the species production rates ω k have to satisfy (2.20).

Boundary conditions
The physical space considered in this paper is a bounded domain Ω ⊂ R 3 with the boundary of class C 2+ν , ν > 0. Generalization to domains with Lipschitz boundaries is possible via a suitable approximation procedure (see Poul [28]).Accordingly, system (1.1 -1.4) has to be supplemented with a suitable set of boundary conditions in order to obtain, at least formally, a mathematically wellposed problem.The concept of the weak solutions introduced below requires the energy flux trough the boundary to be zero.Therefore we suppose that where n stands for the outer normal vector, together with In addition, the impermeability condition (3.3) is supplemented either with the complete slip boundary conditions (Sn) × n| ∂Ω = 0, (3.5) where the latter, combined with (3.4), gives rise the standard no-slip boundary conditions u| ∂Ω = 0. (3.6)

Equation of continuity
Following DiPerna and Lions [8], we shall say that ≥ 0, u represent a renormalized solution of equation (1.1) on a time-space cylinder (0, T ) × Ω provided the integral identity Note that (3.7) already includes the satisfaction of the initial condition (0, •) = 0 in Ω as well as the impermeability condition (3.3).Relation (3.7) anticipates implicitly that all quantities are at least integrable on (0, T ) × Ω, in particular, the distributional div x u is representable by a function.

Balance of momentum
Similarly to the preceding section, the weak formulation of equation (1.2) reads All quantities appearing in (3.8) are supposed to be at least integrable, and S, p obey the constitutive relations (1.5), (2.7), respectively.In particular, the velocity field u must belong to a Sobolev space L p (0, T ; W 1,q (Ω; R 3 )), therefore it is legitimate to require u to satisfy the boundary conditions (3.3), or (3.6) as the case may be, in the sense of traces.

Total energy conservation
In agreement with the conservative boundary conditions introduced in Section 3.2, the total energy of the system is constant of motion.Specifically, where ϑ 0 , Y k,0 denote the initial distribution of ϑ, Y k , respectively.
The essential ingredient of this approach is replacing equation (3.1) by inequality (3.10), or, equivalently, we suppose that the "genuine" entropy production σ associated to a weak solution satisfies rather than (3.2).However, it is a routine matter to check that for any smooth weak solution relation (3.11) reduces to (3.2) and equation (1.3) holds (see Remark 2.3 in [18]).

Species mass conservation
Finally, we introduce a concept of "entropy" solutions for (1.4) requiring the integral identity to be satisfied for any test function ϕ ∈ D([0, T ) × Ω), together with where we have set 4 Main results
Having collected all the preliminary material, we are ready to state the main result of the present paper.
for a.a.x ∈ Ω.Furthermore, assume that the pressure p, the specific internal energy e, and the specific entropy s are functions of the state variables , ϑ, Y The remaining part of the paper is devoted to the proof of Theorem 4.1.In the rest of this section, we discuss briefly the meaning of instantaneous values of the state variables, in particular, the way the initial conditions are fulfilled.

Initial conditions
It follows from (3.7), (3.8), and (3.12) that the quantities , ( u), and Y k , k = 1, . . ., n, are at least weakly continuous with respect to the time variable t ∈ [0, T ].Here a vector valued function Z ∈ L ∞ (0, T ; L 1 (Ω)) is said to be weakly continuous provided the functions Consequently, the instantaneous values of , ( u), and Y k , k = 1, . . ., n, are welldefined and it makes sense to impose the initial conditions for these quantities.In particular, one can assign the initial value to u provided 0 is strictly positive.
On the other hand, in agreement with (3.10), the entropy density s admits only the one-sided limits Therefore the quantities ( s)(τ +), ( s)(τ −) have to be understood as Radon measures on Ω such that ( s Nonetheless, it is still possible to show that the entropy density ( s) attains its initial value at least in the weak sense, that means, for any weak solution in the sense of Definition 4.1.,where we have set as a direct consequence of the total entropy balance expressed through (3.10).Thus (4.3) follows as soon as we are able to show As the entropy production rate is always non-negative, one can use the Riesz representation theorem to rewrite (3.10) in the form After a straightforward manipulation, we deduce from (3.9), (4.6) that ) Here, in addition to (4.1), we have assumed that Clearly, ) and, furthermore, where we have used convexity of the function (y, z) → z 2 y .Thus (4.11) follows.Now, where we have set It follows from the thermodynamics stability hypotheses (2.8), (2.9) that (i) → H ϑ 0 ( , ϑ 0 ) is strictly convex for any fixed ϑ 0 , (ii) ϑ → H ϑ 0 ( , ϑ) − H ϑ 0 ( , ϑ 0 ) attains its global minimum for ϑ = ϑ 0 for any fixed .Consequently, in particular, (4.5) and, consequently, (4.3) hold.Note that we have shown that not only the entropy but also the absolute temperature attains its initial value in the sense of (4.14).

The approximate solutions 5.1 The approximation scheme
The basic strategy adopted in the proof of Theorem 4.1 consists in taking modified constitutive equations for p, e, as well as for the transport coefficients µ, κ, in order to obtain a system that can be solved by the methods developed in Chapter 7 in [13].More specifically, we take where Γ > 8.Moreover, the thermal energy balance (equation (5.10) below) will be augmented by a heat source proportional to δϑ −1/2 in order to ensure strict positivity of the absolute temperature.
The regularized Problem N δ , Problem D δ , with p ≈ p δ , e ≈ e δ , µ ≈ µ δ , κ ≈ κ δ , can be solved by means of an approximation scheme analogous to that developed in [13].More specifically, equation (1.1) is replaced by supplemented with the Neumann boundary condition and the initial condition (5.5) Note that, by virtue of the classical maximum principle, problem (5.3 -5.5) admits a priori estimates in the form exp(− (5.6) in particular, remains positive provided we can control the norm of div x u.
The momentum equation (1.2) is replaced by ) with the initial conditions where both u and u ε,0 satisfy the complete slip boundary conditions The meaning of the extra term ε∇ x u∇ x in (5.7) is that the kinetic energy balance holds even at the level of approximate solutions.In addition to (5.3), (5.7), we introduce the thermal energy balance equation: where ϑ satisfies the Neumann boundary condition together with the initial condition (5.12) Finally, the species balance equations (1.4) are taken in the form (5.15)

The first level approximations
We pursue the approach originated in Chapter 7 in [13], with the necessary modifications to accommodate the extra terms as well as the species balance equation.Our principal concern is to keep the total energy bounded at each level of the approximation procedure, and to meet the natural bounds Let X m = span{w j } m j=1 , be a finite-dimensional space endowed with the standard L 2 −Hilbert structure, where ) is a system of linearly independent functions.In addition, when dealing with Problem N we require the functions w j to satisfy the boundary conditions (3.3), and the set {w j } ∞ j=1 to be dense in Similarly, for Problem D, we take w j ∈ D(Ω; R 3 ), {w j } ∞ j=1 dense in C 1 0 (Ω; R 3 ).At the first level of approximations, equation (5.7) is replaced by a finite system of integral equations: to be satisfied for any w ∈ X m , where u ∈ C([0, T ]; X m ).Simultaneously, the density = [u] is determined in terms of u as the unique solution of problem (5.3 -5.5), while , n, are obtained from the parabolic system (5.10 -5.15).Note that, thanks to the additional (singular) source term in (5.10), the absolute temperature ϑ remains bounded below away from zero uniformly in [0, T ].Indeed equation (5.10) gives rise to where, in accordance with hypotheses (2.8 -2.9), Thus we can apply the classical maximum principle argument in order to conclude that ϑ(t, x) where Θ = Θ T, δ, min Note that, at the level of the Faedo-Galerkin approximations (5.16), boundedness of Θ reduces to showing a bound for u in the space L ∞ (0, T ; L 2 (Ω; R 3 )).
In view of (5.17), we are allowed to divide the thermal energy equation (5.10) on ϑ in order to obtain the approximate entropy balance Finally, the classical maximum principle applied to (5.13), together with hypothesis (2.21), yield (5.20) while (5.13) can be written in the conservative form (5.21) Local in time approximate solutions for fixed m can be found, exactly as in Chapter 7 in [13], by means of a fixed point argument.Their extension to the whole time interval (0, T ) as well as the possibility to take the limit for m → ∞, ε → 0 depends essentially on uniform estimates discussed in the next section.

Uniform estimates
The first rather and rather natural estimate expresses the principle of the total mass conservation.Integrating (5.3) over Ω and taking advantage of the boundary conditions (5.4), we get This is a universal bound independent of m, ε, and δ.The next step is to use the total energy balance established in (5.9), (5.10), namely Since both (5.9) and (5.10) are satisfied even at the level of the Faedo-Galerkin approximations, relation (5.23) can be used in order to deduce uniform estimates provided we can control the singular source term on the right-hand side.
To this end, we integrate the approximate entropy balance (5.18) over Ω and subtract the resulting expression from (5.23) in order to obtain Here, similarly to Section 4, we have set with θ = 1.Writing where we have taken while, by virtue of the thermodynamics stability hypotheses (2.8), (2.9), (see Lemma 2.1 in [17]).Thus relation (5.24) yields estimates independent of m and ε provided we show that the ε−dependent terms can be "absorbed" in the approximate entropy production rate.In order to see this, we compute .
Consequently, integrating by parts we get (5.27)By virtue of the thermodynamics stability hypothesis (2.8), the former term on the right-hand side of (5.27) is non-negative while the latter can be estimated as By virtue of our choice of κ δ specified in (5.2), we have whence the ε−dependent terms in (5.24) are dominated by the modified entropy production rate (5.28) Note that, by virtue of (2.3), (2.20), The uniform bounds provided by (5.22), (5.24) are sufficient for the approximate solutions to be extended on the whole time interval [0, T ] as well as for passing to the limit for m → ∞, ε → 0, respectively.As explained in Chapter 7 in [13], the only reason for splitting the approximation procedure into the ε and δ parts are the refined pressure estimates based on the multipliers ∇ x ∆ −1 [ ν ] (cf.Section 6.2 below).Specifically, we need ν = 1 when the artificial viscosity term ε∆ is present, while uniform estimates with respect to the parameter δ require ν very small.Otherwise, the technical parts of the limits processes m → ∞, ε → 0, δ → 0 are rather similar.For this reason we focus in what follows only on the last and most difficult step letting δ → 0.

The limit for δ → 0
In accordance with our previous agreement, we shall assume that we have already performed the limits m → ∞, ε → 0 in the approximate problems specified in Section 5. Accordingly, we focus on letting δ → 0. Thus we assume that we have a family of approximate solutions { δ , u δ , ϑ δ , Y 1,δ , . . ., Y n,δ } δ>0 solving Problem N or Problem D in the sense of Definition 4.1, where the pressure p δ as well as the transport coefficients µ δ , κ δ have been modified according to (5.1), (5.2).In addition, similarly to (5.23), the total energy balance contains an extra source term proportional to δϑ −1/2 , while the entropy production rate is augmented by δϑ −3/2 .
In accordance with our hypotheses imposed on the data, the expression on the right-hand side of (6.1) is bounded uniformly for δ → 0, the singular term being absorbed in the left-hand side via (5.28).Thus making use of (5.26) we conclude that ess sup ess sup where the constant c is independent of δ → 0. Furthermore, relation (3.13) yields Our next goal is to use estimates (6.2 -6.6) together with the structural hypotheses introduced in Section 2 in order to obtain uniform bounds in suitable function spaces.

Uniform bounds on the pressure
Our next goal is to establish uniform bounds on the pressure as well as of the internal energy in terms of a reflexive space L q , with q > 1.We start with the easier situation when the complete slip boundary conditions (3.3), (3.5) are imposed.Under these circumstances, the quantities where ∆ N stands for the standard Laplacean supplemented with the homogeneous Neumann boundary conditions, represent legitimate test functions in (3.8).Indeed, in accordance (3.7), we have for B, b as in (3.7),where the expression on the right-hand side may be viewed as a bounded linear form on the space L q (0, T ; W 1,p (Ω; R 3 ) for suitable p, q (cf.estimates (6.24), (6.25)).Thus, by virtue of the standard L p −theory for ∆ N , we get Furthermore, in view of (6.27), {∂ t ϕ δ } δ>0 bounded in L q (0, T ; L s (Ω; R 3 )), for certain q > 1, s > 5 (6.32) on condition that for a certain r >> 1 large enough.Consequently, taking δ B( δ ) = ν δ , with ν > 0 suitably small, we are allowed to use the quantities ϕ δ as test functions in the momentum equation (3.8) in order to conclude that uniformly for δ > 0 (the reader may consult Chapter 5 in [15], where the particular case In the case of the no-slip boundary conditions (3.6), the operator ∇ x ∆ −1 N has to be replaced by the so-called Bogovskii operator B (see Bogovskii [2]) -a suitable branch of solutions to the problem Relation (6.33) can be established via a rather lengthy procedure, where, however, the basic steps are the same as above (see [20] for details).

Positivity of the absolute temperature
Positivity of the absolute temperature is necessary for the entropy balance (3.10) to make sense.Recall that we have already established square integrability of ∇ x log(ϑ δ ) in (6.13).However, such a property does not guarantee positivity of ϑ δ itself.Our goal is to show that in particular, ϑ δ > 0 on a set of full measure.
In view of Proposition 6.1, it is enough to show that In order to see (6.35), we use hypothesis (2.9) in order to obtain provided we set S(1) = 0, where S is the function appearing in (2.10).Consequently, (6.35) follows from (6.3), (6.10), in particular, we get (6.34).

The limit in the field equations
The piece of information collected above is sufficient for all terms appearing in (3.1), (3.8), (3.9), (3.10), and (3.12) to belong to the reflexive space L p ((0, T ) × Ω), for a certain p > 1, in particular, all possible asymptotic concentrations are eliminated.Accordingly, we are allowed to pass to the limit at least in the weak topology of the space L 1 .As usual in the literature, the weak limits of compositions C(v δ ) are denoted as C(v).By virtue of (6.6), (6.7), (6.10), (6.14), (6.24), (6.25), we may assume that , and weakly in L 2 (0, T ; W 1,2 (Ω)), (6.40) passing to suitable subsequences as the case may be.
Finally, one can combine (6.14) with (6.44) to conclude that At this point, it is worth-noting that exactly the same procedure can be used to obtain a relation to be used below.
In view of (6.69), the proof of Theorem 4.1 is complete as soon as we show strong (pointwise) convergence of the densities δ → .This is the objective of the last two sections.
The next step is to convert (6.70) to a more familiar relation for any ϕ ∈ D((0, T ) × Ω).To this end, we need the following result in the spirit of Coifman and Meyer [7] (see Proposition 5.1 in [10]).
As already observed by Lions [24], Serre [29], the renormalized equation (3.7) represents a suitable tool for describing propagation of oscillations of { } δ>0 .Following [19] we introduce a family of cut-off functions ) We get provided δ = 1 Ω δ , and u δ are extended to be outside Ω via (6.73).Here, we have set Letting δ → 0 we deduce ∂ t L k ( ) + div x (L k ( )u) + T k ( ) div x u = 0 in D ((0, T ) × R 3 ).(6.77) At this stage, it is convenient to know if the renormalized equation (3.7) holds also for the limit quantities , u, specifically, ∂ t L k ( ) + div x (L k ( )u) + T k ( ) div x u = 0 in D ((0, T ) × R 3 ).(6.78)If this is the case, we obtain T k ( )div x u − T k ( ) div x u dx dt for any 0 ≤ τ 1 ≤ τ 2 ≤ T .Now the crucial observation reads which is a direct consequence of (6.71) as the pressure p is a non-decreasing function of the density (see hypothesis (2.8)).
Letting k → ∞ we conclude that • the validity of the renormalized equation (6.78) for the limit functions , u, • relation (6.80).
These two intimately related questions will be examined in detail in the last section.
A less obvious statement reads as follows (see Proposition 2.4 in [12]).
Then the limit functions , u solve (3.7) in D (Q).