CC BY
24
MSC 97M50, 74J15, 74J30, 74A99
DOI: 10.14529/ mmp150403
MATHEMATICAL MODELLING OF WAVY SURFACE OF LIQUID FILM FALLING DOWN A VERTICAL PLANE AT MODERATE REYNOLDS' NUMBERS
L.A. Prokudina, South Ural State University, Chelyabinsk, Russian Federation, [email protected],
Ye.A. Salamatov, Emerson Process Management, Metran, Chelyabinsk, Russian Federation, [email protected], [email protected]
Development of periodic disturbances on free surface of water film falling down vertical plane for Reynolds' number Re e [5; 10] is investigated. The investigation is implemented in a scope of the nonlinear differential equation for evolution of free surface of falling down liquid film. The equation is solved by a finite differencies method at rectangular uniformly spaced grid. By researching the growth of unit inaccuracy, the conditions on parameters of computation grid for inaccuracies to be not increasing are obtained. As a result, waveforms of water film, time spent to form the regular wave mode and amplitudes of periodic disturbances are calculated. Calculated amplitudes and experimental ones are compared.
Keywords: liquid film; amplitude; waveform; nonlinear evolution of disturbances.
Introduction
Liquid film flows (thin layers of viscous fluids) are implemented in many industries such as chemical, oil, energy and others. Advantages of liquid film flows are energy efficiency, greater product purity etc. These advantages are determined by highly extended contact line between streams moving and interacting each other. The contact line is formed by waves on free surface of liquid film.
Wavy flows of liquid films have been investigated since the middle of the 20th century [1, 2]. Nevertheless they are still interesting [3-10]. Their study is requiered by both practical (due to variety of environments where liquid film flows are realised) and theoretical reasons (since film flows are described by nonlinear partial differential equations).
1. The Model
Consider a thin layer of wavy liquid film falling down a vertical smooth plane under the action of gravity. The liquid is supposed to be viscous, uniform and incompressible. It is interacting with a gas streaming in parallel to the plane. The gas creates constant tangential stress f on the free surface of the liquid film (Fig. 1). If the liquid and the gas are moving in one direction, cocurrent mode is implemented, otherwise countercurrent mode is realised. The equation
dri , d4ri 1 d2tb 1 dri , d2tb 1 , dri , , d2ri , , d4ri , d2ri
= bi^T + + h*JT + + + bzriw^ + + bsri^r+
dt dx4 dx2 dx dxdt dx dx2 dx4 dxdt , ,
dri dri ! f dri\2 , dri d3ri , , 2dri , ,2d2ri , , dri d3ri
+ b9—--— + bio -FT + bli — + bi2ri -FT + bi3ri "rj o + bi4ri^3
dx dt \ dx J dx dx3 dx dx2 dx dx3
(Г ReG 3
describes evolution of liquid film free surface [10-12]. Here b1 =--—, b2 = 40Re3F(F —
5 3 9
т), Ьз = —Re(F — r), b4 = — Re2F, b5 = —2ReF + Rer, b6 = — -Re3Fr + — Re3F2,
24 8 20
3
bj = 3bb bg = 4b^ bg = be, bio = be, bn = 3bb bi2 = —ReF, bi3 = — 4Ree3Fr +
9
-Re3F2, 614 = 6615 Re is the Reynolds number, F is the Froude number, a is the surface 8
tension (unitless), t is the shear stress (unitless). The required function ^(x,t) denotes the displacement of free surface (Fig. 1) from its undisturbed position ¿0 ^ ¿0). We assume t > 0 for countercurrent mode and t < 0 for cocurrent 0ne. If t = 0, free mode is implemented.
2. Searching for Periodical Solution to the Governing Equation
A periodical solution to (1) is searched numerically by finite defferencies method. We consider a rectangular domain Q = {(x,t) : xbegin ^ x ^ xend, 0 ^ t < tend}. Let us cover
Q by a uniform grid (xi,tj). The grid has a step Ax = end — begin along the spacial
Nx
variable x and has a step At along the time variable t. In this case xi = xbegin + iAx, tj = jAt (i = 0,1,... Nx; j = 0,1, 2,...). By ipj denote ^(xi,tj). The derivatives in (1) are approximated as follows
д 4ф
дх4 д 3ф
дх3 д 2ф
дх2
дф дх
Ф1+2 — 41+1 + Щ — 4ф-1 + ф—
Ах4
d4,i
Ах■
Ф1+2 — 3ф|+1 + 3ф-1 — Фi—2 =
2Ах3 Ах3;
ф+1 — 2ф1 + ф— 1
Ах2 ф+ — ф—
2Ах
d2,i Ах ■
d1,i Ах
(2)
(3)
(4)
(5)
J
J
J
J
dtf j dt
d2tf
dxdt
j - tf At '
1 /tf+i - tf+1
At
Ax
j+i- A
Ax
)
(6) (7)
Substituting (2) - (7) in to (1) and modifying it, we obtain
tfj+1 ( XI +
(
1
bd
+
b8tfj
MlA
At AtAx AtAx AtAx
+ tf
j+i i+1
AtAx
b8tf
AtAx AtAx
)
d4,i
dj
Ax2
i Ax
+ bll ^ + bl2
Ax
j + b
d2i + bjtfj j + bio
Ax2
2
2 j
Ax Ax3 ' vri> Ax + b13yriJ Ax2
jA 2 d2 ,i
+ bUtf
Ax2
__+ j
Ax Ax3 At
dj dj d1 i d3
b4 tf+1 - tfj b8tfj tfj+i - tfj b9di ,tfi
At Ax
At Ax
AtAx
Assume
At = nAx3' n> 0'
0 ^ i ^ Nx. (8) (9)
and, multiplying both sides of (8) by ^Ax4, we get
(«j - j j + Pj j = ml + («j - j tf + Pj tfj+1' 0 ^ i ^ N
(10)
where aj = Ax - bgd^, Pj = -b4 - b8tfj, = b^4i + b2d32iAx2 + b3d31iAx3 + b5tf dj1iAx3 +
betf j d2 , iAx2 + b7tf j dj4 , i + b1o (dj , J Ax2 + budj , id\ , i + b12 (j dj , iAx3 + b13 {tfj) d32 , iAx2 + b14tfj d1 , id3 , i-
Thus, in order to search a solution to (1) it is necessary to solve the system of equations (10).
Further, consider how computational inaccuracies behave while transferring from the jth 'time layer' to the (j + 1)th one. We will do that using a method of a unit inaccuracy growth as described in [13].
So, we assume that in any node (i,j) the function tf contains inaccuracy £j (we
suppose ej ^ tfj). But for the other nodes on the jth 'time layer', we assume that the tfij
(i'j), it can be written that tf = tf + a], here tf marks exactly calculated part of the function tfj (further, we denote by 'bar' any variable calculated exactly).
The inaccuracy ej leads to appearing of inaccuracies in the finite differencies (2) - (7) for the nodes (i - 2j), (i - 1j), (i'j), (i + 1'j), (i + 2'j) as follows
(i - 2'j) (i - 1'j) (i'j ) (i + 1'j ) (i + 2'j)
d1 dj d1, i-2 j ej dj + — d1 , i-1+ 2 4 i j ej dj - — d1 ,i+1 2 dj d1,i+2
d2 dj d2, i-2 d2 , i-1 + ei d2,i - 2ej d2, i+1 + £i dj d2,i+2
d3 j ej dj + — d3 , i-2 + 2 4 i 1 -— 3'i-1 2 j i d 3 , i+1 + 2 j ej dj - — d3,i+2 2
d4 d4 , i-2 + ej d4,i-1 - 4ei d34, i + 6£j d4 , i+1 - 4ej d44 ,i+2 + ej
(11)
b
4+1 43 «й ^ №-U + i)№J + i)
XX XXXXX XX
(Oj + 1) (lj + 1) (i-2j + l) (i-lj + 1) (ij + 1) (i+lj + l) (. + 2,3 + 1)
(OJ) (lj) ej (N, - l,j) (NXJ)
Д---Х X X X X X X X
(<-2j) (i-lj) (¿J) (i + lj) (i + 2J)
4 (-X-X-X-X-X-) <-
%begin %end X
Fig. 2. The sketch of Q j
(j + 1)th 'time layer' as depicted on Fig. 2. For example, Equation (10) written for the node (i + 2, j + 1) has the form
j - fe) j + ej+2) + =
- - £J - ■ -■ £J \ П I Yi+2 + b1el + b^l+2el — b11d1,i+2 2 + b14^l+2d\,i+2 2 +
aj -2 " tJi+2) yn+2 T ci+2j i+3 -
+ , i+2 e^
+ j - M+2) j + M+2^+3- (12)
By ej+}, we denote the inaccuracy in the node (i + 2, j + 1).
Taking into account ipj ^ 50 and ej ^ ipj, modify Equation (12) as follows
(Ax + 64) ej+ ~ Vbiej,
j+1 j ei+2 ~ Ki+2ei ,
and
Ki+2 = n ATX • (13)
Similarly, missing manipulations due to their tedious, we find constants of proportionality for the nodes (i + 1, j + 1) (i,j + 1) (i — 1, j + 1) (i — 2, j + 1):
( 4 6 1 6 A x2 63 A x3 _ _ \
Ki+i = n
—4b1 + Ь2Ах2 — -j Ах3 b1b4
Ах + b4 + (Ах + b4)2
V )
(14)
1+ n (6b1 — 2Ь2Ах2) + b4
K =X+ Ах + b, + Ах+Ь~4 Ki+1' (15)
ц(—4b1 + Ь2Ах2 + b3Ах3 J — b4 + b4Ki
2
Ki—1 = —-Ах + * '-' (16)
n (6b1 — 2b2Ах2) + b4Ki—1
Ki—2 =-Ах+b-• (17)
Equations in (10) for the nodes (l,j + 1), 0 ^ l ^ i — 3, have the form (aj — j ■ + ej+1) + ej (ф+ + e^) = nij + (aj — j + ijj, 0 ^ l ^ i — 3^
Вестник ЮУрГУ. Серия «Математическое моделирование 22 и программирование» (Вестник ЮУрГУ ММП). 2015. Т. 8, № 4. С. 30-39
Thus,
b
4 ej+1
e
ej+1 ^_
e ~ Ax + bSl+1
b4
0 < l < i- 3'
j+1
Kl =
Ax + b4
b4
Ax + b4
Kl+1ej' 0 ^ l ^ i - 3' Kl+1' 0 ^ l ^ i - 3.
(18)
eij
(j + 1)th "time layer". Obviously, the inaccuracies on the (j + 1)th 'time layer' will not increase if parameters n and Ax are chosen so that \Kq\ < 10 ^ q ^ i + 2. According to the formula (18), we get \k0\ < • • • < \k1\ < • • • < \Ki-3\ < \n,i-2\. Therefore we may require only
\Kq\ < 1' i - 2 ^ q ^ i + 2. Unequalities (19) can be written with respect to n as follows
Ax + b4
(19)
n<
n<
\b1\ (Ax + b4)2
n1'
b b b 3b1 b4 + 4b1 Ax - b2b4Ax2 + i 4b1 + ^ I Ax3 + -2 Ax4
n2
n<
n<
2(Ax + b4 )3
3b1b4 + 8b1b4Ax + (6b1 - b2b2)Ax2 - ^b4 + ^^ Ax3 - ^2b2 + ^^
(Ax + b4 )4
- (2b2 + Ax4
(4b1 + ¥)
b1b4 + 4b1b2Ax + 6b1b4Ax2 + 4b1 + Ax3 - ^4 + b:ib^i)Ax4-
- + 3b3b^ Ax5 - b3 Ax6
n3'
n4'
n<
(Ax + b4)5
Ax3+
b b3
5b1b4 + 20b1b4Ax + (30b1b4 - 2b2b4)Ax2 - ( 20b1b4 - 8b2b3 - -y4 + (6b1 - 11b2b4 + b3b4)Ax4 + 0b3b2 - 7b2b^ Ax5 + ^^ - 2b2 ) Axe
n5-
Finally, we obtain the next condition
n < min(n1' n2' n3' n4' n5)-
3. Results
Considere an initial periodical disturbance
tf(x' 0) = a0 • cos(kx)
(20)
(21)
Table
Grid parameters
Re k Дж П x 10"3
5 0,07246 0,270 1,11
6 0,08311 0,273 0,94
7 0,09332 0,276 0,78
8 0,10211 0,273 0,62
9 0,10976 0,282 0.49
10 0,11629 0,274 0,37
that is developed 011 free surface of water film falling down a vertical plane. Disturbance (21) has the amplitude a0 = 10"4 ^ 1 and the wave number k (see Table). We consider Re e [5; 10] Mid t = 0.
Put n = -min(ni,n2,V3,V4,V5)- The grid parameters are placed in Table. The length 3
of [xmin; xmax] is multiple of 30 wave lengths of (21). Calculations were carried out untill tend = 1000.
Our calculations demonstrate, that the regular wavy flows are formed 011 free surface of
vertical water film. The waveforms are shown 011 Fig. 3 (we illustrated only one period and
2n
the x axes was scaled as -j^)' Similar waveforms are observed many times in experiments, the reader may see them in [2, 7, 14, 15], for example.
0.13
-0.13
-0.26
Fig. 3. Waveforms on free surface of water film: - Re = 5; 2 - Re = 7; 3 - Re =10
Re
amplitudes depicted 011 Fig. 4 were calculated by formula [1, 2]
^max ^mir ^max + ^mir
(22)
We used the formula [14. 16]
A
¿max - jS)
jS)
(23)
a
to calculate wave amplitudes depicted 011 thickness is denoted.
Fig. 4. Amplitudes (water) calculated by (22): • - [2]; x - our results
ig. 5. By (6) in (23) the mean liquid film
Fig. 5. Amplitudes (water) calculated by (23): ■ - [14]; □ - [16]; x - our results
Fig. 6. Amplitude time dependence for water film: 1 - Re = 5; 2 - Re = 7; 3 - Re = 10
for water film
Fig. 6 shows how amplitudes change while regular wavy flows are forming. On Fig. 6 we marked the times tpi = 823' 2, tp2 = 364' 3, tp3 = 185' 7 when the regular wave flows
Re
tp is required. References
1. Kapit.za P.L. [Wave Flow of Thin Layer of Viscous Fluid]. Zhurnal EksperirnentaVnoy i Teoreticheskoy Fiziki, 1948, vol. 18, no. 1, pp. 3-18. (in Russian)
2. Kapitza P.L., Kapitza S.P. [Wave Flow of Thin Layer of Viscous Fluid]. Zhurnal Eksperimental'noy i Teoreticheskoy Fiziki, 1949, vol. 19, no. 2, pp. 105-120. (in Russian)
3. Vellingiri R. Dynamics of a Liquid Film Sheared by Co-Flowing Turbulent Gas. Int. J. Multiphase Flow, 2013, vol. 56, pp. 93-104. DOI: 10.1016/j.ijmultiphaseflow.2013.05.011
4. Trifonof Y.Y. Stability of the Wavy Film Falling Down a Vertical Plate: The DNS Computations and Floquet Theory. Int. J. Multiphase Flow, 2014, vol. 61, pp. 73-82. DOI: 10.1016/j.ijmultiphaseflow.2014.01.006
5. Dietse G.F., Al-Sibai F., Kneer R. Experimental Study of Flow Separation in Laminar Falling Liquid Films. J. Fluid Mech., 2009, vol. 637, pp. 73-104. DOI: 10.1017/S0022112009008155
6. Dietse G.F., Ruyer-Quil Ch. Wavy Liquid Films in Interaction with a Confined Laminar Gas Flow. J. Fluid Mech., 2013, vol. 722, pp. 348-393. DOI: 10.1017/jfm.2013.98
7. Meza C.E., Balakotaiah V. Modeling and Experimental Studies of Large Amplitude Waves on Vertically Falling Films. Chern. Eng. Sci.„ 2008, vol. 68, no. 19. pp. 4707-4734. DOI: 10.1016/j.ces.2007.12.030
8. Ruyer-Quil Ch., Mannaville P. Improved Modeling of Flows Down Inclined Planes. The European Physical Journal B, 2000, vol. 15, pp. 357-369. DOI: 10.1007/sl00510051137
9. Prokudina L.A., Vyatkin G.P. Self-organization of Perturbations in Fluid Films. Doklady Physics, 2011, vol. 56, no. 8. pp. 444-447. DOI: 10.1134/S1028335811080039
10. Prokudina L.A. Influence of Surface Tension Inhomogeneity on the Wave Flow of a Liquid Film. J. Eng. Phys. Thermophys., 2014, vol. 87, no. 1, pp. 165-173. DOI: 10.1007/sl0891-014-0996-2
11. Prokudina L.A., Vyatkin G.P. Instability of a Nonisothermal Liquid Film. Doklady Physics, 1998, vol. 43, no. 10, pp. 652-654. (in Russian)
12. Prokudina L.A., Salamatov E.A. Shear Stress Influence on Wave Characteristics of Liquid Film. Bulletin of the South Ural State University. Series: Mathematics. Mechanics. Physics, 2012, vol. 7, no. 34, pp. 173-176. (in Russian)
13. Berezin I.S., Zhidkov N.P. Metody vychisleniy [Computing Methods]. Moscow, FIZMATLIT, 1959. 620 p. (in Russian)
14. Rogovaya I.A., Olevskiy V.M., Runova N.G. [Measurement of Parameters of Wavy Liquid Film Flow on Vertical Plate]. Teoreticheskie osnovy khimicheskoi tekhnologii [Theoretical Foundations of Chemical Engineering], 1969, vol. 3, no. 2, pp. 200-208. (in Russian)
15. Lei V.V., Al-Sibai F., Leefken A., Renz U. Local Thickness and Wave Velocity Measurement of Wavy Films with a Chromatic Confocal Imaging Method and a Fluorescence Intensity Technique. Exprements in Fluids, 2005, vol. 39, no. 5, pp. 856-864. DOI: 10.1007/s00348-005-0020-x
16. Penev V., Krylov V.S., Boyadjiev Ch., Vorotilin V.P. Wavy Flow of Thin Liquid Films. Int. J. Heat Mass Transfer, 1972, vol. 15, no. 7, pp. 1395-1406. DOI: 10.1016/0017-9310(72)90019-1
Received January 19, 2015
УДК 532.5 DOI: 10.14529/mmpl50403
МОДЕЛИРОВАНИЕ ВОЛНОВОЙ ПОВЕРХНОСТИ ВЕРТИКАЛЬНОЙ ЖИДКОЙ ПЛЕНКИ ПРИ УМЕРЕННЫХ ЧИСЛАХ РЕЙНОЛЬДСА
Л.А. Прокудина, Е.А. Саламатов
В рамках модельного нелинейного дифференциального уравнения, которое описывает эволюцию свободной поверхности вертикальной жидкой пленки, исследовано развитие периодических возмущений на вертикальной пленке воды в диапазоне чисел Рейнольдса Re £ [5; 10]. Модельное уравнение решалось методом конечных разностей на прямоугольной равномерной сетке. Исследованием роста единичной погрешности получены условия, которым должны удовлетворять параметры расчетной сетки, чтобы погрешности не увеличивались в ходе вычислений. В результате рассчитаны профили волн на свободной поверхности вертикальной пленки воды, время формирования режима стекания с регулярным профилем волны, характер изменения амплитуды периодического возмущения при формировании режима с регулярным профилем волны. Осуществлено сравнение рассчитанных амплитуд с результатами экспериментов.
Ключевые слова: жидкая пленка; амплитуда; профиль волны; нелинейное развитие возмущений.
Литература
1. Капица, П.Л. Волновое течение тонких слоев вязкой жидкости / П.Л. Капица // Журнал экспериментальной и теоретической физики. - 1948. - Т. 18, № 1. - С. 3-18.
2. Капица, П.Л. Волновое течение тонких слоев вязкой жидкости / П.Л. Капица, С.П. Капица // Журнал экспериментальной и теоретической физики. - 1949. - Т. 19, № 2. -С. 105-120.
3. Vellingiri, R. Dynamics of a Liquid Film Sheared by Co-Flowing Turbulent Gas / R. Vellingiri // Int. J. Multiphase Flow. - 2013. - V. 56. - P. 93-104.
4. Trifonof, Y.Y. Stability of the Wavy Film Falling Down a Vertical Plate: The DNS Computations and Floquet Theory / Y.Y. Trifonof // Int. J. Multiphase Flow. - 2014. -V. 61. - P. 73-82.
5. Dietse, G.F. Experimental Study of Flow Separation in Laminar Falling Liquid Films / G.F. Dietse, F. Al-Sibai, R. Kneer //J. Fluid Mech. - 2009. - V. 637. - P. 73-104.
6. Dietse, G.F. Wavy Liquid Films in Interaction with a Confined Laminar Gas Flow / G.F. Dietse, Ch. Ruyer-Quil // J. Fluid Mech. - 2013. - V. 722. - P. 348-393.
7. Meza, C.E. Modeling and Experimental Studies of Large Amplitude Waves on Vertically Falling Films / C.E. Meza, V. Balakotaiah // Chem. Eng. Sci. - 2008. - V. 68, № 19. -P. 4707-4734.
8. Ruyer-Quil, Ch. Improved Modeling of Flows Down Inclined Planes / Ch. Ruyer-Quil, P. Manneville // The European Physical Journal B. - 2000. - Vol. 15. - P. 357-369.
9. Прокудина, Л.А. Самоорганизация возмущений в жидких пленках / Л.А. Прокудина, Г.П. Вяткин // Доклады АН. - 2011. - Т. 439, № 4. - С. 481-484.
10. Прокудина, Л.А. Влияние неоднородности поверхностного натяжения на волновое течение жидкой пленки / Л.А. Прокудина // Инженерно-физический журнал. - 2014. -Т. 87, № 1. - С. 158-166.
11. Прокудина, Л.А. Неустойчивость неизотермической жидкой пленки / Л.А. Прокудина, Г.П. Вяткин // Доклады АН. - 1998. - Т. 362, № 6. - С. 770-772.
12. Прокудина, Л.А. Влияние касательного напряжения на волновые характеристики жидкой пленки / Л.А. Прокудина, Е.А. Саламатов // Вестник ЮУрГУ. Серия: Математика. Механика. Физика. - 2012. - Т. 7, № 34. - С. 173-176.
13. Березин, И.С. Методы вычислений / И.С. Березин, Н.П. Жидков. - М.: ФИЗМАТЛИТ, 1959. - 620 с.
14. Роговая, И.А. Измерение параметров пленочного волнового течения на вертикальной пластине / И.А. Роговая, В.М. Олевский, Н.Г. Рунова // Теоретические основы химической технологии. - 1969. - Т. 3, № 2. - С. 200-208.
15. Local Thickness and Wave Velocity Measurement of Wavy Films with a Chromatic Confocal Imaging Method and a Fluorescence Intensity Technique / V.V. Lei, F. Al-Sibai, A. Leefken, U. Renz // Exprements in Fluids. - 2005. - V. 39, № 5. - P. 856-864.
16. Wavy Flow of Thin Liquid Films / V. Penev, V.S. Krylov, Ch. Boyadjiev, V.P. Vorotilin // Int. J. Heat Mass Transfer. - 1972. - V. 15, № 7. - P. 1395-1406.
Людмила Александровна Прокудина, доктор физико-математических наук, доцент, кафедра «Прикладная математика:», Южно-Уральский государственный уни-вврситбт (г. Челябинск, Российская Федерация), [email protected].
Евгении Александрович Саламатов, К&НДИДсХТ физико-математических наук, Emerson Process Management, ЗАО «ПГ«Метран» (г. Челябинск, Российская Федерация), [email protected], [email protected].
Поступила в редакцию 19 января 2015 г.
