Direct Numerical Simulation of interfacial mass transfer into falling films

International Journal of Heat and Mass Transfer 69 (2014) 343–357

Contents lists available at ScienceDirect

International Journal of Heat and Mass Transfer journal homepage: www.elsevier.com/locate/ijhmt

Direct Numerical Simulation of interfacial mass transfer into falling films Christoph Albert b, Holger Marschall a, Dieter Bothe a,⇑ a b

Mathematical Modeling and Analysis, Center of Smart Interfaces, Technische Universität Darmstadt, Alarich-Weiss-Str. 10, 64287 Darmstadt, Germany DFG International Research Training Group 1529, Mathematical Fluid Dynamics, Technische Universität Darmstadt, Schlossgartenstr. 7, 64289 Darmstadt, Germany

a r t i c l e

i n f o

Article history: Received 23 July 2013 Received in revised form 11 October 2013 Accepted 12 October 2013 Available online 9 November 2013 Keywords: Volume-of-Fluid method Falling film Interfacial mass transfer

a b s t r a c t This contribution is concerned with interfacial mass transfer into falling films, for which detailed insight is achieved by means of Direct Numerical Simulation employing the Volume-of-Fluid method. As for the numerical model for mass transfer across the fluid interface of falling films, we employ the two-scalar approach. We examine the influence of wave regimes and flow patterns on local Sherwood numbers and identify distinct mass transfer enhancement mechanisms taking effect in specific wave regimes. Ó 2013 Elsevier Ltd. All rights reserved.

1. Introduction Both in Process and Chemical Engineering, heat and mass transfer apparatus are in widespread use in which a continuous liquid phase is processed as a continuous thin film falling down flat or structured walls due to gravity. Henceforth we shall refer to such configurations as falling films. Falling films are used in a variety of industrial applications where large rates of heat or mass transfer are desired. Examples are cooling devices [23], gas absorbers [9,43], desalination units [29,45] and falling film microreactors [21,41]. In these falling film apparatus, interfacial waves are known to play a central role for heat and mass transfer enhancement and, hence, process intensification. It has been shown experimentally that the presence of surface waves has an intensifying effect on the observed heat transfer rates from the wall into the film [15], and that the wave effect is even larger for the transport of a chemical species from the gas into the liquid [3]. The observed transfer rates of a chemical species into the liquid depend to a large extent on the flow structure inside the nonlinear waves. The flow structure, in turn, depends on the wavelength of the perturbation from which it develops. As a consequence of these non-trivial interdependencies, detailed experimental studies providing insights into the local sub-processes are not available. Moreover, most numerical simulations are based on simplified model assumptions, the effect of which for a correct prediction of the complex interplay ⇑ Corresponding author. Address: Petersenstr. 17, 64287 Darmstadt, Germany. Tel.: +49 6151 16 5363; fax: +49 6151 16 72022. E-mail address: [email protected] (D. Bothe). 0017-9310/$ - see front matter Ó 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.ijheatmasstransfer.2013.10.025

between hydrodynamics and interfacial mass transfer into falling films would require research in itself. The first numerical simulations of mass transfer into falling films were performed by Wasden and Dukler [39]. They used experimentally measured film thicknesses to compute the local flow fields, and solved a stationary species transport equation for a single wave. Boundary conditions in streamwise direction for the species concentration were computed by a solution of the concentration profile for a flat film. They found that the presence of surface waves results in velocity components normal to the interface, which in turn bring forth mass transfer enhancement. More recent simulations have been performed in [36,30], where the hydrodynamics was not modeled by the full Navier–Stokes equations, but by a simpler model derived by Shkadov [35], where the velocity profile is considered semi-parabolic at each streamwise direction. Both contributions look for an excitation frequency of the film that leads to a wave regime which is optimal in terms of achievable mass transfer rates. In [36], it is found that for low Reynolds numbers (Re), the variation of mass transfer rates exhibits two local maxima, whereas three peaks are found for higher Re. In [30], on the other hand, three peaks are found for low Reynolds numbers, and only two peaks are found for higher Re. In contrast to that, the experimental results in [44] show a single peak for low Reynolds numbers, and two for high Re. In [4], the hydrodynamics of the falling film was captured by Direct Numerical Simulation (DNS), while the species concentration was solved for stationary conditions, employing the method of Wasden and Dukler [39]. They find that recirculation in the large waves causes a decrease in absorption rates, whereas the

344

C. Albert et al. / International Journal of Heat and Mass Transfer 69 (2014) 343–357

Atmosphere

yˆ

ΩG (t) h(t x) Inlet

y 0

Outlet

h0

ΩL (t) 0 x

Wall

x ˆ

g Fig. 1. Computational domain with wavy falling film.

absorption rates increase in the capillary wave region. Here, the effect of the choice of artificial boundary conditions is unclear. Mass transfer simulations into falling films were also performed in [40]. The transient species transport equation was solved over the whole length of the film. They introduced a hydrodynamic parameter (the gradient of the vertical fluctuating velocity at the interface), and find that this parameter can be used to characterize species transfer across the gas–liquid interface. Unfortunately, the range of investigated hydrodynamic parameters was quite small, and local mass transfer rates were only provided for a single simulation. In [20], a coupled heat/mass transfer problem is solved by a finite difference method. Results in the capillary wave regime are compared to results with a smooth film, and it is found that recirculation in the large waves results in enhanced transfer rates compared to the flat film. The present contribution is concerned with the complex interplay of hydrodynamics and interfacial mass transfer into falling films. Since falling films exhibit a transient hydrodynamic nature and complex interfacial structures, their detailed experimental examination is difficult. Hence, we utilize DNS based on the Volume-of-Fluid (VOF) method along with the two-scalar approach of Bothe and Fleckenstein [5] in order to study the interdependency of hydrodynamics of wavy laminar films and mass transfer taking place across the deforming fluid interface. The simple scenario of a falling film on a flat vertical wall shall serve as a prototype case for mass transfer across the interface within (laminar) wavy film regimes. The arising non-trivial flow structures and species concentration patterns are studied, and the pertinent mass transfer enhancement mechanisms are identified. In particular, we examine the influence of wave regimes and flow patterns on local Sherwood numbers. The remainder of this work is organized as follows: in Section 2 the governing equations for both hydrodynamics and interfacial mass transfer into falling films are presented. More precisely, we detail on the sharp interface model for hydrodynamics and interfacial mass transfer. Section 3 briefly describes the employed numerical methods. Section 4 sets forth the computational setup adopted both for validation (Section 5) and results (Section 6). Results are thoroughly discussed in Section 6. Finally, we provide conclusions regarding the obtained results in Section 7.

