Convective evaporation of vertical films
François Boulogne Benjamin Dollet
cond-mat.soft physics.flu-dyn




[
    \begin@twocolumnfalse
        Mag2Pol: A program for the analysis of spherical neutron polarimetry, flipping ratio and integrated intensity data  
    N. Qureshi
    December 30, 2023
====================================================================================================================


        
            Motivated by the evaporation of soap films, which has a significant effect on their lifetime, we performed an experimental study on the evaporation of vertical surfaces with model systems based on hydrogels.
            From the analogy between heat and mass transfer, we adopt a model describing the natural convection in the gas phase due to a density contrast between dry and saturated air.
            Our measurements show a good agreement with this model, both in terms of scaling law with the Grashof number and in terms of order of magnitude.
            We discuss the corrections to take into account, notably the contribution of edge effects, which have a small but visible contribution when lateral and bottom surface areas are not negligible compared to the main evaporating surface area.
        
    \end@twocolumnfalse
]









§ INTRODUCTION


Giant soap films are spectacular structures since they are uncommon in nature <cit.>.
Indeed, soap films are well-known to be fragile and short standing.
To  maintain such large films artificially, the soap solution is constantly injected at the top to balance the effects of drainage and evaporation <cit.>.
These are common in scientific outreach events, as shown in Fig. <ref> where they attract interest from the public due to its gigantic size.
Film evaporation depends on the environmental conditions such as the weather or air currents in the room. Thus, daily fluctuations of these parameters influence the success of the demonstration.







Several scientific studies report effects of evaporation on the soap film stability <cit.>.
Yaminsky et al. evidenced that impurities influence the lifetime of soap film only if placed in an unsaturated atmosphere <cit.>.
In addition, Li et al. pointed out that the nonuniformity of the evaporative flux induces Marangoni effects, which trigger the film burst <cit.>.
Moreover, chemistry can have a significant role on the evaporation kinetics as dense monolayers at interfaces can reduce significantly the evaporation rate <cit.>.
Therefore, the complexity in studying evaporation on soap films arises from the interplay of chemistry, hydrodynamics of liquid films with interfacial effects, and the dynamics of evaporation.


In 1887, Maxwell described the evaporation of a suspended spherical drop as a diffusive process of water in the gas phase.
While this description is successful in many situations such as the evaporation of small sessile drops <cit.>, natural convection has been recently revealed in experimental studies with various liquids <cit.>.

Indeed, water vapor is less dense than dry air.
Above a critical size of the evaporative surface, a convective gas flow takes place, which affects the vapor concentration field.
Typically, the critical size above which convection can be expected is about a few millimeters in room conditions for water.
Models have been developed to describe particular configurations such as the evaporation of pendant drops <cit.> and horizontal circular disks <cit.>. Similar convection effects occur in the dissolution of sessile drops <cit.>.




In this paper, we aim to describe the convective evaporation of vertical films.
First, based on an analogy with heat transfer on flat vertical surfaces, we present the theory describing the natural convective flow in such geometry with a prefactor determined numerically.
Because the evaporation of soap films is particularly challenging to study due to their limited lifetime and to the complex interplay between chemical and hydrodynamical effects, we make the choice to perform experiments on a model system, which consists of hydrogels.
We compare our measurements to the theoretical prediction and we concluded on the significance for future soap films studies.





















§ EXPERIMENTAL METHODS






 §.§ Hydrogels





