INFLUENCE OF A NONLINEAR DEGENERATE DIFFUSION ON AN ADVECTION-DIFFUSION EQUATION IN A DIFFUSE INTERFACE FRAMEWORK ∗

. This work is motivated by the modelling a liquid-vapour flows with phase transition describing the evolution of the coolant within an heat exchanger (e.g. the core of a Pressurized Water Reactor). We investigate an advection-diffusion equation with a degenerate and nonlinear diffusion coefficient. The degeneracy corresponds to a liquid-vapor mixture in the original model whereas the diffusion coefficient is non-degenerate in the pure phase cases. We focus on the influence of the diffusion coefficient on a simple 1D configuration for which some analytical computations can be done, leading to a surprising behavior of the phase transition with respect to the diffusion.


Context
A heat exchanger is a device that facilitates the process of heat exchange between two fluids that are at different temperatures. Heat exchangers are used in many engineering applications, such as refrigeration, heating and air-conditioning systems, power plants, chemical processing systems, automobile radiators, and waste heat recovery units. Nuclear core of a pressurized water reactor is an example of heat exchangers.
In some applications, such as the flow in nuclear reactor cores, convection is characterized by a low Mach number, where the convective velocities are much slower than the speed of sound in the fluid. This has promoted the development of so-called low Mach number models, which filter the sound waves. In the context of pressurized water reactor cores, an asymptotic low Mach number model, called lmnc (for Low Mach Nuclear Core), has been derived and investigated in a series of papers [1,2,4,5,[7][8][9]. The model was derived through an asymptotic expansion performed in the compressible Navier-Stokes equations with an energy source term. It consists of a transport equation upon a thermodynamic variable (here the total enthalpy), of a divergence constraint upon the velocity (with a nonlinear coupling source term which underlines the dilatation property of the flow) and of the momentum equation. The fluid is described by a single equation of state taking into account the phase transition by supposing that, when both vapour and liquid phases are present, they have the same pressure, temperature and chemical potential. The equation of state is piecewise defined with respect to the enthalpy identifying a liquid phase, a gas phase and a mixture of them. The thermal diffusion degenerates in the mixture region.
In the present paper we are interested in studying the influence of the thermal diffusion.

A heat exchanger model
For some bounded domain Ω ⊂ R 3 , the Navier-Stokes-Fourier system with source terms in a 3D conservative formulation reads      ∂ t ρ + ∇ · (ρu) = 0, ∂ t (ρu) + ∇ · (ρu ⊗ u) + ∇p = ρg + ∇ · σ(u), ∂ t (ρh) + ∇ · (ρhu) = Φ + ∇ · (ω∇T ) + σ(u) : ∇u + ∂ t p + u · ∇p. (1) The unknows are the velocity field (t, x) → u(t, x), the total enthalpy (t, x) → h(t, x) and the pressure (t, x) → p(t, x). The specific density ρ and the temperature T are related to the enthalpy h and the pressure p through an equation of state. The heat conductivity ω, assumed to be constant and isotropic for each phase, is also related to the enthalpy h and the pressure p. The power density (t, x) → Φ(t, x) ≥ 0 is a given function of time and space modelling the heating of the coolant fluid (e.g. due to the fission reactions in the nuclear core). Finally, g is the gravity field and σ(u) models viscous effects.
In a low Mach number regime (i.e. when the speed of fluid is much smaller than the speed of sound), an asymptotic expansions with respect to the Mach number can be applied to the compressible model, and this leads to a simpler model where the acoustic waves are filtered out. The "original" pressure p is decomposed in two pressure fields p(t, x) = p * +p(t, x), where p * is the reference state pressure (that is an average pressure constant in time and space) andp is the perturbational pressure (often called "dynamic pressure"). For low Mach number flows, an asymptotic analysis shows thatp(t, We can thus approximate p(t, x) by p * in the computation of thermodynamics quantities and we obtain the so-called LMNC model [2,4,8], that can be written as (2) In this system the unknows are the velocity field u, the enthalpy h and the perturbational pressure (t, x) → p(t, x). The specific density ρ and the temperature T are linked, by an equation of state, to the enthalpy h and to the constant given thermodynamic pressure p * .
In the following the power density Φ is chosen constant and positive. We also neglect the viscous terms.

Diffusion term with phase transition
The fluid can be in liquid phase ℓ, vapour phase g or a mixture of them m.   Figure 2) Since the temperature is a function of the enthalpy h and the (thermodynamic) pressure p, the gradient of the temperature can be written as In the lmnc model, the thermodynamic quantities are evaluated at p = p * , which is constant. Thus, from now letting drop the dependency upon p * , the diffusion term can be written as where λ κ def = ω κ /c p,κ for κ = ℓ or g with ω κ the heat conductivity (assumed to be constant and isotropic for each pure phase κ) and c p,κ def = ∂h ∂T p the isobar heat capacity of the pure phase κ. We note that if the open domain of the mixture is non-empty, a zero diffusion coefficient can be set in this region.
Let us denote By defining λ m = 0 and by virtue of (2), with (3), the obtained equation on enthalpy comes down to It differs strongly on the enthalpy equation written in [13] where the transition liquid-gas is forced at a sharp interface by the addition of a singular source term and a non-degeneracy of the thermal diffusion. The generality of the model proposed here is that the liquid-gas transition can also take place through a sharp interface or through a mixing region (in a diffuse interface framework). In the following, we will see that both can occur and when the liquid-gas transition is a sharp interface, our enthalpy solution coincides with that of the Stefan models written in temperature [12] or enthalpy [13]. This is the main point of this paper and an analytic solution, derived in the 1D steady-state case in the following section, shows this richness of the phase transition.

Steady-state model
In this section, we focus on the steady-state solution, i.e. we seek for h, u andp solution in R 3 + of In the unidimensional case, the first and last equations can be written as with ρ a function of h. To model the configuration described in the Figure 1, the system is closed by the Dirichlet boundary condition corresponding to an inlet condition: the inlet enthalpy h(0) = h e < h s ℓ and the inlet flow rate (ρv)(0) = D e . Then the flow rate is constant on the domain: (ρv)(y) = D e for all y ∈ R + .
In order to focus on a physically relevant solution, we assume that lim y→∞ ∂ y (λ(h)∂ y h) = 0 so that, thanks to (4), We are thus interested in y → h(y) solution of the following non-linear strongly degenerate convectiondiffusion problem: To better understand the role of each parameter, in the following we will consider three cases for the three constant values λ ℓ , λ m and λ g defined in (3): In the last case, which corresponds to the phase transition model, two cases have to be considered. Indeed, it will be seen that the mixture region can disappear for a sufficiently large gas diffusion coefficient λ g .

Steady solution without diffusion
Proposition 2.1 (Steady solution without diffusion). When λ ℓ = λ m = λ g = 0, the steady enthalpy solution of (6)- (7) is Let us denote x s ℓ (resp. x s g ) the liquid/mixture transition (resp. mixture/gas) without diffusion (see Figure 3). The fluid is in pure liquid phase for y ≤ x s ℓ , in pure vapour phase for y ≥ x s g and a mixture at saturation if

Steady solution with a strictly positive diffusion
Without degeneracy of diffusion, the unique solution can be explicited as follows. (3), is piecewise positive constant for all h, the steady enthalpy solution of (6)- (7) is

Proposition 2.2 (Steady solution with a positive diffusion). If λ(h), defined in
where the constants C 1 , C 2 , C 3 , C 4 , depending on y s ℓ and y s g , are for κ = ℓ or g and y s ℓ and y s g satisfies the continuity flux Proof. The non-degenerate diffusion, i.e. λ(h) > 0, leads to smooth solutions in the sense that discontinuities can not occur since order one distributions generated by the diffusive term and discontinuous solutions can not be compensated by order zero distribution. Furthermore, for continuous solutions, the transport term is then a diffuse measure leading to a continuous diffusion flux. Then, the continuous solution induces a mixture region and satisfies the flux continuity (11). Due to the growth of the steady state solution (Φ and D e are strictly positive), three regions, liquid, mixture, gas, are juxtaposed from low to high y values. On each region, an explicit solution solves the equation (6), (12) where the unknows C 1 , C 2 , C 3 , C 4 , C 5 , C 6 , y s ℓ , y s g have to satisfy • the enthalpy continuity conditions at inlet and transition points: • the flux continuity conditions at transition points: • the asymptotic behaviour The conditions (13), (14), (15) and (17) allow to define C 1 to C 6 as functions of the transition points y s ℓ and y s g . These points are thus implicitly defined by (16) and finally we obtain the solution. We note that the transparent boundary condition (17) implies C 6 = 0 so that y → h g (y) is linear and compatible with the semi-infinite domain. Moreover, although y → h ′ g (y) does not involve λ g , the flux continuity (16) does. □ The Figure 4 shows the graph of y → h(y) with λ m > 0, that is a non-degenerate artificial diffusion in the mixture region, leading to a continuous solution. Depending on the increasing ratio λ g /D e when other parameters are fixed, the mixture region reduces and is pushed on the left. A stiffer solution appears at the gas transition. On the last plot (Figure 4c), the mixing region is very sharp and is not included in the large mixing region occurring without diffusion in the whole phases.
The first plot of the Figure 5 is the one of the Figure 4a. In the next figures is plotted a reduced diffusion coefficient λ m . It is then exhibited the convergence of the enthalpy as λ m goes to zero. A discontinuous solution is suggested at the limit. It occurs at the transition to the gas phase. For such diffusions in the liquid and gas phases, the mixture region is slightly affected by the decreasing diffusion λ m and remains at the limit.
We will see in the following section that, when the diffusion is allowed to degenerate, i.e. the diffusive coefficient λ(h) vanishes for values of h in the mixture region, solutions are not necessarily smooth and discontinuous weak solutions must be sought.

Steady solution with degenerate diffusion
The degenerate diffusion can lead to a jump in the diffusion flux at phase transition. For example, as the enthalpy flux is zero in the mixture region, a flux diffusion jump occurs at the interface with the mixture if the space derivative of the enthalpy is not zero at this transition in the liquid (or gas) phase. A Dirac measure follows and has to be compensated by an other Dirac measure in the transport term. A discontinuous solution is then necessary. Due to the growth of the steady state solution (Φ and D e are positive) and to the presence of the liquid phase at inlet, a liquid phase is present for low y values, then follows a transition which can be liquid-mixture (followed by the transition mixture-gas) or liquid-gas. The two possible cases are presented below.
the mixture zone is present and the solution of (6)- (7) satisfies where Proof. As for the nondegenrate case, on each phase or mixture, the integration of the solution gives The five constants C 1 · · · C 5 have to be determine and a particular choice will be made to define a unique viscosity solution. Searching a solution of (6) with mixture, the transition liquid-mixture induces a free dissipative flux in the mixture. A positive jump on dissipative flux can not be compensated by the enthalpy jump on the transport term since the enthalpy is increasing. As a matter of fact, if we denote h the enthalpy jump at the liquid-mixture interface, the jump relation of (6) leads to Since y → h(y) has to be increasing, we should have h ≤ 0. Hence, in the present case, the only possibility is to have both h = 0 and h ′ (y s,− ℓ ) = 0. The inlet Dirichlet boundary condition and the Dirichlet condition h(y s ℓ ) = h s ℓ define the constants C 1 and C 2 parametrized with respect to y s ℓ . Then, the additional Neumann boundary condition at y s ℓ defines implicitly the unique value of y s ℓ , C 1 and C 2 , that is the unique solution in the known liquid phase.
The expressions of the solution in the mixture region follow by the continuity of the enthalpy at the interface with liquid, thus C 3 = h s ℓ . Furthermore, since C 5 = 0 by virtue of the asymptotic behavior (7), The difficulty is then to define the transition point y s g and the constant C 4 . The jump relation at this point is obtained in the same way as for the point y s ℓ , but this time, the two jumps can compensate each other: which has to be lower than the maximal gap h s g − h s ℓ . There is an infinity of solutions to such a problem in the sense of distribution. The transition y s g separating enthalpy lower than h s g to upper can cover an interval, that is, Thus, In order to define an unique solution, we choose a solution corresponding to the smallest mixture region, that is with the smallest value of enthalpy h g (y s g ) = C 4 = h s g in the gas region at the transition. The reader can convince himself (see Figure 5) that this choice corresponds to the viscosity solution by studying the limit of the viscous solution obtained in Proposition 2.2. This implies that and the position y s g follows □ Remark 2.4. Gas diffusion reduces the mixture region for steady solution since where x s κ are the transition points without diffusion defined in (9). We remark also that y s ℓ > x s ℓ if λ ℓ > 0. We point out that if (and only if) the condition is not fulfilled, a transition liquid-gas has to be written (cf. the following section).
The Figure 6 shows some steady solutions when (18) is fulfilled so that the mixture zone is present. The mixture region is always smaller than the one obtained without diffusion, i.e. y s g − y s ℓ ≤ x s g − x s ℓ . More precisely, we observe: • The Figure 6a shows the steady solution when Φ = D e , λ ℓ = 2D e and λ g = 0. We can see the effect of the diffusion only in the liquid phase: the liquid/mixture transition occurs at y s ℓ > x s ℓ and the slope of h ℓ at y s ℓ is zero. • Figure 6b shows the steady solution when Φ = D e , λ ℓ = 0 and λ g = D e . We can see the effect of the diffusion only in the vapour phase: the mixture/vapour transition occurs at y s g < x s g ; the width of the mixture zone is reduced w.r.t. the case without diffusion and there is a jump on the enthalpy h at transition y s g . • Finally Figure 6c shows the steady solution when Φ = D e , 2λ ℓ = D e and λ g = D e . We can see the effect of the diffusion both in the liquid and in the vapour phases: the liquid/mixture transition occurs at y s ℓ > x s ℓ and the slope of h ℓ at y s ℓ is zero; the mixture/vapour transition occurs at y s g < x s g ; the width of the mixture zone is reduced w.r.t. the case without diffusion as observed in the remark 2.4 and there is a jump on the enthalpy h at transition y s g .

Case II: Liquid/gas (no mixture, Stefan problem)
The solution without mixture occurs for a sufficiently large ratio  the mixture zone is not present and the solution of (6)- (7) satisfies where the position y s separating the liquid phase to the gas phase is implicitly defined by h ℓ (y s ) = h s ℓ . Proof. As shown in the proof of Proposition 2.3, there isn't viscosity solutions to (6)-(7) with mixture if On each phase, the integration of the differential equation gives The four constants C 1 . . . C 4 and the transition point y s have to be determined. The inlet Dirichlet boundary condition and the asymptotic condition (7) define the constants C 1 = h e and C 4 = 0.
A particular choice will be made to define an unique viscosity solution. Such a viscosity solution has to satisfy h(y s,− ) = h s ℓ , h(y s,+ ) = h s g . This implies that C 3 = h s g and

The jump relation at the transition point y s reads
which is positive under the assumption (21) and is then compatible with the increasing enthalpy h.
Computing the derivative of the enthalpy in the liquid phase, we obtain an implicit definition of the position y s : We now prove that this value is unique. Let us we denote t def = exp −∞, and f ′′ (t) = −R/t < 0, the function f admits an unique zero for t > 1. □ Figure 7 shows some steady solutions when ℓ so that the mixture zone is not present. • Figure 7a shows the steady solution when Φ = D e , λ ℓ = D e and λ g = 2D e . This corresponds to the smallest choice of λ g such that the jump is equals to h s g − h s ℓ . The mixture is not present, the liquid/vapour transition occurs at y s and the slope of h ℓ at y s ℓ is zero. • Figure 7b shows the steady solution when Φ = D e , λ ℓ = D e and λ g = 4D e . The jump is equals to h s g − h s ℓ , the mixture is not present, the liquid/vapour transition occurs at y s and the slope of h ℓ at y s ℓ is strictly positive. When λ ℓ → 0 the slope tends to +∞: +∞ • Finally Figure 7c shows the steady solution when 2Φ = D e , λ ℓ = 0 and λ g = 5D e . The mixture is not present, the liquid/vapour transition occurs at y s . Since λ ℓ = 0 there is a jump greater to h s g − h s ℓ so that we have a jump in liquid phase that corresponds to a jump on temperature and h ℓ (y s ) < h s ℓ :

Generalization to negative power densities
An enthalpy jump was exhibited at the gas transition point y s g . The gas phase does not play a different role from the liquid one: it is the situation of a heating a fluid (that is Φ > 0) that induces an increasing solution in the direction of flow with a jump localized at the last transition point.
In the case of a cooling heat exchanger (that is Φ < 0), the fluid is injected in gas phase (h e > h s g ): the stationary solution is thus decreasing in the direction of flow and it presents: • a jump at the mixture-liquid transition point in case of • a discontinuity at the gas-liquid transition point otherwise.

Link to the Stefan problems
The solution exhibited in the Proposition 2.5 is exactly the stationary solution of the following stationary Stefan problem, a free boundary problem. That is, find the position y s so that the over-determinate elliptic problem is satisfied: Written with the temperature as unknown, we obtain the classic stationary Stefan problem, with a continuous temperature at the transition and a jump condition on thermal flux, as it can be found in [12], Note that the last relation of (24) describes the jump on heat flux due to latent heat. This problem can be solved even if (18) holds, but the solution would loose a physical meaning due to inconsistent temperature with respect to the considered phase. As a matter of fact, the last relation of (24) gives (18), which induces a temperature above T s in a neighborhood of y s in the liquid phase.
To summarize, the Stefan formulations (23) or (24) have to be used only under (21) contrary to (6)-(7) which is valid under (21) and (18). Under (21), the solutions of the formulations are the same; under (18), the Stefan solution of (23) or (24) isn't physical. That is why the PDE formulation (6)-(7) has to be preferred for numerical solving as it is shown in Appendix A, extended to the non-stationary case.

Conclusion and perspectives
It has been shown the expression of the physical stationary solution for a one dimensional degenerate nonlinear elliptic problem modelling enthalpy subjected to a flow and thermal diffusion with phase change. The expression of the solution depends strongly on the ratio of heat deposit with flow rate to compare with the ratio of flow rate and diffusion coefficient in the gas phase (see (21)). It is notable that discontinuous solutions appear in the context of a modelization with diffuse interface. Furthermore, this mixture region which always exists without thermal diffusion, can disappear depending on the ratio (21). The selected solution is the viscosity solution leading to a minimal mixture region. The higher the gas diffusion coefficient, the smaller the mixing region, until it disappears.
Some unsteady simulations are proposed in the appendix but a more thorough numerical study is to follow for the full model proposed (2) in a future work.

A. Numerical Time-Dependent solution of the constant-density model
The time-dependent 1D model (4) can be written in a non-conservative formulation as In the stationary case, the model does not depend on the expression of the density with respect to the enthalpy. It is remarkable that whatever the density expression, the time-asymptotic solution depends only on the uniform flow rate. In the physical model, the density is assumed to be a continuous function of the enthalpy h, piecewise defined for each phase or mixture. In the transient regime, the law describing the density with respect to the enthalpy affects the unsteady solution. With the EoS describing a mixture at saturation as well as the pure phases, the density is a continuous and piecewise defined function of the enthalpy (and the constant thermodynamic pressure p * ). For simplicity, we choose here to exhibit the transient regime in the case of an artificially uniform in space and constant in time density. Let us consider ρ(h) = r > 0 constant for all h. The system becomes with the boundary condition (rv)(0) = D e > 0 and h(0) = h e < h s ℓ . Thus v(y) = v(0) for all y ∈ R + and (ρv)(y) = rv(0) = D e for all y ∈ R + . The enthalpy equation can be written as The unsteady model, written on the domain (t, x) ∈]0, T [×]0, +∞[, reads , associated to the following boundary and initial conditions, We define δt the time step and h n ≈ h(nδt). It is then essential to numerically solve this evolution problem by a fully implicit scheme, associated to a gradient scheme [6], in order to produce a viscosity solution. All the simulations respects the enthalpy properties at transition (jump relations), namely a jump of enthalpy (when it has to be) at the transition to the gas phase with an enthalpy close to h s g (up to the order of space discretization step) at the first node where enthalpy is upper than h s g . The increasing solution is also observed on the discrete solution avoiding spurious numerical oscillations.
Note that the solution converges, when time goes to infinity, to the analytic solution explicited in the Propositions 2.2 and 2.3. Note also that mixture appears in transient regime even if the asymptotic solution solves the Stefan problem with a transition liquid/gas, as on the Figure 9.
On Figures 8 and 9, the initial condition is in liquid phase, then appears the mixture and finally the gas phase. The parameters of Figure 8 satisfy condition (18) with mixture on time-asymptotic solution.
• Figure 8a shows the solution without diffusion in the gas phase and the solution shows the characteristic transition liquid-mixture with a zero space derivative in the liquid at the motionless transition. • Figure 8b shows the solution without diffusion in the liquid phase and the solution shows the characteristic transition mixture-gas with a growing jump of enthalpy at the moving transition. • Finally Figure 8c shows the solution with diffusion in the liquid and gas phase gathering the two previous behaviors. The parameters of Figure 9 satisfy condition (21) without mixture on the time-asymptotic solution.
• Figure 9a shows the apparition of mixture then gas phase. The jump on the enthalpy is increasing along time and is asymptotically the maximal jump h s g − h s ℓ . This is the critical case where mixture disappears asymptotically with a zero space derivative in the liquid at the transition.
• Figure 9b shows the same phenomena with a bigger diffusion in the gas phase leading to positive slope in the liquid at the transition for sufficient large times.
• Finally Figure 9c shows the particular case of a zero diffusion in the liquid phase leading to an increasing jump at the transition point to the gas phase. For sufficiently large time, the jump exceeds the value h s g − h s ℓ due to the lack of diffusion in the liquid phase which behaves like the mixture.