증발액막의 비선형안정성해석: 2차원 및 3차원 해석

이 재 성, 김 효, 김 덕 찬

서울시립대 화학공학과

Nonlinear Stability For The Evaporating Liquid Film:

Two And Three Dimensional Disturbances

Jae-Sung Lee , Hyo Kim  And Dok Chan Kim

Dept. of Chemical Engineering, Seoul City University

INTRODUCTION

Thin evaporating liquid film provides excellent removal of heat on the hot solid surfaces due to its high transfer surface of heat and mass in comparison with the volume of the film. In contrast with this advantage the evaporation destabilizes the film and thus the stability of this liquid film is of practical importance in several applications, for example the safety analysis of reactor and heat exchanger. Here we consider extremely thin films, which are thick enough to be treated by continuum mechanics but thin enough that gravity effects are negligible and long-range molecular forces are significant in determining their dynamics. A liquid layer placed on a horizontal heated plate can experience a variety of hydrodynamic instabilities due to van der Waals attractions and thermocapillarity induced the surface tension gradient at the interface and vapor recoil. Burelbach et al[1] used an one-sided model of evaporation contained evaporative mass loss, vapor recoil, thermocapillarity and van der Waals attractions so that they derived the long-wave evolution equation for the interface shapes that govern the two-demensional nonlinear stability of films subject to the above coupled mechanisms. We'll extend study linear and nonlinear stabilities in three dimensional layers as well as in two dimensional cases and for the nonlinear stability the Fourier-Spectral method has been employed to calculate the film height.

FORMULATION

We consider an incompressible Newtonian fluid of constant density[eqn: rho], viscosity[eqn: mu], thermal conductivity [eqn: k] and thermal diffusivity[eqn: kappa], lying on a uniformly heated flat surface at a constant temperature [eqn: T_H]. The layer is bounded above by its vapor of temperature[eqn: T_S]. Cartesian coordinates are used to describe the system. In two dimensional case the normal coordinate[eqn: z] is taken to be zero at the solid surface and the interface is located at [eqn: z=h(x,t)], where [eqn: h] is a function of the lateral coordinate [eqn: x] and time[eqn: t] (Fig.1). The unit vectors [eqn: n] and[eqn: t] are the outward normal and tangent vectors, respectively. Evaporation is occured only in a direction normal to the interface. If the film is taken to be a nonionic liquid and thin enough that inertial effects are negligible, we can describe the Navier-Stokes equations augmented with an extra body force,[eqn: del phi] due to van der Waals attractions, such as

[eqn: rho`(`vec v_t `+`vec v CDOT `vec TRIANGLED vec v`) =- vec TRIANGLED p `-` vec TRIANGLED phi`+`mu `{vec TRIANGLED} ^2 `vec v],	(1)

where the potential funtion [eqn: phi] is a function of film thickness [eqn: h], [eqn: p] the pressure and [eqn: vec v] the liquid velocity with components[eqn: u] and[eqn: w] in the [eqn: x]-and[eqn: z]-directions, respectively. The continuity and energy equations are

[eqn: vec TRIANGLED CDOT ``vec v``=``0] ,	(2)

[eqn: T_t ``+``vec v CDOT `vec TRIANGLED T` ``=`` kappa ``{vec TRIANGLED} ^2 `T] .                   (3)

The boundary conditions are described as, at[eqn: z=0], [eqn: vec v``=``0``,``T``=``T_H],	(4)

and at [eqn: z=h(x,t)], jump mass and jump energy balances are

[eqn: J=rho(vec v -{vec v}^{(I`)} ) cdot vec n =rho^(V) ({vec v}^{(V`)} -{vec v}^{(I`)} `) cdot vec n],	(5)

[eqn: &J(L+ 1 over 2 [({vec v}^{(V`)} -{vec v}^{(I`)} ) cdot vec n `]^2 - 1 over 2 [(vec v -{vec v}^{(I`)} ) cdot vec n `]^2 ) +k vec del T cdot vec n `- ~~~~~~~~~~~~~(6)## &`k^{(V`)} vec del T^{(V`)} CDOT vec n `+` 2 mu (VEC VEC D cdot vec n ) cdot (vec v -{vec v}^{(I`)} ) `-` 2 mu^{(V`)} ({VEC VEC D} ^{(V`)} cdot vec n ) cdot ({vec v}^{(V`)} -{vec v}^{(I`)} )`=`0.]

In addtion, at [eqn: z=h(x,t)], the following normal stress, shear stress and no slip conditions are obtained

[eqn: J(vec v `-`{vec v}^{(V`)} ) cdot vec n ``-` ({vec vec T} `-`{vec vec T}^{(V`)} `) cdot vec n cdot vec n```=```2H sigma] ,	(7)

[eqn: J(vec v `-`{vec v}^{(V`)} ) cdot vec t ``-` ({vec vec T} `-`{vec vec T}^{(V`)} `) cdot vec n cdot vec t```=```- vec del sigma cdot vec t] ,	(8)

[eqn: (VEC v `-`{vec v}^{(V`)} ) cdot vec t ``=``0],	(9)