Hydrogels are synthesized in polydimethyl acrylamid (PDMA) <cit.>.
The quantity of the monomer (N,N-dimethylacrylamide denoted DMA) is set by the mass ratio m_ DMA / (m_ DMA + m_w) = 0.1, where m_w is the mass of pure water.
The crosslinker (N,N'-methylene-bis-acrylamide denoted MBA) is added at a molar concentration ratio of [MBA]/[DMA] = 2× 10^-2.
Finally, a quantity of initiators (potassium persulfate and N,N,N',N'-tetramethylethylenediamine) is set to a molar ratio of 1% of the monomer quantity. All the chemical are purchased from Sigma-Aldrich (France).

The solution is poured in cells made of two glass slides separated by a rubber of about 3 mm in thickness, and hold with drawing clips.
The cell is sealed with a piece of parafilm and left in room conditions for 24 h to complete the reaction.
Then, the gel is carefully unmolded and placed in a tank filled with pure water for three days to swell the hydrogel.

The hydrogel is cut with a razor blade in a rectangular shape of the desired size.
The swollen hydrogel has a thickness b=4.5±0.5 mm.



 §.§ Evaporative flux measurements






The full setup is represented in figure <ref>.
Our experiments of controlled evaporation are performed in a box made in polycarbonate with dimensions of 50×50×50 cm^3.
The total evaporative flux of the slab of hydrogel is measured from the time evolution of its weight.
A precision scale (Ohaus Pioneer 210 g) with a precision of 0.1 mg is used with the hook available below the apparatus.
To hang the hydrogel on this hook, we made a clip with two polycarbonate rectangular plates covered with sand papers.
These two parts are placed on the top of the hydrogel and tightened with nylon screws.
This clip is then attached to the hook on the scale with a thin wire.


The scale is interfaced with a Python code using the pyserial library to record the weight every second.
The humidity is regulated with a PID controller based on an Arduino Uno and a humidity sensor (Honeywell HIH-4021-003) positioned far from the evaporating surface.
Dry air is produced by circulating ambient air with an air pump (Tetra APS 300) in a container filled with desiccant made of anhydrous calcium sulfate (Drierite).
Moist air is obtained by bubbling air in water.
The relative humidity is set to R_H=50% in all of our experiments.



 §.§ Parameters and order of magnitudes



Buoyancy is the driving force of the free convection flow.
It is therefore necessary to estimate the density contrast between the ambient air and the air saturated in water vapor.
The air density ρ at a pressure P_0, a temperature T and a relative humidity R_H is given by <cit.>

    ρ = 1/z_mP_0/ R T M_d [1 - f R_H ( 1 - M_w/M_d)  P_s/P_0],

where z_m and f are the compressibility and enhancement factors, R the ideal gas coefficient, M_d and M_w the molar density of dry and saturated air, and P_s the saturated pressure.
At room temperature, we have f = 1, z_m = 1 <cit.> and a saturated pressure P_s = 2.3 kPa <cit.>.
Thus, the density contrast at a relative humidity R_H=0.5 is |ρ_s - ρ_∞| / ρ_∞ = 5 × 10^-3.
In addition, the kinematic viscosity is ν = 1.5× 10^-5 m^2/s and the diffusion coefficient of water vapor in air is D = 1 × 10^-5 m^2/s.
To justify the Boussinesq approximation used in the following Section <ref>, it is useful to estimate the typical relative variations of the viscosity and diffusivity.
According to Tsilingiris <cit.>, the relative dynamical viscosity variation is 1.2× 10^-2 and the relative diffusivity variation is 4× 10^-3, which are both weak.


The Grashof number Gr is defined as the ratio between the buoyant forces driving convection and the viscous forces damping this flow, i.e.

 
    Gr =  g/ν^2ρ_∞ - ρ_s/ρ_∞ L^3,

 
 where L is the characteristic vertical length scale of the surface of evaporation.
This Grashof number will be derived from the equations of hydrodynamics in Section <ref>.

This number predicts the dynamics of evaporation <cit.>.
For small Grashof numbers, the evaporation is diffusive such that the ambient air can be considered as quiescent.
However, for large Grashof numbers, convection takes place.
The transition can be written as a critical length scale L^⋆ defined for Gr=1, i.e.


    L^⋆ = (g/ν^2ρ_∞ - ρ_s/ρ_∞)^-1/3.


For water in room conditions with a relative humidity of R_H=50%, L^⋆ is typically 2 mm.
Thus, for vertical film of few centimeters or more in height, convection is expected to produce an effect on the evaporative flux.

In the next Section, we recall the equations describing the natural convective flow near a vertical surface.
From these equations, we will derive the evaporative flux that can be compared with experiments.



§ THEORY




We consider the convective evaporation of a vertical rectangular surface of width w, height h and thickness b.
The model focuses on the flow due to the largest surfaces of area wh, referred thereafter as plates.
As a first approximation, the flow is assumed bidimensional, independent on z.
The effect of the two vertical edges of area hb will be discussed afterwards.
The velocity field in the gas phase is denoted u = u  e_x + v  e_y.
All the notations are presented in Fig. <ref>.

Equations of mass transfer due to evaporation are analogous to heat transfer.
Therefore, we closely follow the derivations of Schmidt and Beckmann on heat transfer published in 1930 <cit.>, also detailed later by Ostrach <cit.>.
The reader may also refer to textbooks such as <cit.>.
We choose to rewrite in details these derivations with notations adapted to mass transfer for the convenience of the reader.




 §.§ Formulation of the problem


We denote ρ_s the density of saturated air corresponding to a vapor concentration c_s.
We assume that the gas phase obeys the ideal gas law so that the density ρ varies linearly with the vapor concentration c,


    ρ(c) = ρ_0 - Δρc/c_s,


with Δρ = ρ_0 - ρ_s and ρ_0 the density of dry air.

In Cartesian coordinates, the continuity equation is


    ∂ u/∂ x + ∂ v/∂ y = 0.

The momentum equation is written within the Boussinesq approximation, which considers that the gas density is constant and equal to the density ρ_∞ far from the surface, except in the term driving the convection by a buoyancy effect, namely ρ(c) g.
Assuming a stationary flow, we have


    ρ_∞(u ∂ u/∂ x + v ∂ u/∂ y) = - ρ(c)g + μ( ∂^2 u/∂ x^2 + ∂^2 u/∂ y^2),


where μ is the gas viscosity.
The term ∂^2 u/∂ x^2 is negligible compared to ∂^2 u/∂ y^2 due to the boundary layer geometry present in the limit of large Grashof numbers.


The volatile solvent concentration field follows a diffusion-convection equation


    u ∂ c/∂ x + v ∂ c/∂ y =  D∂^2 c/∂ y^2,


where D is the diffusion coefficient of the vapor in the gas phase.



Equations (<ref>-<ref>) must satisfy five boundary conditions.
Two are related to the concentration fields with a saturated vapor concentration near the surface and a concentration c_∞ far from it.
In addition, three hydrodynamical conditions are set, with a no-slip condition on the surface for both components u and v, and a vanishing vertical component u far from the surface.
All these conditions write

    c(x,y=0)   =c_s,
    
            c(x,y→∞)   =c_∞,
            
            u(y=0)   =0,
    
            v(y=0)   =0,
    
            u(y→∞)   =0.








 §.§ Non-dimensionalized equations


To solve the set of equations (<ref>, <ref>, <ref>) with the boundary conditions (<ref>), we non-dimensionalize the system.
First, we define the dimensionless concentration c̃ = (c - c_∞)/(c_s - c_∞) and the streamfunction Ψ as

    u = ∂Ψ/∂ y,      v = - ∂Ψ/∂ x.

Therefore, equation (<ref>) becomes

    ∂Ψ/∂ y∂c̃/∂ x + ∂Ψ/∂ x∂c̃/∂ y =  D∂^2 c̃/∂ y^2,

and the momentum equation (<ref>) is rewritten as

    ∂Ψ/∂ y∂^2 Ψ/∂ x ∂ y  -
        ∂Ψ/∂ x∂^2 Ψ/∂  y^2  =  - ρ_∞ - ρ_s/ρ_∞ g c̃ +  ν∂^3 Ψ/∂ y^3,

where ν = μ/ρ_∞ and with the boundary conditions


    c̃(x,y=0)   =1,
    c̃(x,y→∞)   =0,
            ∂Ψ/∂ y(x, y=0)   =0,
    ∂Ψ/∂ x(x, y=0)   =0,
    ∂Ψ/∂ y(x, y→∞)   =0.


We introduce the self-similar variable

    ξ = ( g/4ν^2ρ_∞ - ρ_s/ρ_∞)^1/4y/x^1/4,

that we can also write ξ = 4^-1/4 Gr_x^1/4 y/x
 with the local Grashof number defined as
 
    Gr_x =  g/ν^2ρ_∞ - ρ_s/ρ_∞ x^3.


For x= h = L^⋆, we recover the definition of the global Grashof number, as defined by Eq. (<ref>).


We rewrite the concentration field c̃(x,y) = ĉ(ξ) and the streamfunction Ψ = 4^3/4ν  Gr_x^1/4 f(ξ) where f is the dimensionless streamfunction.
Equations (<ref>) and (<ref>) become respectively


    
    ĉ” + 3   Sc  f ĉ' = 0,
    
            f”' + 3 f f” - 2 f'^2 + ĉ = 0,



where we introduced the Schmidt number Sc= ν/ D.
The related boundary conditions (<ref>) become


    ĉ(ξ=0)=1,
    ĉ(ξ→∞)=0,
            
            f'(ξ=0)   =0,
    
            f(ξ=0)   =0,
    
            f'(ξ→∞)   =0.





 §.§ Resolution





We solve the system of ordinary differential equations (<ref>) with the conditions (<ref>).
As two of these conditions are specified at ξ→∞, we use a shooting method by guessing the values of ĉ_0' = ĉ'(ξ=0) = ( dĉ/ dξ)|_ξ=0 and
 f̂”(ξ=0).
 A fast convergence toward a solution satisfying ĉ(ξ→∞)=0 and f'(ξ→∞)=0 is obtained with the method described by Nachtsheim and Swigert <cit.>.
Briefly, this method consists in writing and solving differential equations on the shooting parameters to determine the corrections on the starting parameters.
We implemented this algorithm with the Python language, the libraries numpy and scipy <cit.>.
Our implementation has been successfully validated with the data available in <cit.> and our source code is available in Supplementary Materials.

The resolution is performed for Sc = 0.75 and we represent the dimensionless vapor concentration c̃ and the velocity field u in Fig. <ref>.
The thickness of the boundary layer along the plate grows as x^1/4 as indicated by the self-similar variable ξ defined in Eq. (<ref>) and the flow is oriented mainly upward with a horizontal component toward the plate.
In addition, we obtain ĉ_0'≈ -0.51 for our conditions.
The value of this coefficient is used hereafter in Sec. <ref> to estimate the evaporative flux.





 §.§ Evaporative flux
 

As the concentration field surrounding the plate is now determined, we can calculate the evaporative flux.
The local flux is defined as

    j = -  D.∂ c/∂ y|_y=0.


With the definition of the dimensionless concentration ĉ, we obtain, on the surface of the plate of area wh, the local convective flux


    j_ plate = -  D (c_s - c_∞) ĉ_0'  Gr_x^1/4/x.


Note that the numerics is only necessary to determine the prefactor ĉ_0' of Eq. (<ref>) and that the scaling is obtained analytically.


The total evaporative flux due to convection is Q_ plate = ∫ j_ plate  d S, which leads after integration to


    Q_ plate = 2 w  D (c_s - c_∞) (- 4/3ĉ_0')  Gr^1/4,

where Gr is defined by equation (<ref>) for the characteristic length scale L = h.

In the next Section, we discuss the validity of this prediction that we compare to our experimental results obtained by drying vertical hydrogels of different sizes and aspect ratios.





§ DISCUSSION




In Fig. <ref>, we present data on the time variation of the relative weight Δ m for some of the hydrogels used in our experiments.
To obtain the total flux Q, we fit our measurements with a linear law Δ m(t) = -Qt.
This procedure is repeated for hydrogels of different sizes.
The linear fits indicate that the evaporative flux remains constant over the duration of the experiment.
Thus, the variations of the polymer concentration in the hydrogel is negligible enough to ensure a constant saturated vapor concentration.
In a previous work, we also checked that this saturated vapor concentration can be assimilated to the one of pure water <cit.>.





In Fig. <ref>(a), we show experimental measurements of the evaporative flux Q for two hydrogel widths, w=30 mm and w=80 mm and heights ranging between 15 mm and 210 mm.
We plot these fluxes as a function of the Grashof number Gr, which is varied over three orders of magnitude.
As expected, the evaporative flux is systematically larger for wider gels at a given Grashof number, i.e. at a given height h.
For both widths, the measurements show a good agreement with the prediction in Gr^1/4 obtained from Eq. (<ref>), which means that the power law on the Grashof number is well captured by the model.








In Fig. <ref>(b), we compare our measurements with the prediction given by Eq. (<ref>).
By fitting all our measurements Q to the prediction Q_ conv≃ Q_ plate, we obtain Q≈ 1.31  Q_ conv, which is particularly reasonable at this level of theoretical description and for the precision of our experiments.
Nevertheless, in this representation, we can notice for the smallest hydrogel widths that Q_ conv underestimate the measured flux.
We interpret this deviation by the absence of consideration of edge effects in the theoretical description that become visible for our smallest samples.
Edge effects may arise from the two side edges and the bottom edge.
Also, tip effects at the corners introduce a divergence in the local flux, which is not taken into account.
The major difficulty in describing all these effects is that they are not simply additional corrections because edges and tips can have a mutual influence.

Despite these difficulties, we attempt to write some corrections to the prediction given by Eq. (<ref>).
First, we consider that the local flux j_ conv obtained by considering the surfaces of area wh, also holds for the two vertical surfaces of area bh.
Thus, the integral is calculated over a total area 2(w + b)h and the flux due to convection Q_ conv≃ Q_ plate + Q_ side is

    Q_ conv≃ 2 (w+b)  D (c_s - c_∞) (- 4/3ĉ_0')  Gr^1/4.

In Fig. <ref>(b), we observe that this correction reduces the deviation for w=30 mm and provide an estimate closer to the experimental results, Q≈ 1.24  Q_ conv.

A second correction that could be considered is related to the bottom edge.
In our previous study <cit.>, we argued that the boundary layer approximation breaks down at the leading edge of the plate <cit.> over a length scale L^⋆.
This correction must be also related to the evaporative flux over the bottom side of area bw.
Applying these small corrections lead to a better agreement with the measurements but, presently, their validities are particularly difficult to discuss.
Thus, we prefer to do not overinterpret our results as the main trend is well captured.
Further developments, especially from numerical simulations, could help to refine the model on these aspects.





§ CONCLUSION


Motivated by the question of the soap film bursting, we identified that one of the key aspects is related to the evaporation of the soap film.
As shown by the recent literature, evaporation can induce a convective flow in the atmosphere, which depends on the geometry and influences the evaporation kinetics.
Therefore, to understand the effect of evaporation of soap films, it is crucial to develop such models.

In this paper, we followed derivations established for heat transfer to calculate the evaporative flux on a vertical plate.
To verify the validity of this description, we performed model experiments with slabs of hydrogels of various widths and heights.
From the time evolution of the weight, we measured the total evaporation rate that we compared to the theoretical prediction.
We show that our experiments are in good agreement with this model, following the power law Gr^1/4.
Although that the model does not consider edge and tip effects, it provides a reasonable estimation of the prefactor.
We indicate that corrections could lead to a better agreement, especially for small films where edge effects are more significant.
We suggest that a more complete description including these effects would require numerical simulations.

Nevertheless, our results shows that the model presented in this paper provides a good description of the evaporation of vertical rectangular films, provided that the boundary layer remains laminar all along the surface.
In the near future, we would like to investigate the evaporation of soap films in light with this theory to rationalize the mechanisms that relates the bursting of soap films to their evaporation.
We expect that different additional mechanisms can have a contribution to the evaporation of soap films, which are not present in the model system based on hydrogels.
First, the composition of the soap solution can lead to a vapor pressure different from water, eventually with an effect of the increasing concentration of non-volatile molecules. 
The evaporation of thin films may generate temperature gradients, and therefore thermal Marangoni flows in the soap film, in addition to the solutal Marangoni flows, the gravitational drainage, and the marginal regeneration in free standing films.
The film thinning due to evaporation and liquid flow leads the appearance of the so-called black film where van der Waals attraction forces balance the double-layer repulsion forces, and therefore modifies the chemical potential of the liquid <cit.>.




§ ACKNOWLEDGMENTS

We are particularly grateful to Emmanuelle Rio for the initial discussions that motivated this study.
We also thank Lorène Champougny, Jonas Miguet, Christophe Poulard and Frédéric Restagno for discussions.
F.B. thanks Mélanie Decraene for her assistance.
This study has been carried out with a funding support by the ANR (ANR-11-BS04-0030-WAFPI project).




unsrt


