NONLINEAR DIFFUSION EQUATIONS IN FLUID MIXTURES

. The whole set of balance equations for chemically-reacting ﬂuid mixtures is established. The diﬀusion ﬂux relative to the barycentric refer- ence is shown to satisfy a ﬁrst-order, non-linear diﬀerential equation. This in turn means that the diﬀusion ﬂux is given by a balance equation, not by a constitutive assumption at the outset. Next, by way of application, limiting properties of the diﬀerential equation are shown to provide Fick’s law and the Nernst-Planck equation. Moreover, known generalized forces of the literature prove to be obtained by appropriate constitutive assumptions on the stresses and the interaction forces. The entropy inequality is exploited by letting the constitutive functions of any constituent depend on temperature, mass den- sity and their gradients thus accounting for nonlocality eﬀects. Among the results, the generalization of the classical law of mass action is provided. The balance equation for the diﬀusion ﬂux makes the system of equations for diﬀusion hyperbolic, provided heat conduction and viscosity are disregarded. This is ascertained by the analysis of discontinuity waves of order 2 (acceleration waves). The wave speed is derived explicitly in the case of binary mixtures.

1. Introduction. Acoustics in mixtures is deeply influenced by diffusion. Though the subject is well-known, the literature shows different models of diffusion depending on the context at hand. Also in view of the interest of diffusion phenomena in physics, chemistry, and biology, it is then worth investigating the diffusion model equations and their influences on wave propagation.
The simplest and best known model of diffusion is, of course, Fick's law, a relation motivated by the obvious analogy between diffusion and heat conduction (see, e.g., [5]). If j is the diffusion flux and c is the concentration of the pertinent constituent then j = −D∇c (1) provides j as a result of the driving force, the gradient ∇c. The (approximate) balance equation ∂ t c = −∇ · j + ζ (2) and the assumption that the diffusivity D be a constant provide the classical parabolic differential equation ∂ t c = D∆c + ζ, which has the same form of the heat equation for the unknown temperature.
This scheme is much too simple, at least conceptually. The interaction among the constituents and the fact that diffusion describes a motion, relative to the barycentric reference, make the pertinent equations for the concentrations (as well as for the mass fractions) of more than two constituents be a coupled system of equations rather than decoupled parabolic equations. Despite the coarse model (1)-(2), eq. (1) is currently regarded as a reference constitutive equation for the diffusion flux j. It is a standard view that the pertinent diffusion flux is given by a constitutive equation where the gradient of the chemical potential is the essential term. In more general cases, for instance in more complex phenomena involving charged particles, the same view is applied by having recourse to the electrochemical potential (see e.g. [2,13]). In phase field theories, instead, the evolution of concentrations is mainly modelled by kinetic equations in terms of an appropriate free energy function (see, e.g., [24,8]). As a preliminary part, in this paper we re-examine the balance equations of a fluid mixture and to establish an appropriate set of equations governing the evolution of the diffusion fluxes. This in turn shows that the diffusion flux is governed by a balance equation of the rate type (evolution equation), rather than by a constitutive equation, thus eliminating the subjective character of the constitutive assumption.
The purpose of this paper is twofold. First, to prove that, despite the current orthodoxy in the literature of diffusion, the correct modelling of diffusion, in inviscid and non-conducting mixtures, provides a hyperbolic system of differential equations. For definiteness, discontinuity waves of order two are considered and the propagation speeds in two cases are derived. Hence we find how the wavespeeds are affected by the constitutive properties and the concentrations of the constituents. Second, to clarify that (seemingly) conceptually different models appearing in the literature (like, e.g., the Maxwell-Stefan diffusion) can be framed within the general scheme so established. Also, the purpose is to indicate how Fick's law (1) for the diffusion flux, or analogous equations, such as the Nernst-Planck equation, for more complex phenomena, may be obtained in particular limiting conditions.

