Theoretical analysis

Theoretical analysis                           


The first part of this project aims to achieve a theoretical analysis of the fuel layer evaporation in the container. The next sub-parts correspond to the steps that we have followed to determine the equations system solving the problem.

First, we will show the assumptions that we take to model the problem, then we will study some particular cases before considering the evaporation in the whole container.













   Preliminary studies on a particular case

This first study involves to work on evaporation in enclosure where a  pure macro-layer is lied. This enclosure is considered adiabatic except the bottom wall where its temperature is fixed. For the first part, the study is in 2D but the edge effect along the vertical wall is neglected.

The second study involves to work on a particular case with different boundary layer : thermal and mass boundary layer.


Problem modeling

Problem modeling                             


       1- Lists of symbols & subscripts


Here are the list of all symbols and subscripts used in the study.

  • $h$    liquid thickness
  • $L$    container height dimension
  • $l_c$    capillary length
  • $a$    depth
  • $b$    length
  • $S_{flux}=a b$    flux surface ($m^2$)
  • $\rho$    density of the gas mixture ($kg/m^3$)
  • $P$    pressure ($Pa$)
  • $u$    the projection of velocity vector along the x-axis ($m/s$)
  • $\mu$    the dynamic viscosity of the gas mixture  ($Pa.s$)
  • $\nu$    the cinematic viscosity of the gas mixture ($m^2/s$)
  • $v$    the projection of velocity vector along the y-axis ($m/s$)
  • $g$    standard gravity ($m/s^2$)
  • $C_p$    specific heat ($^{-1}.K^{-1}$)
  • $T$    temperature of the mixture ($K$)
  • $\lambda$    the thermal conductivity ($W/m.K$)
  • $D$    diffusion coefficient of the vapour in the air ($m^2/s$)
  • $\rho_S$    mass concentration of vapour at the interface ($kg/m^3$)
  • $\omega$    mass fraction of vapour in the air with $\omega= \frac{\rho_v}{\rho_v+\rho_a}$
  • $\dot{m}$    mass flow rate ($kg/s$)
  • $h_{lg}$    latent heat ($kJ/kg$)
  • $h_c$    thermal exchange coefficient
  • $h_m$    mass exchange coefficient
  • $D_h$    hydraulic dimension


                Subscript :

  • $v$    for the vapour
  • $S$    for the vapour at the interface
  • $in$    for the input condition
  • $w$    for the wall condition
  • $l$   for the liquid


        2- Assumptions related to the vapour- state :

We consider a layer of liquid in a ventilated box. The following figure represents the problem:



                                     Figure: Problem modeling & Simplifications


During this study, we simplify the problem by considering the following assumptions:

  • Bidirectional incompressible turbulent flow.
  • Liquid film thin in relation to the reservoir but sufficient thick for neglect the interfacial resistance.
  • Air charged in vapour is in thermodynamic balance : the phase change happens in saturation conditions.
  • Viscous dissipation and pressure work is neglected.
  • Soret effect ( temperature gradient dependence for the mass flow) and Dufour effect (mass gradient dependence for the heat flux) neglected.
  • Radiation phenomena neglected.
  • Boussinesq approximation considered.

In further, the conservation equations are presented below.




       1- Non-dimensional equations

For a greatest lisibility of the mechanisms, the non dimensional equation can be written as :

$$X=\frac{x}{L} , Y=\frac{y}{L} , U=\frac{u}{U_{in}} , V=\frac{v}{V_{in}}, D_{hydraulic}=L$$ ,

$$ \theta=\frac{T-T_{in}}{T_s-T_{in}} ,  C=\frac{\omega-\omega_{in}}{\omega_s-\omega_{in}}$$

The non dimensional equations become :

$$U\frac{\partial U}{\partial X}+V\frac{\partial U}{\partial Y}=\frac{1}{Re}\left(\frac{\partial^{2}U}{\partial X^{2}}+\frac{\partial^{2}U}{\partial Y^{2}}\right)    (1*)$$

$$U\frac{\partial V}{\partial X}+V\frac{\partial V}{\partial Y}=-\frac{\partial P_{m}}{\partial Y}+\frac{1}{Re}\left(\frac{\partial^{2}V}{\partial X^{2}}\right)      (2*)$$

$$U\frac{\partial \theta}{\partial X}+V\frac{\partial \theta}{\partial Y}=\frac{1}{RePr}\left(\frac{\partial^{2}\theta}{\partial X^{2}}+\frac{\partial^{2}\theta}{\partial Y^{2}}\right)    (3*)$$

$$U\frac{\partial C}{\partial X}+V\frac{\partial C}{\partial Y}=\frac{1}{ReSc}\left(\frac{\partial^{2}C}{\partial X^{2}}+\frac{\partial^{2}C}{\partial Y^{2}}\right)      (4*)$$


        2- Non-dimensional numbers

Some interesting numbers appear in this non-dimensional equation.

  • Reynolds number