2. Sharp-Interface model We consider a continuous thin film of a viscous fluid falling down a vertical flat wall due to gravity. Moreover, we account for the absorption of a weakly soluble gas into the liquid phase. This allows us to neglect any feedback of species concentration

on the hydrodynamics. For this case, the overall resistance to interfacial mass transfer is considered to be inside the liquid film. The Squire Theorem [38], which states that the critical Reynolds number for two dimensional perturbations of parallel flow between fixed walls is lower than the critical Reynolds number for three dimensional perturbations, was extended to the case of falling liquid films in [42]. At the same time, it is known from experiments that the primary instability of the falling film is two-dimensional, see for example [26]. Therefore, we here consider the two-dimensional case only. This also allows for a sufficiently high spatial resolution, which is of significant importance for performing true DNS. ^Þ be the domain occupied by both phases xÞ  ð0; y Let X ¼ ð0; ^ (cf. Fig. 1). The x-axis is aligned with the wall and points in streamwise direction, while the y-axis points in cross-streamwise direction, away from the wall. It is assumed that the position of the interface, separating gas and liquid, is given by a height function h(t, x) over the x-axis. With this notation, the interface position is determined as

RðtÞ ¼ fðx; yÞ 2 X : y ¼ hðt; xÞg:

ð1Þ

The domain occupied by the liquid is given by

XL ðtÞ ¼ fðx; yÞ 2 X : y < hðt; xÞg;

ð2Þ

the domain occupied by the gas is given by

XG ðtÞ ¼ fðx; yÞ 2 X : y > hðt; xÞg:

ð3Þ

2.1. Hydrodynamics As for the mathematical model, we employ a sharp interface continuum mechanical model, i.e. the two-phase Navier–Stokes equations, assuming a sharp jump of material properties at the interface. The fluid interface is a surface of discontinuity separating both fluid phase regions XL and XG (X = XL(t) [ XG(t) [ R(t)), for which we consider two incompressible, Newtonian, immiscible fluids under isothermal conditions and without phase change. The equations for conservation of mass and (linear) momentum in the bulk XL [ XG are

qðut þ u  ruÞ  lDu þ rp ¼ qg r  u ¼ 0;

  1 ; 0

ðx; yÞ 2 X;

ð4Þ

ðx; yÞ 2 X;

ð5Þ

where u = (u, v) is the velocity, p the pressure, and g the gravitational acceleration. The phase dependent material parameters l and q denote the dynamic viscosity and mass density of the fluids. Table 1 Material parameters and dimensionless groups used for Direct Numerical Simulations. Parameter

Symbol

Value

Unit

Density liquid Viscosity liquid Density gas Viscosity gas Surface tension Gravitational acceleration Inlet concentration of oxygen, liquid Inlet concentration of oxygen, gas Henry constant Diffusivity of oxygen, gas Diffusivity of oxygen, liquid, Sc = 15 Diffusivity of oxygen, liquid, Sc = 30 Diffusivity of oxygen, liquid, Sc = 50 Diffusivity of oxygen, liquid, Sc = 100 Diffusivity of oxygen, liquid, Sc = 200 Diffusivity of oxygen, liquid, Sc = 400 Kapitza Number

qL lL qG lG r

9.98744e1 1.068e2 1.3138e3 2.0274e4 7.3638e1 9.81e2 2.37e7 4.186e5 2.70e2 2.0173e1 7.1290e04 3.5645e04 2.1387e04 1.0693e04 5.3467e05 2.6734e05 3.1497e + 003

g cm3 g s1 cm1 g cm3 g s1 cm1 g s2 cm s2 mol cm3 mol cm3

g cijin,L cijin,G H DG DL,15 DL,30 DL,50 DL,100 DL,200 DL,400 Ka

cm2 s1 cm2 s1 cm2 s1 cm2 s1 cm2 s1 cm2 s1 cm2 s1

345

C. Albert et al. / International Journal of Heat and Mass Transfer 69 (2014) 343–357 Table 2 Reynolds, Weber, Froude numbers, mean film thicknesses and excitation frequencies/Strouhal numbers used for Direct Numerical Simulations. Re

We

Fr

h0 in cm

x in Hz

Sr

15 23 31 53

49.79 24.42 14.85 6.08

2.24 2.77 3.21 4.20

1.738e2 2.004e2 2.213e2 2.646e2

15,20,25,30,35 15,20,25,30,35 20,25,30,35,40,45 0,20,30,40,50,60

0.0282, 0.0376, 0.0471, 0.0565, 0.0659 0.0245, 0.0326, 0.0408, 0.0490, 0.0571 0.0296, 0.0369, 0.0443, 0.0517, 0.0591, 0.0665 0, 0.0247, 0.0371, 0.0494, 0.0618, 0.0741

Assuming negligible surface mass density, constant surface tension and no-slip at the gas–liquid interface, the interfacial jump conditions are

sut ¼ 0 at y ¼ hðxÞ; spI  lðru þ ðruÞT Þt  nR ¼ rjnR

ð6Þ at y ¼ hðxÞ:

ð7Þ

Here r is the surface tension constant, j the interfacial curvature, and nR the outer unit normal. The notation

sWtðt; xÞ :¼ lim ðWðt; x þ hnR Þ  Wðt; x  hnR ÞÞ for x 2 RðtÞ h!0þ

ð8Þ

stands for the jump of a physical quantity W across the interface. For a detailed introduction to the continuum mechanical modeling of two-phase flows see, e.g., [19,37]. Since we do not consider phase change, the evolution equation for the interface reads as

ht þ uhx  v ¼ 0:

ð9Þ

The time-averaged mean value of h, as determined by the liquid flow rate at the inlet, is denoted as h0. In order to complement the mathematical problem formulation, boundary and initial conditions are to be provided. We prescribe a no-slip boundary condition at the wall and homogeneous Neumann boundary conditions at the upper boundary, i.e.

u ¼ 0 at y ¼ 0;

ð10Þ

@u ^: ¼ 0 at y ¼ y @y

ð11Þ

At an artificial outflow boundary, a homogeneous Neumann condition can cause numerical problems when negative streamwise velocities hit the outlet. In Albert et al. [1], this problem was met by means of a damping zone near the inlet, which exponentially damps deviations from the steady state. In the presented work, an atmospheric boundary condition is employed, that ensures

numerical stability while distorting only a smaller part of the flow field. It reads as

 @u ¼ 0; if limuðxÞ P 0; @x x¼^x x!^x

ð12Þ

ujx¼^x ¼ 0; if limuðxÞ < 0; x!^x  @ v  ¼ 0: @x x¼^x

ð13Þ ð14Þ

Since the pioneering work of Kapitza and Kapitza [22], it has been customary to apply a perturbation with a fixed frequency at the inlet, in order to control the film wave regime that develops. As in [17], we implement this by prescribing a time-dependent, non-homogeneous Dirichlet boundary condition at the inlet, i.e.

u ¼ ð1 þ  sinð2pxtÞÞuN at x ¼ 0; 0 6 y 6 h0 ;   L 3=2u ^: u ¼ ð1 þ  sinð2pxtÞÞ at x ¼ 0; h0 < y 6 y 0

ð15Þ

We call x the frequency of the excitation and  its magnitude. Here,

L ¼ u

1 h0

Z

2

h0

uN ðyÞ dy ¼

0

gh0 q0 ¼ 3mL mL

ð16Þ

is the mean velocity of the steady state solution. The parameter q0 is the time averaged volumetric flow rate of the liquid per unit width  L. The at the inlet. The free surface velocity amounts to 3/2u expression

uN ¼

qL g 2lL

ð2h0 y  y2 Þ 0

! at x ¼ 0;

0 6 y 6 h0

ð17Þ

is the well-known Nusselt profile, a steady state solution of the onephase Navier–Stokes equations with a free boundary, named after [27].

Fig. 2. Validation study – Comparison of solutions after 1.4 s and 1.5 s.

346

C. Albert et al. / International Journal of Heat and Mass Transfer 69 (2014) 343–357

The boundary conditions for (9) are:

h ¼ h0 at x ¼ 0;

ð18Þ ð19Þ

The initial conditions at t = 0 are u = uN and h = h0. Employing dimensional analysis, the hydrodynamics of a gravity-driven wavy falling film can be completely described by the following dimensionless groups [16]:

Ka ¼

3

 L h0 u

¼

mL

gh0 ; 3m2L

r 4=3 1=3 L L g

qm

or Fi ¼ Ka1=3 ;

Fr ¼

L u 1=2

ðgh0 Þ

;

ð20Þ

ð21Þ

;

h0 ^x

L u

ð25Þ

:

As for the mathematical modeling of interfacial species transfer, let us consider a two-phase system containing the transfer species i being transferred from the gaseous into the liquid phase. In addition to the already exposed hydrodynamic assumptions, species i is assumed to be dilute (passive scalar transport) and physico-chemically inert (no reaction, no adsorption onto the interface, i.e. no surface activity). Furthermore, we impose local chemical equilibrium at the interface such that Henry’s law holds. Together with conservation of species mass, these assumptions lead to the transport equation for species i, i.e.

@ t ci þ r  ðci u þ ji Þ ¼ 0; We ¼

3

1=3

Ka

Re5=3

;

Mo ¼

l ¼ K1 F ; qr 4 Lg 3 L

ð22Þ

ð23Þ

namely the film Reynolds, Kapitza (alternatively the equivalent Film, Weber, Morton or Fluidkennzahl) and the Froude number, respectively. Additionally, to account for the film geometry one introduces

C¼

xh0

2.2. Interfacial mass transfer

@h ¼ 0 at x ¼ ^x: @x

Re ¼

Sr ¼

ð24Þ

and in the case of forced waves the frequency of disturbance is characterized via the Strouhal number, being defined as

ðx; yÞ 2 X:

ð26Þ

The jump conditions in this simplified modeling framework (cf. [5]) are

sji t  nR ¼ 0 at y ¼ hðxÞ;

ð27Þ

ci jR;L ¼ Hci jR;G

ð28Þ

at y ¼ hðxÞ:

where we have assumed an ideally mixed gaseous phase for the latter condition. Herein, H denotes the Henry constant and cijR,L is the liquid-sided interfacial limit of the species concentration. Given an ideally diluted system (as considered here), the Fickian law of diffusion with constant diffusivities Di > 0 is suitable to model the molar flux density ji, i.e.,

ji ¼ Di rci :

Fig. 3. Validation – Species concentration averaged over the film height.

ð29Þ

C. Albert et al. / International Journal of Heat and Mass Transfer 69 (2014) 343–357

347

Fig. 4. Film height (dashed line) and streamwise velocity, relative to wave velocity, at the interface (solid line); Re = 15. Velocities lower than the lower dotted line indicate backflow (Section 6.1.3).

Fig. 5. Film height (dashed line) and streamwise velocity, relative to wave velocity, at the interface (solid line); Re = 23. Velocities larger than the upper dotted line indicate a vortex (Section 6.1.2), Velocities lower than the lower dotted line indicate backflow (Section 6.1.3).

Fig. 6. Film height (dashed line) and streamwise velocity, relative to wave velocity, at the interface (solid line); Re = 31. Velocities larger than the upper dotted line indicate a vortex (Section 6.1.2), Velocities lower than the lower dotted line indicate backflow (Section 6.1.3).

As for the boundary conditions at the inlet, fixed concentrations are prescribed according to

ci ¼ cijin;L at x ¼ 0;

y < h0 ;

ci ¼ cijin;G at x ¼ 0;

y > h0 :

ð30Þ

Here, cijin,L/G is the species concentration at the inlet inside the liquid and gas phase, respectively. At the other boundaries, homogeneous Neumann conditions are prescribed, i.e.

@c ¼ 0 at y ¼ 0; @y @c ^; ¼ 0 at y ¼ y @y

ð31Þ

@c ¼ 0 at x ¼ ^x: @x As initial condition, the concentration of species i is set to cijin,G in the gas and to cijin,L in the liquid phase.

348

C. Albert et al. / International Journal of Heat and Mass Transfer 69 (2014) 343–357

Fig. 7. Film height (dashed line) and streamwise velocity, relative to wave velocity, at the interface (solid line); Re = 53. Velocities larger than the upper dotted line indicate a vortex (Section 6.1.2), Velocities lower than the lower dotted line indicate backflow (Section 6.1.3).

Interfacial species transfer in falling films is commonly described by means of two dimensionless groups, namely by the Schmidt and Sherwood number,

Sc ¼

m Di

;

R ^x kint ¼

0

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 kloc ðxÞ 1 þ ð@ x hðxÞÞ dx q ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ; R ^x 2 1 þ ð@ x hðxÞÞ dx 0

ð34Þ

ð32Þ with the local species transfer coefficient being defined as

kh0 Sh ¼ ; Di

ð33Þ

where k = kint represents the overall (integral) species transfer coefficient, which is determined by integration over the whole interfacial area as

kloc ðxÞ :¼

  1 @cðxÞ cR ðxÞ  c1 :  @nR R Di

Fig. 8. Local backflow in the wave troughs, Re = 53 x = 50 Hz. Velocity is shown in the (fixed) reference frame of the wall.

Fig. 9. Film height (dashed line) and local Sherwood number (solid line) at Re = 15, Sc = 50.

ð35Þ

C. Albert et al. / International Journal of Heat and Mass Transfer 69 (2014) 343–357

349

Fig. 10. Film height (dashed line) and local Sherwood number (solid line) at Re = 23, Sc = 50.

Fig. 11. Film height (dashed line) and local Sherwood number (solid line) at Re = 31, Sc = 50.

Fig. 12. Film height (dashed line) and local Sherwood number (solid line) at Re = 53, Sc = 50.

3. Numerical solution method 3.1. Hydrodynamics For the Direct Numerical Simulation of the hydrodynamics of falling films the in-house code Free Surface 3D (FS3D), originally developed by Rieber [32], has been adapted to the investigation of falling liquid films in [1].

The underlying numerical approach is based on the Finite Volume Method (FVM). The position of the interface is captured by the Volume-of-Fluid (VOF) method [10,18]. Fluxes of volume fraction are determined geometrically with the Piecewise Linear Interface Calculation (PLIC) algorithm of Rider and Kothe [31]. Surface tension is treated by the Continuum-Surface-Force (CSF) model [7], in a so-called balanced force implementation with height functions, cf. [14].

350

C. Albert et al. / International Journal of Heat and Mass Transfer 69 (2014) 343–357

Table 3 pffiffiffiffiffi Variation of Shint = Sc (absorption rate) with frequency, Reynolds and Schmidt number.

x = 15 Hz

x = 20 Hz

x = 25 Hz

x = 30 Hz

x = 35 Hz

15 30 50

0.11341 0.16022 0.20006

0.11408 0.16054 0.19952

0.11398 0.16083 0.20017

0.11413 0.16045 0.19892

0.11386 0.15996 0.19657

23 23 23

15 30 50

0.18997 0.26046 0.31087

0.18910 0.26074 0.31105

0.19101 0.26203 0.31215

0.19087 0.26183 0.31098

0.1897 0.26118 0.31065

31 31 31

15 30 50

0.25511 0.34408 0.39921

0.25680 0.34508 0.40005

0.25668 0.34555 0.40075

0.25773 0.34614 0.40138

53 53 53

15 30 50

Re

Sc

15 15 15

x = 0 Hz

0.38829 0.43363 0.43761

0.39639 0.49920 0.54140

0.40144 0.50855 0.55430

For an in-depth description of the numerical approach, and a thorough validation of the solver with respect to the hydrodynamics of falling films, see [1].

For species transport we employ a two-field representation. For details, the interested reader is referred to [6] and, more recently, [5]. For the sake of completeness, we briefly recall central aspects of the approach.

x = 45 Hz

0.25783 0.34606 0.40035

0.25749 0.34402 0.39514

0.40430 0.51211 0.55763

x = 50 Hz

x = 60 Hz

0.40597 0.51574 0.56088

0.40698 0.51455 0.55758

The two-scalar approach to interfacial mass transfer accounts for the amount of species in each phase separately. For instance,

ci;L ðxÞ ¼ 3.2. Interfacial mass transfer

x = 40 Hz



ci ðxÞ if x 2 XL 0

otherwise;

ð36Þ

denotes the concentration of species i in the liquid phase. This exhibits the advantage of retaining information on the one-sided limit of species concentration within each interfacial cell. Moreover, the same methods adopted for the VOF-transport of the volumetric phase fraction field can be used for species advection.

Fig. 13. Streamlines (upper picture) and concentration fields for Re = 15, Sc = 50, x = 20 Hz. Streamlines are shown in the (moving) reference frame of the wave. Species concentration in 1e7 mol cm3.

C. Albert et al. / International Journal of Heat and Mass Transfer 69 (2014) 343–357

351

Fig. 14. Streamlines (upper picture) and concentration fields for Re = 15, Sc = 50, x = 30 Hz. Streamlines are shown in the (moving) reference frame of the wave. Species concentration in 1e7 mol cm3.

In order to accurately capture the steep concentration gradient at the interface, a subgrid-scale model is employed. Herein, we exploit a specific class of functions describing the species distribution near the interfacial boundary [6,5], namely

 cðdÞ ¼ cR ðxÞerfc

 d ; dd ðxÞ

ð37Þ

where cR denotes the interfacial species concentration and d the distance in normal direction to the interface. The local thickness (in d-direction) of the concentration boundary layer is denoted by dd(x). Note that this class of functions can be derived under the assumptions of a velocity that is constant and tangential to the interface, and a species diffusion which acts only in direction normal to the interface. With this approach, the concentration gradient at the interface reads as

 @c  2 1 ðxÞ ¼ pffiffiffiffi c ðxÞ; @nR R p dd ðxÞ R

in water is small, such that ScL  400. Since the thickness of the 1 diffusive boundary layer scales with Sc2 , interfacial mass transfer problems at high Schmidt numbers require a high spatial resolution. In order to evaluate the limit up to which the problem can be sufficiently resolved at acceptable computational costs, we performed numerical simulations with artificially increased diffusion constants corresponding to Schmidt numbers of 15, 30, 50, and 100. In Table 1, a list of all used material parameters can be found. Specified Reynolds numbers are accounted for by setting h0 according to (20). The size of the computational domain was fixed at ^ x ¼ 512h0 and yˆ = 4h0. We used quadratic computational cells with an edge length of h0/16, if not specified otherwise. The time evolution of the film was simulated for 1.5 s. For a list of simulated Reynolds numbers and excitation frequencies, see Table 2. The relative amplitude of the excitation was set to  = 0.05 (cf. (15)).

ð38Þ

which provides a nonlinear flux correction to save spatial resolution. Then, interfacial species flux can be computed from a concentration value c at a known distance to the interface and by fitting of the parameter dd(x) in (37) to the local conditions. The fitting is performed by Newton’s method. 4. Computational setup We consider a water film at 18°C, surrounded by an atmosphere of pure oxygen of the same temperature. The diffusivity of oxygen

5. Validation For a flat film flow, i.e. in the stationary case, the position of the interface and the velocity field is known, and we can numerically solve the resulting advection–diffusion equation for mass transfer inside the liquid phase with high numerical accuracy. We validate our method against this theoretical case, since sufficiently resolved measurements for species concentration in falling films are not available. In order to obtain reference solutions for validation at high resolutions, we ignore diffusion in streamwise direction, i.e.

352

C. Albert et al. / International Journal of Heat and Mass Transfer 69 (2014) 343–357

Fig. 15. Streamlines (upper picture) and concentration fields for Re = 53, Sc = 50, x = 30 Hz. Streamlines are shown in the (moving) reference frame of the wave. Species concentration in 1e7 mol cm3.

x-direction. Since the problem is strongly convection-dominated, the absolute error compared to the full diffusion problem is negligible. Accordingly, the system of governing equations simplifies to become

qL g ð2h0 y  y2 Þ@ x ci ¼ @ yy ci ; 2lL ci ¼ cijin;L ; @ci ¼ 0; @y

x ¼ 0; y ¼ 0;

ci ¼ Hcijin;G ;

y ¼ h0 :

0 < x < ^x;

0 < y < h0 ;

ð39Þ ð40Þ ð41Þ ð42Þ

This can be viewed as a one-dimensional heat equation, with diffusivity depending on the y-coordinate (cf. Eq. (39)). In this case, (40) provides the initial condition. The system is complemented by a homogeneous Neumann and a non-homogeneous Dirichlet condition – cf. Eqs. (41) and (42), respectively. Using second order Finite Difference Method (FDM) for discretization in y-direction, the system can be integrated in x-direction using a standard Ordinary Differential Equation (ODE) solver. Since the ODE-problem is still very stiff, the ode15s-solver of MatlabÒ (v. 2011 b) was employed to integrate the equation up to x ¼ ^x. The reference solution is calculated with a resolution of 4096 equidistant points in y-direction. This resolution is sufficient to achieve grid independent results. The step size in x-direction is adaptively adjusted by the ODE-solver.

For validation, we used FS3D to obtain the solution for a flat film at Re = 53. The onset of surface waves was avoided by setting  = 0. As for the required mesh density, this is to be considered the most challenging case, since the computational cell sizes are the larger, the higher the Reynolds numbers are chosen (cf. Section 4). We consider the species concentration averaged over the film height as the validation quantity, i.e. 

c ðt; xÞ ¼

1 hðt; xÞ

Z

hðt;xÞ

cðt; x; yÞ dy:

ð43Þ

0

Simulations at the three resolutions h0/8, h0/16, and h0/24 were performed. In order to prove that stationarity has been reached after  1.5 s, we consider the averaged concentrations c for different Schmidt numbers at 1.4 s and 1.5 s for medium resolution, cf. Fig. 2. The numerical solutions can be assumed stationary, since the corresponding lines for constant Sc cannot be distinguished anymore. The same holds for other spatial resolutions considered herein.  For validation, the averaged concentration values c are compared at different axial positions, cf. Fig. 3. Given the number of cases examined, the performed DNS reached the limit of acceptable computational costs at the intermediate resolution h0/16 (taking computation times of up to several weeks). Evidently, the numerical simulations at high Sc numbers show significant deviations from the reference solutions. With increasing Schmidt numbers the deviations become more pronounced. For this reason, we shall discuss only simulation results obtained for Sc up to 50 in the remainder, while disregarding results beyond this limit, as they

C. Albert et al. / International Journal of Heat and Mass Transfer 69 (2014) 343–357

353

Fig. 16. Streamlines (upper picture) and concentration fields for Re = 53, Sc = 50, x = 50 Hz. Streamlines are shown in the (moving) reference frame of the wave. Species concentration in 1e7 mol cm3.

cannot be considered reliable. Notably, however, DNS of interfacial mass transfer into falling films at these moderate Schmidt numbers already shows relevant features of high Schmidt number cases as will be shown in the remainder.

6. Results & discussion 6.1. Hydrodynamics For models of wavy falling films, aiming at the prediction of mass transfer enhancement due to surface waves, it is essential to identify wave regions that show distinct interfacial mass transfer characteristics, which in a second step have to be modeled, i.e. understood. The inner flow structure of the waves exhibits complex phenomena which are expected to have a significant influence on the observed species transfer rates. We are especially interested in the effect of two hydrodynamic phenomena that can be observed on falling films (cf. Fig. 4–7) within distinct wave regimes, depending on Reynolds number and frequency of excitation. Therefore, we first discuss the hydrodynamics of falling films by means of characteristic wave regimes, vortices and backflow. Then, the influence of these hydrodynamic phenomena on interfacial mass transfer into falling films is revealed. We identify distinct mass transfer enhancement mechanisms taking effect in the considered wave regimes. Our discussion is set out in terms of local Sherwood number profiles (i.e., local mass transfer rates) and the

typical concentration patterns arising thereof. In the results section, the interface position is obtained as a height function by summing up the volume fraction field in y-direction.

6.1.1. Observed wave regimes As depicted in Fig. 6 for x = 20 Hz, low excitation frequencies are known to lead to waves with long wavelength, that exhibit a typical roll wave/capillary wave structure. In this wave regime, there are large wave humps (roll waves), preceded by several smaller wave humps (capillary waves). An increase in excitation frequency first leads to a decrease of the number of capillary waves as shown in Fig. 6 for x = 30 Hz, and finally their disappearance (x = 40 Hz). This means that in this wave regime only large wave humps with an approximately sinusoidal form remain. Once waves in these wave regimes are fully developed, they stay unchanged over a considerable distance downstream. The velocity field inside these waves is stationary in a reference frame moving downstream with phase velocity; see, however, Section 6.2.3. In the experiment or in three-dimensional simulations, these time periodic waves are susceptible to spanwise perturbations and disintegrate into horseshoe-shaped solitary waves sufficiently far downstream, see [28]. Furthermore, there exists an intermediate wave regime, where genuine time-periodic waves do not appear (Fig. 6 for x = 35 Hz). In general, both wave amplitude and wave length decrease with excitation frequency. The notion of ‘‘low’’ and ‘‘high’’ frequencies depends on the Reynolds number. The higher the Reynolds number, the higher the necessary

354

C. Albert et al. / International Journal of Heat and Mass Transfer 69 (2014) 343–357

frequency for a certain wave regime. This can be seen by rewriting the Strouhal number defined in (25) in terms of Reynolds number as

Sr ¼ x

 2 3mL 3 13 Re : g

ð44Þ

In the performed simulations, capillary waves vanished for a common threshold, namely for Strouhal numbers exceeding a value of about 0.05. 6.1.2. Vortices It is known that for sufficiently large wave humps a vortex can occur in the coordinate system of the wave. This phenomenon has been studied in detail by Malamataris and Balakotaiah [24], and direct experimental evidence was presented in [2], in the slighty different, yet related case of liquid flow down the outside of an inclined tube. In Fig. 4–7, solid lines show the streamwise velocity component at the interface in the coordinate system of the moving wave, urel,R = uR  uwave. Here, uR = u|y=h and uwave is the phase velocity of the wave. The phase velocity of the wave uwave was determined from two downstream distances, which a wave peak traveled in Dt = 0.01 s. The presence or absence of a vortex is indicated by an auxiliary line drawn at uR  uwave = 0. For positive

values, the fluid at the interface moves faster than the wave, and a vortex is present in the wave. We observe vortices in the large wave humps only for films at sufficiently large Re. Since the ratio between maximal and residual film thickness increases with Re, this is in agreement with [8,25], who predict that this vortex occurs when the film thickness ratio exceeds a threshold value of about 2.5. Accordingly, vortices can not be observed in the crests of capillary waves.

6.1.3. Backflow Backflow denotes a localized region of the flow field where negative streamwise velocities, opposing the direction of gravity, occur near the interface in the (fixed) reference frame of the wall – cf. Fig. 8 for velocity arrows in the wave trough at Re = 53 and x = 50 Hz. In Fig. 4–7, a line is drawn at urel,R =  uwave. When interfacial velocities fall below this threshold, it means that uR is negative in the wall-fixed coordinate frame, indicating the occurrence of local backflow. Detailed numerical and experimental studies of the backflow phenomenon were performed in [12,11]. There, the occurrence of backflow was linked to the variation of interfacial curvature in the following way: Due to the small thickness of the film, the bulk

Fig. 17. Concentration fields for a wave at Re = 53, Sc = 50, x = 20 Hz on different downstream positions. Species concentration in 1e7 mol cm3.

C. Albert et al. / International Journal of Heat and Mass Transfer 69 (2014) 343–357

pressure in the liquid is almost entirely determined by the law of Young–Laplace,

pL  rj þ pG ; where the atmospheric pressure pG is nearly constant. Curvature in the capillary troughs is negative, while curvature in the wave humps is positive. The resulting sharp increase of pressure leads to a pressure gradient that exerts a force in negative x-direction, leading to a deceleration of the liquid. For low Re, backflow is only observed for waves of a low frequency, associated with the presence of capillary waves. This is in accordance with [12], where flows of low Reynolds numbers were investigated, and backflow was only found in the presence of capillary waves, for which the steepest pressure gradients occur. For large Re, however, backflow is also observed in the high frequency regime, where no capillary waves are present. Here, the gradient of interfacial curvature between two large wave humps is large enough to induce backflow. This has also been recognized recently in [13]. 6.2. Species transfer 6.2.1. Local Sherwood number profiles (local transfer rates) We define the local Sherwood number as

Shloc ðxÞ ¼ 

h0 @c ðxÞ; cR  cL;0 @nR

ð45Þ

i.e. as an appropriately normalized derivative in normal direction, where cR is the species concentration at the interface. The integral Sherwood number is defined as the line integral of the local Sherwood number over the free surface, averaged over the interface length, i.e.

R ^x Shint ¼

0

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 Shloc ðxÞ 1 þ ð@ x hðxÞÞ dx q ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi : R ^x 2 1 þ ð@ x hðxÞÞ dx 0

ð46Þ

In Fig. 9–12, local species transfer rates are shown over a representative part of the film, and are indicated by solid lines. For ease of

355

comparison, film thickness is shown by a dashed line. Integral Sherwood numbers for all simulated cases are given in Table 3. For all considered Reynolds numbers and excitation frequencies, a local maximum occurs in the wave troughs, irrespective of the type of adjacent wave humps. This observation is affirmed by measurements reported in [34,33] for heat transfer coefficients. The occurrence of a local maximum in the wave troughs is, in particular, independent of the presence of backflow. This behavior can be explained by the fact that urel,R exhibits a local minimum that coincides with the local minimum of film height. The quantity urel,R can be interpreted as a surface renewal rate, since it means that species concentration is removed with a very high velocity from the interface in the wave troughs. For roll waves without vortex as well as in capillary waves, a local minimum of species transfer rate coincides with the wave crest. These local minima are met by local maxima of urel,R, meaning that the velocity of fluid at that point is close to wave velocity, resulting in a smaller surface renewal rate. In larger wave humps – either due to a higher Reynolds number, or due to a smaller excitation frequency – the maximum of urel,R gradually approaches 0, until a vortex occurs. In this case two critical points appear – saddle points, such as stagnation or separation points. These points are denoted as A for the critical point upstream of the wave crest and B for the critical point at the downstream position (cf. streamline plots provided in Figs. 15 and 16). The local minimum of species transfer rate migrates from the wave crest to the downstream critical point B, where an affluent transport of species occurs from two positions: from the capillary wave region downstream, and from the part of the interface between both critical points. The upstream point A, on the other hand, experiences an efflux of liquid. If the local maximum in the nearest wave trough is sufficiently strong, this leads to the development of a plateau as can be seen in Fig. 12 for x = 50 Hz. However, if capillary waves are sufficiently close to this point, another local maximum can occur at this position; cf. Fig. 12 for x = 30/40 Hz.

Fig. 18. Streamlines for a wave at Re = 53, Sc = 50, x = 20 Hz on different downstream positions. Streamlines are shown in the (moving) reference frame of the wave.

356

C. Albert et al. / International Journal of Heat and Mass Transfer 69 (2014) 343–357

6.2.2. Concentration patterns In Fig. 13–16, species concentration fields and local flow fields in the moving reference frame are shown for four representative wave regimes: small and large Reynolds numbers, each with and without capillary waves. Since no relevant depletion of species in the gas phase takes place, only the concentration field inside the liquid is shown, resulting in step-like patterns visible around the interface. In Fig. 13 and 14, we see monotonous profiles in the concentration boundary layers; the concentration pattern is mainly determined by diffusive transport from the interface. No roll-up of species concentration can be observed in the roll waves, and species concentration decreases monotonously with distance from the interface. In the wave trough, the deceleration of liquid due to the adverse pressure gradient leads to a concentration boundary layer that is significantly thinner than the one in the large wave hump. In Fig. 14 we see a wave in the sinusoidal regime at a low Reynolds number, where no backflow occurs. The species distribution is very similar to the one in Fig. 13. In Fig. 15 and 16, waves with a vortex are considered. Species concentration is highly influenced by the presence of the vortex. A finger of high concentration starts at the downstream critical point and is rolled up, following the streamlines. In this regime, species concentration does not monotonously decrease with distance from the interface. In [24] it was pointed out that the vortex is not seen in the stationary frame of reference. This is true for the velocity field. In the presence of an advected scalar like species concentration, however, the vortex becomes clearly visible. 6.2.3. Time periodic behaviour inside the large wave humps We find that the assumption of stationarity inside the large wave humps is not valid anymore for large Reynolds numbers. The wave hump at Re = 53 and x = 20 Hz is the largest one considered in this work. In a time periodic manner, a second finger of high concentration develops and vanishes downstream of the large vortex. In Fig. 17, the concentration field inside a roll wave is shown at different times as it travels downstream. The time difference between the figures is 0.0125 s. In the first picture, only one finger of high concentration is visible, but there is a noticeable dent in the 8.25e7 mol cm3 isoline downstream of the vortex. In the second picture, the concentration boundary layer separates from the wave front, and in the remaining figures it can be seen how it is swept towards the vortex. In the last picture, after 0.05 s, the finger has vanished, and the concentration field is again very close to the situation at t = 0 s. This means that the second finger oscillates with the excitation frequency of 20 Hz. It was checked whether the phenomenon might be a purely numerical artifact by repeating the simulation with a higher resoh0 lution of 24 , without any significant effect on the results. Streamlines of the velocity field are shown in Fig. 18; they do, however, give no clear indication about the reason for this behaviour. The small visible oscillation of the vortex is caused by oscillations of the wave velocity. It might be that this behaviour is simply a remnant of the wave inception near the inlet. 7. Summary & conclusion The inner flow structure of waves in falling films exhibits complex phenomena which are expected to have a significant influence on the observed species transfer rates. Hence, for the development of physically sound models to predict interfacial mass transfer into wavy falling films, detailed knowledge about the influence of wave regimes and flow patterns on local mass transfer rates are of utmost relevance.

By means of Direct Numerical Simulations, we have identified distinct mass transfer enhancement mechanisms taking effect in specific wave regimes. Local Sherwood number profiles have been related to local hydrodynamic phenomena such as vortices and backflows occurring in falling films. The presence of the former has been found to be irrelevant with respect to mass transfer enhancement, while the occurrence of the latter determines, to a significant extent, both the local profiles of mass transfer rates into the liquid film and the concentration patterns developing inside the waves. This finding is believed to be of immediate relevance for the development of interfacial mass transfer models for wavy falling films. The gained insights from our Direct Numerical Simulations provide the ground for further research aiming at the disclosure and analysis of specific transport phenomena in wavy falling films – preferably along with results from high resolution experimental techniques accompanying the numerical results. In particular, the observed migration of the local species transfer minimum from the wave peak to a downstream stagnation point, and the possible development of a local maximum at the upstream stagnation point has been reported here for the first time. Moreover, the hitherto unknown phenomenon of a time-dependent concentration double finger in very large wave humps was described. Since the double finger is observed in one simulation only, however, a final conclusion regarding the hydrodynamic cause of this could not be drawn. Certainly, for a comprehensive understanding of interfacial mass transfer into falling films, these observations need further attention and should be experimentally confirmed.

References [1] C. Albert, H. Raach, D. Bothe, Influence of surface tension models on the hydrodynamics of wavy laminar falling films in Volume of Fluid-simulations, Int. J. Multiphase Flow 43 (2012) 66–71. [2] S. Alekseenko, V. Antipin, A. Bobylev, D. Markovich, Application of PIV to velocity measurements in a liquid film flowing down an inclined cylinder, Exp. Fluids 43 (2007) 197–207. [3] S.V. Alekseenko, V.E. Nakoryakov, B.G. Pokusaev, Wave effect on the transfer processes in liquid films, Chem. Eng. Commun. 141–142 (1996) 359–385. [4] S. Bo, X. Ma, H. Chen, Z. Lan, Numerical simulation on vapor absorption by wavy lithium bromide aqueous solution films, Heat Mass Transfer 47 (2011) 1611–1619. [5] D. Bothe, S. Fleckenstein, A Volume-of-Fluid-based method for mass transfer processes at fluid particles, Chem. Eng. Sci. 101 (2013) 283–302. [6] D. Bothe, M. Kröger, H. Warnecke, A VOF-based conservative method for the simulation of reactive mass transfer from rising bubbles, Fluid Dyn. Mater. Process. 7 (2011) 303–316. [7] J.U. Brackbill, D.B. Kothe, C. Zemach, A continuum method for modeling surface tension, J. Comput. Phys. 100 (1992) 335–354. [8] N. Brauner, D.M. Maron, Modeling of wavy flow in inclined thin films, Chem. Eng. Sci. 38 (1983) 775–788. [9] M.V. Dam, J.P. Corriou, N. Midoux, A.S. Lamine, C. Roizard, Modeling and measurement of sulfur dioxide absorption rate in a laminar falling film reactor, Chem. Eng. Sci. 54 (1999) 5311–5318. [10] R. DeBar, Fundamentals of the KRAKEN code, Lawrence Livermore Laboratory, UCIR-760, 1974. [11] G.F. Dietze, F. Al-Sibai, R. Kneer, Experimental study of flow separation in laminar falling liquid films, J. Fluid Mech. 637 (2009) 73–104. [12] G.F. Dietze, A. Leefken, R. Kneer, Investigation of the backflow phenomenon in falling liquid films, J. Fluid Mech. 595 (2008) 435–459. [13] E.O. Doro, C.K. Aidun, Interfacial waves and the dynamics of backflow in falling liquid films, J. Fluid Mech. 726 (2013) 261–284. [14] M.M. Francois, S.J. Cummins, E.D. Dendy, D.B. Kothe, J.M. Sicilian, M.W. Williams, A balanced-force algorithm for continuous and sharp interfacial surface tension models within a volume tracking framework, J. Comput. Phys. 213 (2006) 141–173. [15] D.P. Frisk, E. Davis, The enhancement of heat transfer by waves in stratified gas–liquid flow, Int. J. Heat Mass Transfer 15 (1972) 1537–1552. [16] G.D. Fulford, The flow of liquids in thin films, Adv. Chem. Eng. 5 (1964) 151– 236. [17] D. Gao, N.B. Morley, V. Dhir, Numerical simulation of wavy falling film flow using VOF method, J. Comput. Phys. 192 (2003) 624–642. [18] C.W. Hirt, B.D. Nichols, Volume of fluid (VOF) method for the dynamics of free boundaries, J. Comput. Phys. 39 (1981) 201–225. [19] M. Ishii, T. Hibiki, Thermo-fluid dynamics of two-phase flow, Springer, 2011.

C. Albert et al. / International Journal of Heat and Mass Transfer 69 (2014) 343–357 [20] M.A. Islam, A. Miyara, T. Setoguchi, Numerical investigation of steam absorption in falling film of LiBr aqueous solution with solitary waves, Int. J. Refrig. 32 (2009) 1597–1603. [21] K. Jöhnisch, M. Baerns, V. Hessel, W. Ehrfeld, V. Haverkamp, H. Löwe, C. Wille, A. Guber, Direct fluorination of toluene using elemental fluorine in gas/liquid microreactors, J. Fluorine Chem. 105 (2000) 117–128. [22] P. Kapitza, S. Kapitza, Wave flow of thin layers of a viscous fluid: III. Experimental study of undulatory flow conditions, Zh. Eksp. Teor. Fiz. 19 (1949). D. ter Haar (Ed.), Collected Papers of P.L. Kapitza, Macmillan, New York, 1964, vol. II, pp. 690–709. [23] S. Kulankara, K.E. Herold, Theory of heat/mass transfer additives in absorption chillers, HVAC& R Res. 6 (2000) 369–380. [24] N.A. Malamataris, V. Balakotaiah, Flow structure underneath the large amplitude waves of a vertically falling film, AlChE J. 54 (2008) 1725–1740. [25] D. Moalem Maron, N. Brauner, G.F. Hewitt, Flow patterns in wavy thin films: numerical simulation, Int. Commun. Heat Mass Transfer 16 (1989) 655– 666. [26] T. Nosoko, P.N. Yoshimura, T. Nagata, K. Oyakawa, Characteristics of twodimensional waves on a falling liquid film, Chem. Eng. Sci. 51 (1996) 725–732. [27] W. Nusselt, Die Oberflächenkondensation des Wasserdampfes, Z. Ver. Dtsch. Ing. 60 (1916) 541–546. [28] C.D. Park, T. Nosoko, Three-dimensional wave dynamics on a falling film and associated mass transfer, AlChE J. 49 (2003) 2715–2727. [29] H. Raach, J. Mitrovic, Simulation of heat and mass transfer in a multi-effect distillation plant for seawater desalination, Desalination 204 (2007) 416–422. [30] A. Rastaturin, E. Demekhin, E. Kalaidin, Optimal regimes of heat-mass transfer in a falling film, J. Non-Equilib. Thermodyn. 31 (2006) 1–10. [31] W.J. Rider, D.B. Kothe, Reconstructing Volume Tracking, J. Comput. Phys. 141 (1998) 112–152. [32] M. Rieber, Numerische Modellierung der Dynamik freier Grenzflächen in Zweiphasenströmungen, Ph.D. thesis, University of Stuttgart, 2004.

357

[33] A. Schagen, M. Modigell, Local film thickness and temperature distribution measurement in wavy liquid films with a laser-induced luminescence technique, Exp. Fluids 43 (2007) 209–221. [34] A. Schagen, M. Modigell, G. Dietze, R. Kneer, Simultaneous measurement of local film thickness and temperature distribution in wavy liquid films using a luminescence technique, Int. J. Heat Mass Transfer 49 (2006) 5049–5061. [35] V. Shkadov, Wave flow regimes of a thin layer of viscous fluid subject to gravity, Fluid Dyn. 2 (1967) 29–34. [36] G.M. Sisoev, O.K. Matar, C.J. Lawrence, Absorption of gas into a wavy falling film, Chem. Eng. Sci. 60 (2005) 827–838. [37] J.C. Slattery, L. Sagis, E.S. Oh, Interfacial transport phenomena, Springer, 2007. [38] H.B. Squire, On the stability for three-dimensional disturbances of viscous fluid flow between parallel walls, Proc. R. Soc. London. Ser. A 142 (1933) 621–628. [39] F.K. Wasden, A.E. Dukler, A numerical study of mass transfer in free falling wavy films, AlChE J. 36 (1990) 1379–1390. [40] Z. Xu, B. Khoo, N. Wijeysundera, Mass transfer across the falling film: simulations and experiments, Chem. Eng. Sci. 63 (2008) 2559–2575. [41] K.K. Yeong, A. Gavriilidis, R. Zapf, V. Hessel, Catalyst preparation and deactivation issues for nitrobenzene hydrogenation in a microstructured falling film reactor, Catal. Today 81 (2003) 641–651. [42] C.S. Yih, Stability of two-dimensional parallel flows for three-dimensional disturbances, Quart. Appl. Math 12 (1955). [43] S.M. Yih, C.C. Kuo, Design and testing of a new type of falling film gas–liquid contacting device, AlChE J. 34 (1988) 499–501. [44] P.N. Yoshimura, T. Nosoko, T. Nagata, Enhancement of mass transfer into a falling laminar liquid film by two-dimensional surface waves – some experimental observations and modeling, Chem. Eng. Sci. 51 (1996) 1231– 1240. [45] L. Zhang, H. Zheng, Y. Wu, Experimental study on a horizontal tube falling film evaporation and closed circulation solar desalination system, Renewable Energy 28 (2003) 1187–1199.