Notation.
Let Ω be a time-dependent region occupied by a mixture of n constituents and denote by x the position vector of a point relative to a chosen origin. The functions describing the evolution of the mixture have Ω × R as their common space-time domain. We denote by ∂ t the derivative with respect to t ∈ R and by ∇ the gradient operator in Ω. The subscripts α, β = 1, 2, ..., n label the quantities pertaining to the constituents of the mixture. Hence ρ α is the mass density, v α the velocity, ε α the internal energy density, θ α the absolute temperature, η α the entropy density, on Ω × R, of the αth constituent. Throughout α and β are shortands for sums over (the constituents) 1, 2, ..., n.

2.
Balance equations for single constituents. The mass density ρ α and the velocity v α are required to satisfy the continuity equations where τ α is the mass produced, per unit time and unit volume, of the αth constituent. The overall conservation of mass requires that The mass density of the mixture, ρ, and the barycentric velocity, v, are defined by Hence summation of (3) over α provides Let a superposed dot stand for the derivative along the barycentric motion,ψ = ∂ t ψ + v · ∇ψ for any function on Ω × R. Hence we can write (6) in the forṁ The mass fraction ω α and the diffusion velocity u α are defined by Since α ω α = 1 we may regard ω 1 , ..., ω n−1 and ρ as n variables in 1-1 correspondence with ρ 1 , ρ 2 , ..., ρ n . Substitution of ρ α with ρω α in (3) and use of (7) provides It is worth remarking that in continuum mechanics the mass fraction ρ α /ρ is called the concentration and denoted by c α . In physics and chemistry, more often c α is the molar concentration, M α being the molar mass of the αth constituent. For ease of comparison with the main literature on diffusion, we adopt the notation (10) and hence hereafter c α is the molar concentration. Still within the Eulerian description, namely in dealing with functions on Ω × R, we write the equations of motion in the form where T α is the stress tensor, f α the body force, and m α the growth of linear momentum, due to the action of the other constituents. Also, we write the balance of energy in the form where q α is the heat flux, r α is the external heat supply, and e α is the growth of energy. Due to the internal character, the growths {m α } and {e α } are subject to The balance of entropy is written in the form where φ α is the entropy flux and s α the entropy supply. Also, σ α is the entropy growth. The second law of thermodynamics is expressed by saying that has to hold for all admissible fields compatible with the balance equations (3), (11), [26,3].

ANGELO MORRO
The growths {m α }, {e α }, {σ α } are non-objective in that they are related to the particular reference chosen for the Eulerian description. Objective quantities are provided in a while within the Lagrangian description.
For any function φ α on Ω × R we denote byφ α the αth peculiar time derivative that is the derivative following the motion of the αth constituent, By means of (3) we have As a consequence, eqs. (3) and (11) becomè Hence, applying again (16) we find that whereε By means of (4), the constraints (13) provide Equations (22) involve only objective vectors and scalars.
3. Balance equations for the whole mixture. The third principle listed by Truesdell [26] states that the motion of the mixture is governed by the same equations as is a single body. Hence, the balance of mass, linear momentum, and energy are required to take the formρ In addition, the balance of moment of momentum results in the symmetry of T. It is known (see, e.g., [26,3,17]) that eqs. (23)-(25) hold with ρ and v given by (5) whereas and The subscript I denotes the inner quantities, namely the sum, or the average, of the corresponding peculiar terms, Equation (24) can be written in the form which shows the nonlinear effect of diffusion velocities {u α } on the motion of the barycentric reference. Now,ε We can expressu α asu As a consequence, owing to (25) and (30) we obtain Some rearrangements provide the identities Consequently we find that 4. Evolution equations for the diffusion fluxes. We now determine the evolution equation for the fluxes {h α } (see [17,15]). By (17) and (18) we have Since Substitution of ρv from (24) giveṡ Equation (33) governs the evolution of h α relative to the barycentric reference frame. It is rather cumbersome and nonlinear. Yet it is exact, i.e., free from any approximations or constitutive assumptions. The occurrence of (33) is quite exceptional in the literature; it was derived first in [17] and re-visited in [19,15]. It is of interest to consider some "generalized equations of motion" mainly developed as approximations to the Boltzmann equation [14,23]. They show features that look intermediate between the equation of motion (18) and the balance equation (33).