The Reynolds number is characteristic of pure flow conditions, inertial effect compared with momentum diffusion:


  • Prandlt number

Another one is the Prandlt number. It appears in thermal equation and it compares the momentum diffusion with the thermal diffusion

$$Pr=\frac{\frac{\mu}{\rho_{in}}}{\frac{\lambda}{\rho C_p}}$$

  • Schimdt number

The Schimdt number, quiet close with the Prandt number but it appear in mass equation. It compares the movement quantity diffusion with the molecular diffusivity.


In addition the expression of the thermal flows through the liquid film has to be written clearly because the study of the heat transfers can explain the energy transfer between the liquid phase and the vapour phase.

Firstly, the sensible heat flux can be modeled as:

$$ q''_s=-\lambda_{liq} \left(\frac{\partial T}{\partial y}\right)_{y=h}          (5) $$

That is the energy transfer through the thickness of the liquid  resulting from temperature gradient. We can suppose in this case the convection is neglected.


Then the local energy balance on the interface is :

$$q''_{evap}=q''_{fluid}+q''_ {gas}$$

that is to say :

$$h_{lg}h_m(\rho_{in}-\rho_{S})=h_c(T_{in}-T_S)+\lambda_l \left(\frac{\partial T}{\partial y}\right)_{y=h}$$

  • Nusselt number & Sherwood number

After all that, the Nusselt number and the Sherwood number are interesting to be studied. These numbers compare the convective effect and the diffusion effect.

The local Nusselt number compares the thermal convection flux with the thermal diffusion flux and it is defined as :

$$ Nu=\left(\frac{\partial \theta}{\partial Y}\right)_{Y=h/L}$$

Then we can deduce the global Nusselt Number as :



To finish this first part, the Sherwood number compares mass convective flux with molecular diffusion.

$Sh=\frac{h_m*D_h}{D}$ with $h_m=\frac{\dot{m}}{\rho(\omega_{in}-\omega_{S})}$

The total flux of vapour mass is $\dot{m}=\rho_l v_l=\rho_S v_S=\omega_S\rho v_S-\rho D(\frac{\partial \omega}{\partial y})_{y=h}$.



This work was made according to the thesis "Transfert couplé de chaleur et de masse par convection mixte avec changement de phase dans un canal" fulfilled by Othmane OULAID.

Diffusion dominant position

Diffusion dominant position                


In this part, the context is the same as in the precedent part  but the air charged in vapour is fixed. Then the molecular diffusion is the only mechanism which transports the vapour into the air.

For this, we will simplify the developed equations in the previous section depending on this specific case.


         1- Mass balance equation 

$$ \frac{d}{dt}(hDl\rho)=-\rho_S v_S S_{flux}$$

In this equation the liquid volume variation is equal to the vapour flux out going from the surface.


          2- Interface mass balance equation 