and the linearized constitutive equation can be derived from Palmer[2]

[eqn: J ``=``LEFT( {alpha rho ^(V) L } OVER T_S ^3/2 right)left(M_w over { 2pi R_g }right)^1/2 (T^{(I`)} - T_S )],	(10)

where [eqn: J(x,t)] is the mass flux due to evaporation,[eqn: T^(I) (x,t)] the temperature of the interface,[eqn: vec vec T`(=`-p vec vec I ``+ 2 mu vec vec D )] stress tensor in the liquid, [eqn: vec vec I] the identity tensor,[eqn: vec vec D] the deformation tensor,[eqn: L] the latent heat of evaporation, [eqn: M_w] the molecular weight, [eqn: R_g] the universal gas constant, and [eqn: alpha] the accomodation coefficient. The superscripts [eqn: I] and [eqn: V]  denote the interface and vapor phase, respectively.  The surface tension is described by

[eqn: sigma ``=``sigma_0 `-` gamma `` (T^{(I`)} `-`T_S `)] ,	(11)

where [eqn: sigma_o] is the surface tension at the reference temperature[eqn: T_S], taken[eqn: T_S] the saturation temperature at the system pressure and[eqn: gamma(=-d sigma/dT)] has positive value for common liquid.

For nondimensionalization of the governing equations and B.C's we need scale factors, that is, the initial film thickness[eqn: h_o] for length scale, while [eqn: nu/h_o] for a viscous velocity scale,[eqn: h_o ^2 /nu] for a viscous time scale and[eqn: rho nu^2 / h_o ^2] for a pressure, where[eqn: nu(=rho/mu)] is the kinematic viscosity of the liquid. The temperature difference [eqn: T-T_S] is scaled on the[eqn: TRIANGLE T =T_H -T_S]. Evaporating time scale[eqn: t_E] is[eqn: rho h_o ^2 L/k TRIANGLE T]. Mass flux is scaled on its initial value for a linear temperature profile across an initially flat film[eqn: k TRIANGLE T/d_o L]. And we get the potential funciton[eqn: phi] followed by Ruckenstein & Jain[3] as[eqn: phi=A h^-3] which is independent of [eqn: z].

BASIC STATE

The basic state is a flat evaporating interface. For slow evaporation, the evaporation number [eqn: E], the ratio of the viscous to evaporative time scales, is a small parameter. Since the velocity field in the basic state is zero, time rescaled on the evaporative scale to retain the effect of mass loss. The dependent variables are expanded in powers of [eqn: E], and substituted in the temperature equation and boundary conditions. The resulting leading order basic-state solution is

[eqn: bar h```=``-K``+``(``K^2``+``2 `K ``+``1`` -``2`E`t``)`^1/2],	(12)

[eqn: bar T``=``1``-``(``K^2``+``2 `K ``+``1`` -``2`E`t``)`^-1/2`` z] ,	(13)

[eqn: bar J``=``(``K^2``+``2 `K ``+``1`` -``2`E`t``)`^-1/2],	(14)

where [eqn: bar h , bar T ,and { bar J}] are the nondimensional time-dependent film thickness, temperature profile, and surface mass flux, respectively, and

[eqn: K`=` LEFT({k T_S `^{3/2} } OVER {alpha` d_o `rho^(v)`L`^2 }RIGHT)``left({2`pi`R_g }over{M_W } right)^{1/2}].                                          (15)

The parameter [eqn: K] measures the degree of nonequilibrium at the evaporating interface. [eqn: K=0] corresponds to the quasi equilibrium limit where the interfacial temperature is constant and equal to the saturation value,[eqn: T=0.][eqn: K^-1 =0] corresponds to the nonvolatile case in which the evaporative mass flux [eqn: J] is zero.

LONG-WAVE APPROXIMATION

The characteristic length of the interfacial disturbance is much larger than the thickness of the liquid layer. So we consider long-wave disturbance periodic in [eqn: x] and define as a small parameter [eqn: k], wavenumber. We  rescale the governing equations and boundary conditions consistent with lubrication theory by writing [eqn: xi= kx, zeta=z,tau=kt]. The leading order system contains the dimensionless velocity components,[eqn: u_o] and[eqn: w_o], as well as dimensionless temperature and pressure,[eqn: T_o] and[eqn: p_o]. From the kinematic condition, we obtain a nonlinear evolution equation for the film thickness:

[eqn: h_t' + epsilon & `(h+K)^-1 +(h^3 `h_XXX`)_X`+## &([h^-1 + OMEGA (h+K)^-3 h^3 +K PSI (h+K)^-2 h^2 ] h_X`)_X =0 ,]                    (16)

where[eqn: X ``=``(`A`/`S`)^1/2 `x ``,``t'``=``(`A^2 `/`S`)``t],[eqn: epsilon ``=``E`S`/`A`^2 ``,`` OMEGA ``=``E`^2 `D^-1 /A],[eqn: PSI ``=``M P^-1 /A].

It is assumed that [eqn: Ak, Dk^-3 ,Ek^-1 ,Mk] and [eqn: Sk^3] are all [eqn: O(1)]. The parameter [eqn: S``=``{ sigma _o d_o } /{3 rho nu`^2 }],[eqn: M ``=``{gamma `TRIANGLE T `d_o }/{2 `rho `nu `k}] and [eqn: D``=``3 rho^(v) / 2 rho] describe nondimensional surface tension, Marangoni number and the ratio of the vapor to liquid densities in that order.

LINEAR STABILITY

To take into consideration of the stability of the liquid film, first we employ  linear stability theory and perturb the basic state[eqn: bar h] by a small amount[eqn: h'] which is treated by normal modes of the form[eqn: h'(X,t'`)=H(t'`)exp^{i`k`X}]. Then we can get a resulting  ordinary differential equation for the normal mode amplitude [eqn: H] with initial condition [eqn: H=H_o] at [eqn: t'=0]. To define the effective growth rate [eqn: Gr(t')] the exponent in resulted [eqn: h'(x,t')] is divided by time[eqn: t']. We find a range of unstable wavenumbers [eqn: 0<k`^2 < k_c `^2] where cutoff wavenumber [eqn: k_c] is given by [eqn: Gr(t')=0].

NONLINEAR STABILITY

The nonlinear development of the instability can be studied by posing an initial value problem in a spatially-periodic domain an integrating the evolution equation numerically. As an initial condition, we use a  small-amplitude disturbances;

[eqn: h`(x,0)``=``1- delta ` cos(`k_M `x`)],	(17)

where[eqn: delta] is 0.1 and[eqn: k_M `(=1/ SQRT { 2})] is the fastest growing wavenumber for the isothermal linear analysis. In order to calculate[eqn: h], the Hamming modified predictor-corrector method is used, while the spatial derivatives are computed using the Fourier spectral method.

CONCLUSIONS

We studied the evolution equation contained the effects of mass loss, vapor recoil, van der Waals attractions and Marangoni instability to calculate the film instability  in 2- and 3-dimensional cases. The linear stability in case of[eqn: K=0](Fig 2.) indicates that the maximum growth rate increases with time while the range of unstable wavenumbers gets wider. Fig 3. shows the effects of [eqn: K] for the evaporating film when thermocapillary effects are neglected. As[eqn: K] increases, the liquid approaches near the nonvolatile case and the rupture time becomes longer.

REFERENCES

1. Burelbach, J.P., Bankoff, S.G., and Davis, S.H.: J.Fluid Mech., 195, 463(1988)

2. Palmer, H.J.: J.Fluid Mech., 75, 487(1976)

3. Ruckenstein, E. and Jain, R.K.: J. Chem. Soc. Faraday Trans., II 70, 132(1974)