Generalized equation of motion; Maxwell-Stefan diffusion.
In the present notation, eq. (44) of [12] reads (34) The stress T T T α is said to represent the shearing force on α and then we can identify −p α 1 + T T T α with T α of (18). It seems that chemical reactions are disregarded. Hence the remaining terms stand for the growth of linear momentum, The restriction (13) for {m α } holds provided only that The quantities D αβ are usually referred to as Maxwell-Stefan diffusion coefficients while x α is the αth mole fraction that is related to ω α by M being the average molar mass. Hence m α is taken as due to the temperature gradient ∇θ and the velocity differences v β − v α = u β − u α .
Other models take ρ α f α and m α so that occurs as the force on the αth constituent, in addition to the term involving the temperature gradient (see (92) of [12]). Here c α is the αth molar concentration, c t the total molar concentration. The motivation for the force (36) seems to be the assumption that any deviation from equilibrium leads to a diffusion flux and that, at equilibrium, the chemical potential gradient ∂ xα µ α ∇x α , regarded as driving force, is counterbalanced by molecular friction that is proportional to the velocity differences [21]. It is apparent that eq. (33) requires only that we select the growth m α whereas no subjective choice is required about the so called driving force.

5.
Exploitation of the entropy inequality. Look now at the entropy inequality (15). It is convenient to express φ α in the form k α being called the extra-entropy flux, while Hence, owing to (16), inequality (15) becomes Upon replacing ∇ · q α − ρ α r α from (20) we have By means of the free energy we can write the entropy inequality in the form The dependence on ∇ρ α is motivated also by applications to the modelling of phase transitions. Moreover, we let the interaction termsε α andp α , of the αth constituent, depend also on quantities pertaining to the other constituents. In addition, we assume that η α , τ α , T α , and q α are continuous functions while ψ α and k α are differentiable. Hence (39) becomes As a consequence, (39) reduces to For any function g α on Ω × R the identitỳ holds. Hence some rearrangements give and replaceρ α with −ρ α ∇ · v α + τ α . The inequality can then be written in the form This relation suggests that we take the extra-entropy fluxes k α as In addition replace L α with D α + W α , D α ∈ Sym being the stretching and W α ∈ Skw the spin. The arbitrariness and linearity of W α implies that T α + ρ α ∇ρ α ⊗ ∂ ∇ρα ψ α ∈ Sym, α = 1, 2, ..., n.
(43) If ψ α is independent of ∇ρ α then p α = ρ 2 α ∂ ρα ψ α , which is the classical form of the (partial) pressure. We can then view (43) as the general form of the pressure when the dependence on ∇ρ α is allowed.
Inequality (42) becomes The viscous and heat-conducting constituents are modelled by letting for any constituent α. The compatibility of wave motion with the requirements of the entropy inequality has been investigated for a number of dissipative models as, e.g., those in [10,16,11]; see also [25]. Here waves are investigated in the particular case of inviscid and non-conducting constituents. Inequality (42) reduces to Again we allow for separate non-negative contributions to the inequality namely If we regard µ α = ψ α + p α ρ α as the αth chemical potential then inequality (45) may be viewed as the condition governing the reactions between constituents at different temperatures.
Inequality (44) is a requirement on the admissible e α , m α , τ α . Often models are considered for diffusion in non-reacting mixtures with a single temperature, θ α = θ. In such cases (44) reduces to which is the thermodynamic requirement on the momentum growths.
The symmetry of M and the skew-symmetry of N implies that α m α = 0.
To check the compatibility with the second law we observe that, if the constituents have a single temperature and chemical reactions do not occur, then inequality (42) becomes If, as is the case, m α and q α are independent of D α then the inequality requires that Now, We let M αβ ≥ 0, α, β = 1, 2, ..., n. Inequality (47) is then satisfied provided with ξ α possibly being a tensor. It then follows that a model of the Maxwell-Stephan type for the linear momentum growths is compatible with thermodynamics provided we take the heat fluxes in the form (48).