At the surface there is :

  • a production of vapour $\rho_S v_S S_{flux}$,
  • a transport by convection and diffusion (Fick's law)

The convection velocity is $v_c=(1-\omega)v_a+\omega v_v$ with $v_a$ air velocity and $v_v$ vapour velocity. The convection velocity in the case of fixed air is $v_c=\omega v_v$

At the surface $v_{cS}=\omega_S v_S$.


The output vapour is considered on the interface as : $\rho_S v_{cS} S=\rho_S \omega_S v_S S$.

The balance equation corresponds to "production=diffusion+evaporation",  then :

$$\rho_S v_S S=-D_S \rho_S (\frac{\partial \omega}{\partial y})_{y=h} S+\rho_S\omega_S v_S S$$

We then develop this expression and we obtain:

$$v_S=\frac{D_S (\frac{\partial \omega}{\partial y})_{y=h}}{\omega_S-1}$$

Afterwards the Spalding parameter $b=\frac{\omega}{\omega_S-1}$ then $v_S=D_S(\frac{db}{dy})_{y=h}  (8)$

After calculations, it's simplified to:

$$j_S=\rho_S \frac {D_S}{L-h}ln(1+(b_L-b_S))$$

and becomes:

$$j_S=\rho_S \frac {D_S}{L-h}ln(1+\frac{(\omega_S-\omega_L)}{1-\omega_S})$$


         3- Conclusion

We can observe on this evaporation model that the flux depends on the boundary conditions and thermodynamic balance.

The model is available for a mono-constituent. For a particular case we can deduce the thermodynamic balance parameters $T_{Sat}=T_{wall}$ and $\rho_{Sat}$ with tables.

The comparison inside the logarithm from the last formula ( the sign of the difference between $b_L$ and $b_S$)  gives us the sign of the flux:

  • If the mass concentration is more important far from the interface than the mass concentration at the interface, then flux is negative : there is condensation.
  • Else if $b_S>b_L$ there is evaporation. Else $b_L=b_S$ the system is in total balance state.

Moreover, we can deduce a characteristic time of the problem:


To have an order of magnitude, here there's a numerical application in the water case : 

$T_w=300K, \rho_S=2.38*10^{-2}kg/m^3, L=50 cm, D_S=5*10^{-5} m^2/s$

we obtain $\tau=60$ hours.

Aeration dominant position

Aeration dominant position                


        Evaporation leaded by the heat science

If the air velocity is very important, we can simplify the problem, it may still be necessary to consider that the convection transport is very effective. Then, in that case, it is the conduction inside the fluid which could master the mechanism. In other words, the convection takes the vapour faster than the vapour is produced.

At the interface, we can made a thermal flux balance. As we know, the interface is in a balance state :

$$ 0=q_s+q_L$$

$$0=-\lambda_{liq} S_{flux} \frac{\partial T}{\partial y}\mid_{y=h}+\dot{m_{vap}}h_{fg}$$

with that we can write

$\dot{m}=\rho_l v_l S_{flux}=\rho_s v_s S_{flux} $

and $v_s$ is a boundary condition for the air at the interface liquid-vapour. This is brought by the mass conservation.

Now we can write $v_l=-\frac{d h}{d t}$

with $h$ the effective length of the liquid film.

We have $q_s=q_l$

afterwards : $$\frac{\lambda_l \Delta T}{h}=\rho_l \frac{d h}{d t} h_{fg}$$

We can deduce : $$h^2(t)={h_0}^2-\frac{2\lambda_l \Delta T}{\rho_l h_{fg}} t$$

with a characteristic time $$ \tau=\frac{h_0 \rho_l h{fg}}{2 \lambda_l \Delta T}$$


           Numerical Application :

for $\Delta=T_{wall}-T_{sat}=1 K, h_0=100\mu m, h_{fg}=2400 kJ/kg, \lambda_l=0,67 W/Km, \rho_l=1000 kg/m^3$.

We obtain:

$$\tau=17,9 s$$






Fuel evaporation in the container

Fuel evaporation                                       

in the container


Evaporation rate model

Evaporation rate model                      


  In the previous parts, the general equations and the global context have been presented, as well as two examples of application in particular situations :

  • pure diffusion 
  • pure aeration

In this part, the industrial case is studied. This work requires both a theoretical and a numerical study in order to identify the parameters which give the mass concentration profile throughout the container. The aerodynamic equations is intrinsically linked with thermal and mass coefficient. In further, the link between the thermal and the mass mechanism is developed. As referencing with the first part, the Sherwood and the Nusselt number appear in the equations in their dimension form $h_c$ and $h_m$.


       1- Mass conservation equation

$$ab\rho v= ab\rho_s v_s   (1)$$


     2- Spalding mass transport equation

$$\rho v \frac{\partial b}{\partial y}=\frac{\partial}{\partial y}\left(\rho D\frac{\partial b}{\partial y}\right)   (2)$$

with the spalding number $b=\frac{\omega}{\omega_S-1}$


         3- Heat transport equation 

$$\rho v \frac{\partial  C_p T}{\partial y}=\frac{\partial}{\partial y}\left(\lambda\frac{\partial   T}{\partial y}\right)   (3)$$

$$\lambda_s\left(\frac{d T}{d y}\right)_s=h_c(T_\infty-T_s)$$

\[h_c(T_\infty-T_s)=\rho_s v_s L_v\]


          4-  Boundary conditions

The boundary conditions are expressed as follows:

  • $b=b_{\infty}$ for $y=\delta_\omega$
  • $b=b_{s}$ for $y=h$
  • $v_s= D_s \left(\frac{d b}{dy}\right)_s$
  • $T=T_\infty$ for $y=\delta_T$
  • $T=T_s$ for $y=h$


Therefore, we can simplify the previous expressions some simplifications, the result is :

$W_v(x)=b\rho_sD_s ln(b_{\infty}-b_s+1)\frac{1}{\delta_{\omega}(x)-\delta}$

$$W_v(x)=\frac{b\lambda}{C_p} ln\left( 1+\frac{W_v(x)}{b}\frac{C_p}{h_c(x)}\right)\frac{1}{\delta_T(x)-\delta}$$


The result is close to the droplet result of J.Reveillon ( Rouen University). The second equation is entirely dependent on thermal conditions. That is why we can have an expression of the mass flux of vapour with the values of the thermal boundary layer $\delta_T$ and the thermal convection coefficient $h_c$. Then we can deduce that the mass flux influences the mass boundary layer.


Parametric study

Parametric study                                    


In this section, we will use the expressions of different parameters that we've found to compute these parameters in multiple flight's cases.


                  1- Physical quantities

We need the values of the diffusion coefficient in the air and the saturation vapour pressure to can compute our parameters. 

The following table summarizes the values of these physical quantities in different flight's cases:




                 2- Parametric study

Now we can calculate the non dimensional numbers involved in this problem. We summarize all the computed values in the next table for different flight's cases: