A theoretical model for double diffusive phenomena in cloudy convection

Using classical rheological principles, a model is proposed to depict the molecular diffusion in a moistsaturated dissipative atmosphere: due to the saturation condition existing between water vapor and liquid water in the medium, the equations are those of a double diffusive phenomenon with Dufour effect. The double diffusivity is important because of the huge diffusivity difference between the liquid phase and the gaseous phase. Reduced equations are constructed and are then applied to describe the linear free convection of a thin cloudy layer bounded by two free surfaces. The problem is solved with respect to two destabilizing parameters, a Rayleigh number Ra and a moist Rayleigh numberRh. Two instabilities may occur: (i) oscillatory modes, which exist for sufficiently large values of the Rayleigh number: these modes generalize the static instability of the medium; (ii) stationary modes, which mainly occur when the moist Rayleigh number is negative. These modes are due to the molecular diffusion, and exist even when the medium is statically stable: the corresponding motions describe, in the moist-saturated air, configurations such as “fleecy clouds”. Growth rates are determined at the instability threshold for the two modes of instability occurring in the process. The case of vanishing moisture concentration is considered: the oscillatory unstable case appears as a singular perturbation (due to the moisture) of the stationary unstable state of the Rayleigh-B énard convection in pure fluid, and, more generally, as the dynamical perturbation of the static instability. The convective behaviour of a cloud in the air at rest is then examined: the instability of the cloud is mainly due to moisture, while the instability of the surrounding air is mainly due to heating.


Introduction
The role of the moisture on the dynamics of atmospheric air has been studied several years ago, from various points of view.For instance, Dudis (1972), Einaudi and Lalas (1973), Durran and Klemp (1982a, b;1983), considered its influence on the Brunt-Väisälä frequency and inviscid flows (such as mountain lee waves).Einaudi and Lalas and later, Durran and Klemp, focused their attention on media with a nonconstant distribution of moisture with respect to the altitude.Deardorff (1976Deardorff ( , 1980)), and Betts (1982) constructed reduced schemes using the so-called "liquid-water moist potential temperature": such schemes were used by Bougeault (1981a, b), in order to study the moist atmospheric turbulence.From another point of view, Kuo (1961Kuo ( , 1965)), Ogura (1963) studied convective instability in dissipative media.More recently, Bretherton and Smolarkiewicz (1987, 1988, 1989) considered the motion, the appearance and the disappearance of clouds under moist convection effects.Some consequences of the diffusion taking into account the microstructure of the medium have also been exhibited by Kambe andTakaki (1975), andMerceret (1977).In these works, the double diffusive characteristics of the medium are either absent (when dissipation is neglected in the medium), or not really taken into account.
In fact, it has been admitted a long time ago, that, since the diffusivities of dry air and water vapor are almost equal, double diffusion has no important effect on the atmospheric flows (see, for instance, Huppert and Turner, 1981).This property is valid in unsaturated atmosphere, where the molecular diffusion equation takes the classical form of the binary diffusion of one constituent (say, the water vapor) in the whole mixture.The diffusion, in this case, follows the well-known Fick's law.In saturated atmosphere, the diffusivity of the liquid phase also plays a role in the phenomenon.This diffusivity is not of the same order of magnitude as that of the gaseous phase, and the laws governing the phenomenon are not yet well known.In practice, two models of molecular diffusion may be found in the literature: (i) a binary Fick's law is chosen for the diffusion of the total water in the total mixture {air + total water} (e.g. in Bretherton's papers); (ii) a binary Fick's law is chosen for the diffusion of water vapor in the gaseous phase, i.e. {water vapor + dry air}, the diffusion of the liquid water being neglected: this law (see Hijikata and Mori, 1973) is also adopted in works dealing with industrial thermodynamics, and a variant may be found in Deardorff's papers and Bougeault's papers.Hence, the first question set in the present paper is that of the choice of a law of molecular diffusion: by applying the basic principles of thermodynamics and of rheology, we construct a diffusion law stating precisely the validity of Bougeault's formulae.A simplified form of this law is also discussed in the case of weak moisture.
Another singularity arises in the analysis of the phenomenon of moist convection.Since the water concentration, say q v , is small, the Schmidt number S * , characterizing the molecular diffusion, plays an effective role in the equations through the combination 1/S = q v /S * only; this number is very small when S * is of order unity (in practice, in this mixture, S * ≈ 0.721).However, the factor 1/S multiplies a laplacian, so that even for small 1/S, the molecular diffusion acts through some boundary layers, and, in a stability problem, such boundary layers may originate specific instabilities.
The following paper is organized as follows: in Sect.2, the equations of motion are set up using a classical mixture theory; the basic equations are written using a model of continuous medium, with the diffusive dissipation being given using the principles of thermodynamics.The major simplifying assumption is that the condensed water is treated as an aerosol (its pressure is neglected in the whole mixture, so that we are led to describe it as a polytropic gas whose adiabatic constant is zero); the additional assumption of negligible diffusion of the liquid water is not needed.A second basic assumption is that of a saturated medium.These assumptions are classical in non-precipitating cloud theories.Because the molecular diffusion in saturated mixtures is not very well known, the molecular diffusion law is especially examined.It is shown that the phenomenon follows a generalized Fick's law with Dufour effect.In fact, the problem of ternary convection may be finally reduced to a convection in a binary mixture, and this reduced problem is that of a diffusion with Dufour effect.However, some properties of the ternary nature of the mixture remain present in the model, as will be seen later.
In the Sect.3, the linearized gravity wave equation governing the motions is derived.The study is made using the asymptotic frame of the Boussinesq approximation already used in previous papers (Bois, 1991(Bois, , 1994)).The use of this approximation first implies the existence of a particular solution describing a static state of the medium, depicted with the help of unique variable (the altitude referred to the "dry" atmospheric scale).Second, the characteristic vertical scale of the perturbation motion is small before the atmospheric scale, and the characteristic time of the motion is also scaled with the help of static data (namely t 0 = c 0 /g, where c 0 de-notes the speed of sound in the medium at the rest); this formulation follows the classical analysis of Spiegel and Veronis (1960).An advantage of using the asymptotic formulation of the Boussinesq approximation is to allow one to write the equations with the help of two scales (the fast scale, which is the scale of the convection, and the scale of the stratification of the medium, or slow scale); their ratio is the small Boussinesq parameter, say ε.After linearizing the equations around the static state, we obtain a wave equation generalizing to the considered medium the classical gravity wave equation used, for instance, in the Rayleigh-Bénard convection.The analogy with the equations of thermohaline convection is also noted.
In Sect.4, we consider shallow convection in the medium at rest.Since the considered problems are dissipative, two diffusion factors remain present in the reduced problem.In order to solve the singular perturbation occurring in the phenomenon, we first assume that the Schmidt number S * is small, of the same order as q v .We introduce the reduced Schmidt number S, which is, itself, of order unity.After solving the corresponding instability problem (so-called auxiliary problem) it is simpler to obtain the solution of the real problem by examining the behaviour of the solution of the auxiliary problem when S goes to infinity.By applying this procedure, a moist Rayleigh-Bénard problem is solved, and its solution is described with respect to two destabilizing parameters, namely the Rayleigh number R a , and a moist Rayleigh number R h , itself proportional to the total water gradient.Finally, the following properties are derived: first, because of the permanent exchange of mass between the two liquid phases, a stationary instability mainly occurs for a negative moist Rayleigh number, even when R a is, itself, negative.The corresponding cloud configurations appear in a first approximation, as motions of the dry air and the aerosol only.These motions describe fleecy cloud configurations.Such motions are partially located in the statically stable region of the (R a , R h ) plane, in the same manner as the salt fingers found in thermohaline convection (see Baines and Gill, 1969).Second, an oscillatory instability appears in an angular region of the (R a , R h ) plane, itself located between the stationary stable region and the stationary unstable region, in such a manner that the instability threshold may be either stationary or oscillatory, except in the vicinity of the polycritical point, where the two instabilities may simultaneously exist.The oscillatory instability is related to the static instability of the medium, while the stationary instability is due to the molecular diffusion.Starting from the results of this problem, we then deduce those of the corresponding problem when the Schmidt number is of order unity.The case of vanishing moisture is also examined, and it is shown that the classical Rayleigh-Bénard convection appears when the moisture vanishes, this state appearing as a limit case of the oscillatory unstable state already exhibited.The growth rate of the unstable modes is determined in the neighbourhood of the instability thresholds: this determination is straightforward for oscillatory modes.For stationary modes, because the bifurcation is singular (all wavelengths are unstable), this determination allows one to select the preferential instabilities occurring in the medium.To end with, we match the instability conditions in the cloud with the corresponding conditions in unsaturated air (clear air), so that we may conclude about the instability of a cloud surrounded by clear air.

Thermodynamics of moist saturated air
The medium is modelled as a ternary mixture: dry air, water vapor and condensed water.The constituents are identified by the subscripts g (dry air), v (water vapor), L (liquid water).The densities are denoted ρ g , ρ v , ρ L .The density of the medium being ρ = ρ g + ρ v + ρ L , the concentrations are denoted by q α = ρ α /ρ, and satisfy the relation q g + q v + q L = 1.We will also use the total water concentration q w = q v + q L (more generally, the subscript w will denote "total water").The velocity of each constituent is denoted u α , and the barycentric velocity of the fluid is u = q g u g + q v u v + q L u L .The diffusion velocities v α are then defined by v a = u α − u, and satisfy the relation (1) The specific variables defining the system satisfy the thermostatic Gibbs-Duhem equation where e denotes the internal energy per mass unit, and the g α 's denote the free enthalpies per mass unit of the constituents.For a reversible infinitesimal transformation of the medium, because of the saturation hypothesis, the free enthalpies g L and g v are equal.Hence, setting dq w = dq L + dq v and using dq g + dq w ≡ 0, Eq. ( 2) reads After Eq. ( 3) we have e = e(s, ρ, q w ), T = (∂e/∂s)(s, ρ, q w ) and p = −(∂e/∂υ)(s, ρ, qw), where υ = 1/ρ.Eliminating s between these relations, we have in the saturated mixture e = e(T , ρ, q w ), p = p(T , ρ, q w ). (4) In the sequel we assume that the mixture is an ideal mixture of polytropic gases.Moreover, following a classical approximation for the aerosols, we neglect the pressure of the liquid phase so that this phase is also a polytropic gas (with a zero polytropic constant).

The Clausius-Clapeyron equation
The saturation equation is the Clausius-Clapeyron equation.We introduce the latent heat of vaporization of the liquid, namely L v = h v − h L , (h α is the specific enthalpy of the phase α).A form of this equation, valid when the liquid phase is an aerosol (see, for instance, Zemansky, 1968), may be written as It is of interest to examine some consequences of the Clausius-Clapeyron formula, when the partial pressure of the liquid phase is neglected.The dry air, the water vapor and the whole mixture follow the equations of state, respectively, where R g and R v are two constants.The third relation ( 6) expresses Dalton's law.As said above, the partial pressure p L is neglected in Dalton's law.Let us now differentiate Dalton's law.After some calculation and using the Clausius-Clapeyron equation, we obtain the relation The relation (7) expresses the total water variations with respect to the variations of vapor, pressure and temperature.Now, denoting c p v and c p L the heat capacities of the water vapor and the liquid water, respectively, the latent heat of vaporization can be written as L v 0 and 0 being two constants.L v 0 is the latent heat at the temperature 0 .Introducing Eq. ( 8) in the Clausius-Clapeyron equation (written with the variables q v , ρ and T ), we obtain the relations where q 0 , ρ 0 , 0 are some constants.By eliminating ρ in (9) with the help of the third of (6), we deduce, after some calculation, the relations Since q v /q g denotes the concentration of water vapor in the gaseous phase, the pressure p s0 G(T ) is the saturation pressure of the water vapor for the corresponding vapor concentration, with the constant p s0 being the saturation pressure of the water vapor at the temperature 0 .The function Q s (p, T ) defined in Eqs. ( 10) is the saturation concentration of the water vapor in the gaseous mixture.Both Eqs. ( 9) and (10) may be taken as "equations of state" and used in order to eliminate secondary variables in the equations of motion.Extracting q v /q g as functions of p, ρ, and T from Eqs. (10) and the third of ( 6), an explicit expression of the equation of state p = p(q v , ρ, T ) could be obtained by inserting these expressions in Dalton's law.However, although these expressions theoretically could allow one to eliminate some unknowns in the equations of motion, it is simpler to keep the original variables in these equations.

The diffusion equations
We now consider a moving system.The barycentric motion of the system is described using the material derivative d/dt = ∂/∂t + u • ∇.First, consider the balance of mass for the dry air: this equation, (so-called molecular diffusion equation) reads The balance of mass for the whole mixture in the barycentric motion reads The diffusion velocities are themselves related to the gradients of concentrations by Fick's law.For the mixture considered here, this Fick's law (see Appendix A) reads where C is a positive constant coefficient.Because q v is small, an approximation of Eq. ( 13) is ρ g v g = C∇q v ; such equation has been used in Kubicki and Bois (2000).However, it will be seen later that, in convection problems, as we further investigate, this approximation is too rough if we look for a realistic solution.Hence, even for small q v , we will keep the whole relation ( 13) in what follows.The relation (13) inserted in Eq. ( 11) provides the new equation Now consider the first law of thermodynamics: for a moving mixture of dissipative fluids, this law, written in terms of enthalpy in place of the internal energy, after some calculation, provides the equation The coefficient c p denotes the heat capacity of the mixture with constant pressure and concentrations.On the right-hand side, the dissipation is the sum of the viscous dissipation , the thermal dissipation k T /ρ (k is the thermal conductivity), and the diffusive dissipation in the third term: the h α 's are the partial enthalpies of the constituents.
According with classical hypothesis of small diffusion velocities (Bowen, 1976), the thermodynamical potentials e, h, g are related to the partial potentials e α , h α , g α by the approximate relations e = q α e α , h = q α h α , g = q α g α .Note that (∂h/∂q α ) p,T = h α .Finally, using Eq. ( 5) and Eq. ( 14), we replace dq g by −(dq v + dq L ), and h v − h L by L v : the Eq. ( 15) reads, after some rearrangement, The motion of the medium is described by using the classical variables u, p, ρ, T , and two concentrations, say q v and q g .It is convenient to write the diffusion equations with the same variables: hence, replacing dq g by −(dq v + dq L ) in Eq. ( 14) and using Eq. ( 7) to eliminate dq L , Eq. ( 14) also reads The two Eqs.( 16) and ( 17) are two symmetrical forms of the diffusion equations: assuming that dp/dt is known in these equations (it is an approximate consequence of the Boussinesq equations, see Sect. 3) the real unknowns, in Eqs. ( 16) and ( 17), are dT /dt and dq v /dt.These equations can be solved in two independent forms dealing with linear combinations of these quantities.The molecular diffusion and the heat conduction appear in laplacians figuring in the righthand sides of these equations and cannot separately be considered; this property characterizes a double diffusive phenomenon.
Finally, the equations of motion for a saturated mixture are: the third Eq. ( 6), and Eqs. ( 9), ( 12), ( 16), ( 17), to which we must add the balance of momentum where µ 0 denotes the dynamic viscosity.Finally, we have a complete set of 8 scalar equations for the 8 unknowns u, p, ρ, T , q v , q g .

Validity conditions and nondimensional equations
The Boussinesq approximation refers to the motion to a static state of the medium, in which the variables depend only on one variable, namely the altitude scaled by the so-called atmospheric height H = p 0 /(ρ 0 g), where p 0 and ρ 0 denote characteristic values of pressure and density.
The general validity conditions of the Boussinesq approximation (Bois, 1991) first imply that the characteristic scale L of the motion is small before H : L/H = ε 1, and second, that the characteristic velocity scale U of the motion is related to ε by the equality U √ ρ 0 /p 0 = ε: that is equivalent to choosing a characteristic time of the motion t 0 = H √ ρ 0 /p 0 .We don't discuss the validity of these conditions, which are assumed to be true.We note, however, that, when rewritten for the significant nondimensional parameters of the problem (defined further, see relations ( 25), F = U/ √ gL is the Froude number of the flow), these conditions imply the relations The variables are now scaled by characteristic values: U, L, ρ 0 , 0 , c p g .The scaling pressure is p 0 = ρ 0 c p g 0 , and we have R g = (γ − 1)c p g /γ , γ being the adiabatic constant of the dry air.Rewritten with nondimensional variables the former equations read The Eq. ( 12) remains unchanged.We have set in the Eqs.( 20)-( 24) 0 and β being defined in Eqs. ( 9) and ( 10).The two diffusion parameters are the Prandtl number P r and the Schmidt number S * , which match the thermal dissipation (Prandtl) and the molecular dissipation (Schmidt) with the viscous dissipation (Reynolds number R e ).We have used the relation ( 8) and the relations The first two equations of ( 26) are valid for polytropic gases, the third equation results from Eq. ( 8).Since we are concerned only with shallow convection, it is possible to find equilibrium solutions to the system Eqs.( 20)-( 24).Hence, we consider an equilibrium state (u = 0), having a prescribed temperature distribution and a prescribed water vapor concentration, say Eqs. ( 20)-( 24) are satisfied (up to and including order ε) by pressure, density and dry air concentrations denoting with primes the derivatives d/dζ , these quantities satisfy the static equations The variable z is the physical vertical variable (say z ) scaled by L, while ζ (= εz) is the same variable z scaled by H .Note that Eq. ( 28) is nothing but the static form of the Eq. ( 7).

Boussinesq asymptotic expansion
Following the Boussinesq procedure (Bois 1991), we now expand the variables with respect to ε in the following form: We expand Eqs.( 12) and ( 20)-( 24) with respect to ε. Discarding the static terms, we obtain the following perturbation equations We have set in Eqs. ( 32)-( 35) Simplified forms of Eqs. ( 30)-( 35) may be looked at in two ways: (i) by assuming that the water concentration is weak (see below) and, (ii) by approximating the values of the coefficients by their values at a given level (shallow convection, see Sect. 4).

The reduced problem
In order to look for simplifying assumptions about the system (30)-( 35), we first note that the magnitude of both water vapor concentration and liquid water concentration are usually small, of the same order: their order of magnitude (namely the parameter q 0 introduced in Eqs. ( 9) and ( 23)) is, in practice, about 10 −3 .Hence, we expand the system with respect to q 0 : keeping only the leading terms of the expansion, many coefficients in Eqs. ( 30)-( 35) take simpler approximate forms.However, this simplification induces some singularities.The most important of these is that the double diffusive character of the system disappears when we neglect terms of order q 0 (the laplacians qv and qg figuring in Eqs. ( 32) and (33) become negligible) .In order to avoid this singularity, we now assume that S * is also small, of the same order q 0 as the water concentrations.Thus, we set where S is a modified Schmidt number, which is itself assumed of order unity.Since the Schmidt number appears only in the combination (1/S * ) qv , the assumption S = O(1) allows one to keep the laplacian in the simplified equations, and, further, when we let S go to infinity in the solutions of realistic problems, we can directly examine the influence of the singular perturbation induced by the double diffusivity when S * = O(1).For small q 0 , taking into account the assumption (37), the Eqs.( 32) and ( 33) read In Eqs. ( 38) and (39) we have approximated the variables by constant values when their expansions for small q 0 don't lead to singularities: for instance, for small q 0 , C p 0 and C g 0 are approximated by 1.The system [( 30), ( 31), ( 34), ( 35), ( 38), (39)] may be yet simplified by extracting the unknown qg (or, equivalently, the unknown qv ) from the algebraic Eq. ( 35), as a linear combination of the other variables.The remaining system is the reduced system associated with the problem.Note that, since we have allowed in the system the coefficients depending on the variable ζ , the equations are valid for deep convection, as well as for shallow convection.

Linearized equation
The preceding equations are now applied to the study of shallow free convection in a thin layer (a cloud).The level z = 0 is the mean altitude of the convecting cloud.The thickness of the layer is chosen as scaling length L in the equations.
Because of the weak thickness of the layer, and since we are interested only in terms of order 0 with respect to ε, all slowly varying coefficients of the Eqs.( 31)-( 36) and ( 39)-( 40) are approximated by their values at z = 0, namely Taking into account the assumption of small q 0 , we set qv = q 0 qv , qL = q 0 qL , The state of rest (namely ũ = 0, ρ = p = T = qv = 0) is an exact solution of the reduced system (30)-( 31)-( 34)-( 35)-( 38)-( 39).The linearized equations associated with the system, approximated by their first term with respect to q 0 , then read We have used, in the preceding equations, the notations Because qg figures in Eqs. ( 43)-( 47) only through the combination q 0 qg , the limit of Eqs. ( 43)-( 47) for vanishing q 0 is singular (8 equations for seven unknowns only).Because all boundary conditions are zero in a free convection problem, this singularity may be avoided by assuming that the variables ũ, p, ρ, T , qv are themselves of order q 0 , while qg is of order 1.Hence, we set The first approximation of Eq. ( 45) now takes the nondegenerate form Hence, after Eq. ( 44) we also have Finally, extracting q * v from Eq. (51), the system ( 42)-( 47) reduces to the following The system ( 52)-( 53) clearly exhibits the particularities of the problem: first, the problem is a double diffusive problem for a fictitious binary mixture whose equation of state is the third equation listed in Eq. ( 52).Second, due to the form of the right-hand side of the first equation listed in Eq. ( 53), the diffusion equations in this medium involve the Dufour effect (the time derivative ∂q * g /∂t doesn't depend on q * g but on T * only !).The significance of the first diffusion Eq. ( 53) may be understood by the help of the static Eq. ( 28): effectively, taking into account the shallow convection assumption and the smallness of q 0 , the leading approximation of this equation reads so that − 0 N 2 is nothing but the total water gradient in the medium at rest.Since Q w0 is of order q 0 , a rigorous study would impose to cancel − 0 N 2 (at the order 0 with respect to q 0 ) in the system (52)-( 53).In fact, it is equivalent to solve this system for arbitrary values of − 0 N 2 and, furthermore, to consider small values of this parameter; this procedure is adopted in the sequel.

Solution of the linear Rayleigh-Bénard problem with two free surfaces
The Rayleigh-Bénard problem with free surfaces is the simplest case of boundary conditions associated with the system (52)-( 53).In spite of its academic character, it allows one to simply exhibit analytical solutions (in fact, with the help of the only dispersion equation).The linear stability of the system may be studied by a normal modes analysis, by looking for solutions in the form u * = U(x, y, z)e ηt , with η = ξ +iω.Since we don't envisage boundary conditions depending on the horizontal variables, all modes are assumed in the form: where U, V , W etc., depend on z only, and where k, h, are the components of a horizontal wave number vector.We search the values of N 2 and such that all possible exponents η have a negative real part.First, inserting the variables (55) in the Eqs.( 52)-( 53), the system reduces, for the unknown W , to an equation of sixth order with constant coefficients.This equation itself takes a simpler form, by introducing the characteristic parameters R a is the Rayleigh number, R h is the moist Rayleigh number, τ is the ratio of the diffusivities of the medium.Moreover, we adopt the notations Rewritten with these notations, the equation satisfied by W reads The boundary conditions are free surfaces conditions at z = 0 and z = 1.If we require continuity of the temperature and the concentration at the boundaries, the other boundary conditions are standard (see Drazin and Reid, 1981), so that we assume The boundary conditions (59) are satisfied by particular solutions of (58) of the form W = W sin nπ z, n ≥ 1, where W = const.The dispersion equation takes the algebraic form where we have set Equation ( 60) is a third-degree equation with respect to σ , similar (but, however, simpler) to the equation from Baines and Gill (1969) for thermohaline convection.The thresholds are reached for Re(σ ) = 0, so that these thresholds are associated with stationary solutions or purely oscillatory solutions of Eq. ( 60): the procedure is standard (see Drazin and Reid, 1981) and leads to the following conclusions: (i) stationary solutions: Eq. ( 58) possesses stationary solutions if the condition is satisfied independently of n; this property characterizes a singular bifurcation, because all harmonics of a solution appear at the same threshold as the fundamental mode.Hence, the growth rates of the unstable solutions must effectively be computed in order to select the mode which effectively appears in an unstable process (see Sect. 4.3); (ii) oscillatory solutions: two conditions must be satisfied: first, denoting the pulsation of the solution (σ = i ), exists only if the inequality is satisfied.Second, in the half plane defined by the condition (63), 2 is, itself, positive, only if is satisfied.In the (R a , R h ) plane, both inequalities are simultaneously satisfied inside the angular sector U AY of Fig. 1.A is the polycritical point of the problem.The formulae ( 62) and ( 63) allow one to determine the coordinates (R a a , R ha ) of A, namely R a a = 27π 4 4 Since the left-hand side of Eq. ( 54) is of order q 0 , we see that the only realistic values of the number R h /(P r R 2 e ) are small of the same order.Hence, the solutions exhibited here are available only in the region of the (R a , R h ) plane, such as R h /(P r R 2 e ) = O(q 0 ) 1.Note that this condition is not very restricting, because, in general, the Reynolds number is large.The above analysis also shows that the least value of R h in an oscillatory unstable state is that reached at the polycritical point A. Hence, in order to satisfy the condition R h /(P r R 2 e ) 1, it is of interest to study the possible positions of A when the intensity q 0 of the moisture (even remaining small) varies.We have plotted, in Fig. 2, the way followed by A in the (R a , R h ) plane when τ varies from 0 to + ∞.Eliminating τ between the two formulae (65), the equation of this curve may be easily derived, namely }/(2µ).For infinite τ, A is located at the point C = 27π 4 /4(≈ 657) of the R a -axis: this point is just the threshold of the stationary instability in a pure fluid.For smaller values, and, in particular, for vanishing τ , the trajectory is asymptote to the straight line of equation R a a = µR ha + 27π 4 µ/[4(µ − 1)].Hence, for vanishing q 0 (large values of τ ), A is located near the R aaxis, and for large values of q 0 (small values of τ ), A is far from the R a -axis.However, the inequality R h /(P r R 2 e )) 1 remains satisfied due to the presence of the Reynolds number in this condition: in the considered problem, the Reynolds number characterizes the thickness of the unstable layer, and this thickness should be very small (practically less than a few meters), in order that the inequality is not satisfied.In In order to exhibit the different regions more explicitly, the scales are not respected in this figure.The realistic instability is located in the shaded region R h = O(q 0 R 2 e ), so that the polycritical point may be either inside or outside the instability domain.The straight line R a = R h delineates the statically stable region.The numerical values chosen for determine the polycritical point A are (units M.K.S.A.): L 0 = 2600 10 3 J/kg, 0 = 288 • K, c p v = 1004 J/(kg.• K), R v = 464 J/(kg.• K), R g = 287 J/(kg.• K), S = S * /q 0 = 721, P r = 0, 76, q 0 = 10 −3 .practical problems, because the point A is located in a realistic region of the (R a , R h ) plane, oscillatory instability may also exist.

Stable and unstable regions and growth rates at the instability thresholds
The preceding analysis, although allowing one to determine the instability thresholds, doesn't place in evidence the stable and the unstable regions of the (R a , R h ) plane.In order to localize these stable and unstable regions, because of the indetermination led by the equality (62), it is convenient to examine the real part of σ in the neighbourhood of the thresholds: (i) Stationary threshold XA: as a preliminary remark, we note that a procedure followed by Baines and Gill (1969) may be applied here: setting Eq. ( 60) takes the canonical form Equation ( 67) is an algebraic equation of the third degree for θ .By examining the sum and the product of the roots, it is easily verified that, when the inequality is satisfied, this equation possesses one real positive root only.Hence, the half plane defined by Eq. ( 68) is an unstable region with unstable direct modes: one mode exactly for a given value of (R a , R h ), i.e. for fixed R a , R h , and K.In order to determine the growth rates and the wavelengths of those unstable modes (and the most unstable modes themselves), we linearize Eq. ( 67) near the threshold defined by Eq. ( 68): hence, we set where R h0 and R a 0 satisfy the relation and where R h and θ are assumed small.After linearizing Eq. ( 67), and taking into account the relation (69), we obtain the following relation between R h and θ After the second equation listed in (65), because R h0 < R ha along the half-straight line XA, the denominator of Eq. ( 71) is negative along XA.Hence, for negative R h , θ is real and positive: the stable region is located over XA, and the unstable region is located under XA (Fig. 1).Moreover, Eq. ( 71) allows one to find the most unstable mode in the following manner: first, we deduce from Eq. ( 71) by reinserting the variables σ and R h with the help of Eq. ( 66) By studying the variations of σ with respect to n, with K being fixed, we deduce from Eq. ( 72) (see Appendix B) that the greatest value of σ is obtained when n = 1; hence, the most unstable mode is always the fundamental mode.Furthermore, restricting the study to the case n = 1, we deduce, after some calculation, the wave number K max for which the maximum of σ , say σ max , is reached: K max is the only root of the algebraic equation (the calculations are given in the Appendix B).The corresponding σ max is given by the formula (72).K max is a function of R h0 (or, equivalently R a 0 ); hence, it varies along the whole straight line XA, up to A. For instance, for the numerical values used in Fig. 1, the curves are drawn on Fig. 3, where Fig. 3a shows that K max is always greater than the value π/ √ 2(≈ 2.221) of the wave number of a regular Bénard convection.K max decreases when R a 0 increases, until the value is reached at the critical value of R a 0 .Figure 3b shows that the growth rate of the fundamental mode is much larger than the first harmonic (n = 2), and that, for n ≥ 3, these growth rates become negligible.Of course the formula (72) becomes invalid in the neighbourhood of the critical point: near this point we must use a second order expansion, which leads to an expression of σ max proportional to ( R h ) 1/2 (see later).
(ii) Oscillatory threshold AU .Since the instability near the threshold AU is a regular problem, the instability always occurs for the fundamental mode n = 1, and the corresponding wave number K c = π/ √ 2: this situation is analogous to that of the classical Rayleigh-Bénard convection.However, the preceding analysis is useful in order to determine the growth rates of the oscillating unstable modes; we set where (R a 0 , R h0 ) are the coordinates of a point of AU , and R a denotes a small variation of R a from an arbitrary value R a 0 .We set, moreover, where φ 0 is real, and we assume that θ is small.Since iφ 0 is an exact imaginary root of Eq. ( 67), R a 0 and R h0 satisfy the relations Now, linearizing Eq. ( 67) for the unknown θ and using Eq. ( 76) and the first of ( 77), we obtain after some calculation After Eq. ( 79), for positive R a , ψ is always positive along the threshold AU .In fact, ψ is nothing but some scaling of the real growth rate Re( σ ); hence, using Eq. ( 66) to reintroduce the variables σ, R a , R h , K, and noticing that, for the oscillatory instability, We have drawn on Fig. 4, for the same numerical values as in the preceding figures, the curve | σ/ R a | = F (R a 0 ) along the half straight line U A. Note that φ 0 2 can be extracted, in terms of R a 0 , from the relations (77).In the same manner as in Fig. 3, the neighbourhood of the critical point is singular, but, after Eq. ( 79), only the phase becomes singular at this point.
(iii) Neighbourhood of the critical point.Near the critical point, the linear expansion of the dispersion equation is invalid, for the stationary case, as well as for the oscillatory case.Since Re( σ ) is small in this neighbourhood, we now set where R a , R h , ψ, φ are assumed small.Since the critical point is a stationary threshold and an oscillatory threshold simultaneously, the following relations ( are satisfied.We first rewrite Eq. ( 67) by separating its real part and its imaginary part, taking into account the relations (82).Hence, we obtain the two following real equations According to Eq. ( 84), we separately consider the cases φ = 0 and φ = 0. We obtain the following conclusions: (a) if φ = 0 : φ 2 may be extracted from Eq. ( 84) and inserted in Eq. ( 83).This equation is then an equation of the third degree for the unknown ψ.Since we look for solutions of Eq. ( 83) only for small R h and R a , approximated values of the roots can be looked for.We obtain one positive root if the condition τ P r R h − [ 0 (µ − 1) + τ + τ P r ] R a is negative; this condition corresponds to a point (R a , R h ) located beyond the oscillatory threshold V U of the (R a , R h ) plane (see Fig. 1).In this half plane, the condition φ 2 ≥ 0 is satisfied if the inequality ( 0 µ + τ ) R h − 0 R a ≥ 0 is satisfied.Finally, the corresponding modes exist if the figurative point (R a , R h ) is located in the oscillatory region previously determined on Fig. 1.The values of ψ, φ are (b) if φ = 0 : ψ is directly given by Eq. ( 83), which possesses a solution only if the inequality ( 0 µ + τ ) R h − 0 R a ≤ 0 is satisfied; this corresponds to a figurative point inside the stationary unstable region of Fig. 1.In this case the solution reads After Eqs. ( 85) and ( 87), the order of magnitude of the growth rates is O( 0 R a − ( 0 µ + τ ) R h ) 1/2 in the stationary region, and 0(τ [ 0 (µ − 1) + τ + τ P r ] R a − τ 2 P r R h ) in the oscillatory region.The three expansions (85), ( 86) and ( 87) may be asymptotically matched without difficulty, either with Eq. ( 71) (for the stationary instability), or with Eq. ( 79) (for the oscillatory instability).The behaviour of the oscillations is given, in the latter case, by Eq. ( 86).Finally, in the (R a , R h ) plane, the stable region is the region inside the angular domain XAU (Fig. 1).The threshold XA leads to direct instabilities (stationary onset).
The threshold AU leads to unstable oscillating modes (oscillatory onset).The growth rates of these modes are those calculated above.

The fleecy configurations of clouds
It is of interest, in order to understand the role of the diffusion in the medium, to compare the instability regions in the (R a , R h ) plane, to the static instability regions.In fact, it is well-known (see Durran andKlemp, 1982a or Bois 1994) that the moist Brunt-Väisälä frequency of the medium, say N m , is given by the formula 88) so that the neutral line of static stability, in the (R a , R h ) plane, is (at order q 0 ) the straight line R h = R a (the dotted line of Fig. 1).Hence, the stationary instability is partially located in the statically stable region (the angular sector XOT of Fig. 1).In this region, the instability regime which sets in is entirely analogous to the well-known "salt fingers regime", existing in the thermohaline convective instability (see Nield, 1967 or Huppert andTurner, 1981).Indeed, as said formerly, this regime is a very slow motion of dry air and liquid water only; effectively, the variables qg and qL defined in Eq. ( 29) are of order q 0 , while the variable qv , (after Eq. 49) is of order q 0 2 .Hence, rather than "moisture fingers", the solution depicts, in fact, the fleecy appearance of clouds in an air at the rest.

The case S
As we have seen above, the case S * = O(1) may be deduced from the preceding analysis by examining the behaviour of its solution when S → ∞; in fact, this case corresponds to τ → ∞ in the former formulae.When τ → ∞, the point A of the (R a , R h ) instability moves along the curve CD of this plane, towards C (see Fig. 2 and Fig. 5b).This point becomes exactly the point C when τ = +∞.Since the instability thresholds are delineated by straight lines joining A and O itself, the corresponding positions of these thresholds may be very easily followed.More precisely, for infinite τ the inequalities ( 63) and (64) become while the equality (62) now reads After Eq. ( 90), at order 0 with respect to q 0 , the stationary instability threshold becomes the R a -axis itself: the stable region is the half plane R h > 0, the unstable region is the half plane R h < 0. From Eq. ( 89), the oscillatory unstable region is the angle U CY of the (R a , R h ) plane (see Fig. 5c), where C is the point R a = R a 0 = 27π 4 /4 of the R a -axis.In fact, for vanishing moisture, i.e. q 0 → 0, the R a -axis for R a > R a 0 is not a singular bifurcation line.At order 0 with respect to q 0 (Fig. 5c) the thickness of the realistic region of the flow is zero, so that this region reduces to the R a -axis itself; stationary instability is, in this case, irrelevant.On the contrary, oscillatory instability is realistic: the stable part of the R a -axis is R a < R a 0 , the unstable part is R a > R a 0 .Moreover, since the unstable region corresponds to the boundary 2 = 0 of the "oscillatory" instability, this instability is, in fact, stationary.We recognize the instability scheme of the classical Rayleigh-Bénard convection, where the instability is due to the only Rayleigh number.This behaviour is regular, since, for vanishing q 0 , the fluid becomes a pure fluid.At order 1 with respect to q 0 (Fig. 5b), the realistic region of the flow is a very narrow strip around the R a -axis, and the critical point A is distinct of C. The case τ → ∞ also allows one to understand the singularity occurring if we directly assume {q 0 1, S * = O(1)} in Eqs. ( 20) to ( 24): in this latter case Eq. ( 53) is replaced by and the problem is now governed by the system [( 52), ( 91)], which reduces to the equation Let us look for stationary solutions of Eq. ( 92): such solutions (except the state of rest) exist only if R h = 0, i.e if − 0 N 2 = 0.The first equation ( 91) degenerates in this case, so that the system (52)-( 91) is now a system of 6 scalar equations for seven unknowns.However, we can note that, in this case, the changes of variables defined by (49) are superfluous.The original variables ũ, p, ρ, T are solutions of the system (42), (43), and the two equations The first equation listed in Eq. ( 93) is the leading approximation of Eq. ( 45), and the second is the stationary version of the second equation in Eq. ( 91).Equations ( 42), ( 43), ( 93) constitute a complete, regular set for ũ, p, ρ, T , which are solutions of an exact stationary Rayleigh-Bénard problem.
Figure 5 C Fig. 5. Successive positions of the oscillatory instability threshold for decreasing q 0 .When q 0 → 0 (and, hence, τ → +∞), the thickness of the realistic region decreases so that, for vanishing q 0 , the oscillatory instability region progressively goes to the stationary instability region of the classical Rayleigh-Bénard convection (in order to exhibit the different regions more explicitly, the scales are not respected on the figure).
Finally, in this case, a stationary solution of the full problem exists only if R a ≥ 27π 4 /4, instead of the full R a -axis, as shown above.Contrary to the stationary solutions, the oscillatory solutions are not singular in this case.Moreover, if we let the heat conductivity go to zero (i.e. for P r → ∞), taking into account the definitions (56), the oscillatory threshold defined by Eq. ( 89) becomes −N 2 = −( − 0 N 2 ).Furthermore, let R e go to infinity and the frequencies of the oscillating modes go to zero (for instance, after Eq. 57), so that these modes become stationary, and the condition 2 > 0 has no sense.The oscillatory instability threshold is, in this case, the whole straight line U V of Fig. 1, i.e. exactly the static threshold.Hence, the oscillatory instability is the regular dynamical degeneracy of the static instability of the medium.

Convection in unsaturated air
Because a cloud is always confined between layers of clear air, it is of interest to match the results of the preceding analysis with those of the corresponding analysis in unsaturated air: it is simple to repeat the analysis in this case.In such a medium, Eq. ( 3) remains valid, taking into account that, now, q v = q w .The diffusion Eq. ( 11) also remains valid, but Eq. ( 16) is replaced by The scaled form of the diffusion equations now read Since there is not a change of phase in the medium, the static Eq. ( 28) is replaced by the equation The Boussinesq analysis of the problem does not require one to use the parameter S instead of S * .The reduced linearized Eqs. ( 42) and ( 43) remain valid, while Eqs.( 45), ( 46), ( 47) are replaced by Equations ( 99) and ( 100) explicitly show the disappearance of double diffusivity when S * = P r , the operator applied to T and qv being, in this case, the same in both Eqs. ( 99) and ( 100).The singularity of Eq. ( 98) is the same as that of Eq. ( 45), and may be similarly treated.Finally, defining the parameters and looking for solutions in the form (55), we obtain for W the equation Equation ( 102) with boundary conditions (Eq.59) may be discussed as follows: (i) stationary solutions: such solutions exist if the condition is satisfied.In the (R a , R h ) plane, the inequality (103) defines a half plane bounded by the straight line X Y (see Fig. 6a), with the unstable region mainly corresponding to positive values of R a and positive values of R h ; (ii) oscillatory solutions: denoting the pulsation (σ = i ) two conditions must be satisfied: first, 2 being given, solutions exist only if the inequality is satisfied.Second, in the half plane defined by Eq. ( 104), 2 is, itself, positive only if the inequality is satisfied.In the (R a , R h ) plane, both conditions ( 104) and ( 105) are simultaneously satisfied inside the angular sector U A W of Fig. 6.The point A , which also lies on the straight line X Y , is the polycritical point.Figure 6 shows the behaviour of the system when the parameters S * and P r are almost equal (here: S * = 0.72, P r = 0.76).In Fig. 6a the variables are R a and R h The coordinates of A are very large (A is far from the origin), and the angle U A W is very sharp.Moreover, 2 is small inside the domain U A W .In Fig. 6b the variables are R a and R h (itself defined in Eq. ( 56): R h is a relevant variable if we consider a problem where both unsaturated and saturated fluids coexist, such as the convection of a cloud in unsaturated atmosphere.For instance, if we consider the convective instability of a cloud in order to be able to match together stable and unstable regions of both saturated and unsaturated cases, as it occurs in the coupled instability of two layers of fluid.In the latter case, three values of R e are considered: (i) R e = 1 (stable region inside the angle XF 1 Y 1 ); (ii) R e = 150 (stable region inside the angle XF 2 Y 2 ); (iii) R e = 300 (stable region inside the sector XF 3 A 3 Y 3 ): in that case, the instability of the clear air may be oscillatory (threshold surrounded by two layers of unsaturated air, the two stability diagrams of the clear air layers and of the cloud itself may be superimposed (Fig. 6b); the stability is confined in the region XF Y or XF A Y of this figure, according to the values of the Reynolds number.The figure also shows that the instability of the cloud is mainly due to moisture, while the heating first destabilizes the surrounding air.

Concluding remarks
We have exhibited, in this paper, some singularities of the instability phenomena related to the double diffusive structure of the moist-saturated air.The most important conclusion concerns the law of molecular diffusion in the medium: following Onsager's assumptions, generalized expression of Fick's law of diffusion is given (Eq.13).Furthermore, with the magnitude q 0 of the water concentration assumed to be small, the cases where this law may be simplified or not simplified are also studied.With the parameter q 0 being taken as a small parameter, we have made an asymptotic expansion of the equations with respect to that parameter.It is interesting to note that the method we used in order to derive Fick's law (13) also indicates how one would arrive at other laws of diffusion (in particular viscoplastic diffusion: such a law is not necessarily irrelevant if we consider that, very often, the motion of clouds seem like rigid, solid motions).
From a physical point of view, the main conclusions are the following: (i) There exist stationary, unstable states, analogous to the salt fingers of the thermohaline convection: in the present case, these states describe fleecy clouds (rolls or cells), and mainly involve motions of the dry phase and the liquid phase of the mixture.These stationary states are due to the combined influences of molecular diffusion in the system and change of phase.
(ii) Oscillatory instability may occur, mainly because of the heating (Rayleigh number).The classical destabilizing influence of the Rayleigh number in a pure fluid is the limit of the oscillatory instability for vanishing moisture.Moreover, it is this instability which generalizes in a dissipative medium the eventual statically (nondissipative) instability of the medium.
(iii) By matching the results of a saturated instability with those of the corresponding unsaturated instability, we have shown, as a application, that stationary instability of a cloud can develop in a stable unsaturated atmosphere, mainly because of the moisture gradient, while the surrounding air becomes unstable, mainly because of the temperature gradient.
From a mathematical point of view, although the dispersion equation is of the sixth degree (instead of the eighth degree) this problem is more singular than the analogous problem of the thermohaline convection in the oceans.Physically, this singularity expresses that the wavelengths of unstable modes strongly depend on the values of the Rayleigh number and of the moist Rayleigh number in convecting modes.Some assumptions necessary for the modelling have been made in the paper: the first, which is the assumption of a small Schmidt number S * , is used in order to find Eq.( 58).This assumption is only an artifice related to the mathematical procedure.The Schmidt number S is a Schmidt number defined with respect to the total water density (instead of the total density of the mixture).
A second mathematical assumption is that of the "Boussinesq free surfaces" bounding the medium.This assumption facilitates the determination of the instability thresholds, but also allows one to qualitatively estimate the behaviours of the solutions of problems involving other boundary conditions.The last assumption, which is that of a small q 0 , is straightforward: this assumption is classically used, supplemented by the assumption of a large 0 , in such a manner that 0 2 q 0 remains of order unity (Einaudi and Lalas, 1973); such a new assumption would not change the basis of our analysis.

Appendix A Derivation of the molecular diffusion
Eq. ( 13) The diffusion velocities in a fluid mixture are related to the gradients of concentrations by phenomenological relations (in classical mixtures: Fick's law).The general way to obtain such relations consists of calculating the rate of the entropy production (the dissipation in the medium): by assuming that this dissipation is a quadratic positive form (it is the so-called Onsager hypothesis), we obtain the most usual correspondences between the variables.In the present medium, denoting the dissipation by , we obtain after some calculation where τ ij denotes the stress tensor, D ij denotes the deformation rate tensor, ∇ T denotes a gradient at constant temperature, q denotes the heat flux through the medium, and q' is related to q by the formula q = q − ρ g v g (h g − h v ). (A2) The relation (A1), which is not straightforward, is derived in Bois (2002).The three dissipations, v , t , d , are the viscous dissipation, the thermal dissipation, and the dissipation by molecular diffusion.The second law of thermodynamics stipules that must be positive for any thermodynamical process applied to the medium.The three dissipations v , t , d are written, in Eq. (A1), using independent panels of variables, so that they must separately be positive: the Onsager hypothesis corresponds to the simplest case where this condition is satisfied.After this assumption, the dissipation appears as a positive quadratic form, so that the variables figuring in the quantities v , t , d are related by linear correspondences: first, writing v as quadratic form leads to the classical Navier-Stokes equations for the whole mixture; second, the examination of t leads to the Fourier law, but this law does not affect q but q'; hence, q'= −k∇T (k is the thermal conductivity) ; third, writing d as a quadratic form of its arguments provides the generalized Fick's law, namely ρ g v g = −D∇ T (g g − g v ) = −D ∂ ∂q g (g g − g v ) T ,q v ∇q g − D ∂ ∂q v (g g − g v ) T ,q g ∇q v ,(A3) where ∇ T denotes a gradient taken at constant temperature, and D is a scalar coefficient.The last expression (A3) results from the Clausius-Clapeyron relation (9).Furthermore, the Fourier's law joined to (A3) yields

Fig. 1 .
Fig.1.Linear instability diagram of moist-saturated air in the (R a , R h ) plane.In order to exhibit the different regions more explicitly, the scales are not respected in this figure.The realistic instability is located in the shaded region R h = O(q 0 R 2 e ), so that the polycritical point may be either inside or outside the instability domain.The straight line R a = R h delineates the statically stable region.The numerical values chosen for determine the polycritical point A are (units M.K.S.A.): L 0 = 2600 10 3 J/kg, 0 = 288 • K, c p v = 1004 J/(kg.• K), R v = 464 J/(kg.• K), R g = 287 J/(kg.• K), S = S * /q 0 = 721, P r = 0, 76, q 0 = 10 −3 .

Fig. 3 .
Figure 3Fig.3. Wave number K max (a) and growth rates (b) of the three first stationary unstable modes.The numerical values are the same as in Fig. 1 (the scales are not respected on the figure for R a 0 ).The x-axis is taken along the stationary threshold, the variable being the Rayleigh number R a 0 .The figure is invalid in the vicinity of the asymptotes and, in particular, R a a (shaded region), because the linearized formula (71) becomes invalid.It is also invalid when R a > R a a , because this region is already unstable.When R a 0 → −∞, | σ/ R h | → 0 very slowly (as |R a 0 | −1/2 ).

Fig. 4 .
Figure 4Fig.4.Growth rate of the oscillatory most unstable mode.The xaxis is taken along the oscillatory threshold, the variable being the moist Rayleigh number R h0 .The numerical values are the same as in Fig.1.When R a 0 → ∞, Re( σ )/ R h → 0 as R a 0 −1 .

Fig. 6 .
Fig. 6.Instability diagram in clear air: (a) in the (R a , R h ) plane; (b) in the (R a , R h ) plane.The numerical values are the same as in Fig. 1.In Fig.4athe stability parameter R h is defined with respect to the total moisture of the medium.Because the Schmidt number is very near the Prandlt number, the double diffusion is weak (almost absent) in this figure, the oscillatory instability being confined in a very sharp angular sector U A W .In Fig.6b, the parameter R h = −R h + 0 R a + P r R e