5.2.
A representation of m α in reacting mixtures. An expression for m α is now determined in the case where each constituent has the same temperature. By To establish a thermodynamically-consistent model we let Hence (13) requires that α w α = 0.
Moreover, substitution of m α in the left-hand side of (49) and v α = v + u α gives As a consequence we have A simple model arises if we restrict attention to binary mixtures. Take (50) as an equality so that This implies that 6. Discontinuity waves. We now investigate some propagation properties of mixtures undergoing diffusion. This investigation is of interest in many respects. Here, diffusion proves to be modelled by a hyperbolic system of equations. The intrinsic nonlinearities play an essential role in the propagation modes. The evolution of the constituents of the mixture may be described by the 2n fields {ρ α }, {v α } on Ω × R. Diffusion may be described more conveniently by the fields ρ, {ω α }, v, {h α }, obviously subject to α ω α = 1, α h α = 0. The pertinent differential equations arė The system is closed if we let τ α and m α be functions of ρ, {ω α }, {h α }. It is natural though to let τ α depend also on the temperature θ. In this sense we might close the system by appending eq. (31) with ε I as a function also of θ. Yet, for definiteness, we regard the temperature as a parameter. This is allowed by letting and making the assumption that the entropy η α is constant (adiabatic approximation). Of course the constancy of {η α } is reasonable if heat conduction and viscosity are negligible. This assumption in turn is consistent with the next analysis of discontinuity waves which would be ruled out by the Navier-Stokes-Fourier model. We now look at singular surfaces of order 2 [27,4] and let ρ, {ω α }, v, {h α } be continuous everywhere whereas their time and space derivatives suffer jump discontinuities at a surface σ. Though there are differences among the velocities of the constituents, each velocity v α , and hence v, is continuous in space and time.
For simplicity we restrict attention to plane discontinuity surfaces. Hence σ ⊂ Ω is a plane moving in the direction of the unit vector n. Let Ω − , Ω + be the subregions behind and ahead of σ so that Ω = Ω − ∪ σ ∪ Ω + and n is directed from Ω − to Ω + . Moreover, by we let c α be continuous so that By means of (54) we obtain the following nonzero jumps, Hence, by (26) we find that We now compute the jump relations associated with (52)-(53). In view of (54)-(58) we obtain The occurrence of {u β } in (61) and (62) imply that the propagation modes are affected by the velocity fields ahead of the front σ. For definiteness and simplicity we make the assumption that h α (t, x) = 0, x ∈ Ω + , which means that the mixture is not diffusing ahead of the wave. Equations (59) and (60) remain unchanged whereas (61) and (62) simplify to where κ α := ∂ ρα p α .
Hence two speeds U + , U − may occur, As expected, if κ 1 = κ 2 = κ then U 2 + = U 2 − = κ. It is worth mentioning that a similar analysis of harmonic waves, in binary mixtures, is developed in [18,20] for non-reacting constituents subject to a non-uniform common temperature. 7. Limiting properties of the diffusion fluxes. The fluxes {h α } are determined by the evolution equations (33) which, in view of (24), can be written in the forṁ By (26), T involves all of the fluxes {h α } in the nonlinear form α h α ⊗ h α /ρ α . Hence (70) is a system of nonlinear first-order differential equations for {h α }.
The mixture is allowed to be chemically reacting and hence m α − τ α v is given by (56). Restrict attention to binary mixtures so that (51) holds. By (10), Since v β − v α = u β − u α we then have Additional assumptions are now made. First, diffusion occurs in stationary conditions,ḣ α = 0. Second, the mixture is so dilute that we can regard the solute as diffusing in a steady solvent and then L = 0. Third, the nonlinear terms h α ⊗ u α are disregarded. Equation (70) then provides ζ h 1 = ∇ · T 1 − ω 1 ∇ · (T 1 + T 2 ) + ρ 1 (f 1 − f ).