ВЕСТНИК ПЕРМСКОГО УНИВЕРСИТЕТА
2011
УДК 536.25
Серия: Физика
Вып. 2 (17)
Исследование устойчивости потока жидкости над насыщенной пористой средой
Д. Т. Байдина, Т. П. Любимова
Пермский государственный университет, 614990, Пермь, ул. Букирева, 15
Изучается устойчивость стационарного потока жидкости над насыщенной пористой средой, применительно к задаче устойчивости воды над дном, покрытым растительностью. Рассматриваемая двухслойная система слой чистой жидкости - пористый слой, насыщенный той же жидкостью - описывается на основе единого подхода (с помощью системы уравнений движения в насыщенной пористой среде), а наличие двух разных слоев моделируется введением переменных (зависящих от вертикальной координаты) проницаемости и пористости. Поверхность, ограничивающая систему снизу, предполагается твердой, верхняя граница - свободной и недеформируемой. Получены нейтральные кривые устойчивости и зависимости минимального числа Рейнольдса Re(k) и волнового числа наиболее опасных возмущений от отношения толщины пористого слоя к полной толщине. Показано, что зависимость Rem (d) немонотонна.
Ключевые слова: устойчивость, пористая среда.
1. Введение
При движении потока жидкости над пористой средой с достаточно высокой проницаемостью и пористостью поток частично увлекает жидкость в пористой среде, однако, из-за больших сил сопротивления возникает резкий перепад касательной компоненты скорости, что способствует развитию неустойчивости, аналогичной неустойчивости Кельвина-Г ельмгольца.
В работе исследуется устойчивость потока жидкости над насыщенной пористой средой применительно к задаче устойчивости потока воды над дном, покрытым растительностью. Данная задача актуальна в связи с приложениями к проблеме вымывания загрязняющих веществ из слоя придонной растительности. Оценки дают для растений, расположенных на расстоянии 10 см, значение Дарси, используемое в задаче, порядка 10~4.
2. Постановка задачи. Определяющие уравнения
Рассмотрим двухслойную систему, состоящую из слоя вязкой несжимаемой жидкости и расположенного под ним пористого слоя, насыщенного той же жидкостью (рис.1). Поверхность, ограничивающая систему снизу, предполагается твердой, верхняя граница - свободной и недеформируемой.
Рис. 1. Геометрия задачи
Движение в вязкой несжимаемой жидкости описывается уравнениями
дУ Vp
+ У -VV =-----^+vAV + g,
д р
divV = 0, (1)
где V - скорость жидкости, р - плотность, р -давление, g - ускорение силы тяжести.
Для описания движения в пористой среде наиболее часто применяются модели Дарси и Бринкмана [1]. При использовании модели Дарси предполагается, что скорость фильтрации пропорциональна градиенту давления:
К^7
v=------Vp,
V
(2)
где V - скорость фильтрации, К - проницаемость, х - динамическая вязкость.
Скорость фильтрации пропорциональна средней скорости жидкости в порах: V = фV , где ф -
пористость.
© Байдина Д.Т, Любимова Т. П., 2011
Уравнения, описывающие движение жидкости в пористой среде в рамках модели Дарси, в которых использована средняя скорость движения жидкости в порах V, имеют вид
p
ґ dV Л
---+V■VV
V d У
И
= - Vp--------—V + pg,
K (3)
div(^V) = 0;
если использована скорость фильтрации v , то следующий:
P<P~l (f + VM^lv)) = -Vp-tv + pg, (4)
div v = 0.
В модели Бринкмана учитывается вязкая диффузия импульса. В рамках этой модели уравнения, описывающие движение жидкости в пористой среде, применительно к средней скорости движения жидкости в порах записываются в виде
p{jV + У ■Vy j = -VP +VAV + Pg, (5)
div(<pV) = 0.
Здесь V - эффективная вязкость.
В настоящей работе движение жидкости в двухслойной системе слой чистой жидкости - пористый слой, насыщенный той же жидкостью, описывается на основе единого подхода (с помощью системы уравнений движения в насыщенной пористой среде), а наличие двух разных слоев моделируется введением переменных (зависящих от вертикальной координаты) проницаемости и пористости. Движение жидкости в насыщенной пористой среде описывается в рамках модели Бринкмана, с учетом средней скорости движения жидкости в порах. Эффективная вязкость V полагается равной v . Предположение о том, что V = v можно применять в случае высокой пористости среды [8], что справедливо в рассматриваемом в статье случае.
Выберем следующие единицы измерения: для длины - полную толщину h, для скорости - характерную скорость U = (gh2 sina)/v - удвоенную максимальную скорость стекания слоя однородной жидкости по наклонной плоскости, для времени - h2 / v , для давления - pgh sin а , где а - угол наклона слоя к горизонтали. Уравнения в безразмерной форме имеют вид
ЯI/
-----+ Re V ■VV =-Vp - Da-1- V + ДV - (sina)-1 у,
дt к
div(—V) = 0.
(6)
ординаты, К0 - характерное значение проницаемости.
На твердой нижней границе учитывается условие прилипания:
V = 0, (7)
на верхней, свободной границе - условия отсутствия нормальной компоненты скорости и касательных напряжений:
'дК„ дК ,
„ , , „ (8)
=
+ —- | = 0 . дz дx
3. Стационарное течение
Ограничимся рассмотрением плоской задачи д / ду = 0. Будем искать стационарное решение задачи (6)-(8), соответствующее плоскопараллельному течению:
V = К( ^), 0,0). (9)
Уравнения, описывающие стационарное решение (9), в проекциях на направления х (вдоль слоя) и г (поперек слоя) имеют вид
дР0
dx
= -Da-1 —У +
d% dz 2
- +1,
др0
= -ctga .
Dz
Граничные условия:
z = 0 : V0 = 0,
z = l: DVo = 0.
dz
(10)
(11)
(12)
(13)
В задаче (10)-(13) используются три параметра: число Дарси Ба, отношение пористости к проницаемости ф/к и угол наклона слоя к горизонтали а .
Интегрируя уравнение (11) по г , получаем
р0 =-х ctgа + С(х). (14)
Из граничного условия для давления на свободной поверхности
Ро(1) = 0 (15)
(давление на свободной поверхности принято за начало отсчета) следует
С (х) = ^а. (16)
Окончательно получаем
Ро = (1 -z)ctga .
(17)
Из выражения (17) следует, что dp0 / dx = 0, для V получается уравнение
д2У
-Da1 — VO + 2
к dz
с граничными условиями
+1 = 0
Здесь Ба = К0 /к - число Дарси, Яе = ПН /у -число Рейнольдса, к = К / К0 - безразмерная проницаемость среды, зависящая от вертикальной ко-
z = O: V = O
z = 1:
У
dz
-= 0.
(1В)
(19)
(20)
к
В качестве функции ф/к использовалось выражение
ф_ (Л - г) (Л-1)
— = ,--- - , . (21)
к ^1 + а-(Л - г)2 -^1 + а-(Л-1)2
Выражение (21) содержит 2 параметра: параметр а , определяющий ширину переходного слоя между слоями чистой жидкости и пористой среды, насыщенной жидкостью, и параметр Л - отношение толщины пористого слоя к полной толщине. На рис. 2, а-б функция ф / к изображена для нескольких значений параметров а и с!.
Рис. 2а. Зависимость функции ф/ к от вертикальной координаты для й = 0.7 : 1) а = 50,2) а = 100, 3) а = 150
Уравнение (22) решалось методом скалярной прогонки. При вычислении скорости на верхней границе из граничного условия (20) производная по г аппроксимировалась односторонней разностью. Основные расчеты проводились для фиксированных значений параметров Ба и а: Ба = 10-4, а = 100. Параметр Л варьировался.
На рис. 3, а-б изображены полученные численно профили скорости стационарного течения при различных значениях Л .
Рис. 2б. Зависимость функции ф/ к от вертикальной координаты для а = 100:
1) й = 0.5 , 2) й = 0.7, 3) й = 0.9
Задача (18)-(20) решалась численно, методом конечных разностей. Конечно-разностный аналог уравнения (18) записывался в виде
Ккм Ж°г + ^ 1 - Ба-1 \ф\ У0І = -1. (22)
б
Рис. 3. Профиль скорости стационарного течения для Л = 0.35 (а) и 0.7 (б).
4. Устойчивость стационарного течения
Рассмотрим устойчивость полученного стационарного течения. Запишем скорость и давление в виде сумм стационарных полей и возмущений
V = V0+V', р = р0 + р . Считая возмущения малыми, запишем линеаризованные уравнения малых возмущений: дV'
— + Яе (у-УУ + V-Уу ) =
-1 фу к
Аі\(фУ ') = 0.
= -Ур' - Ба- V +АУ’, к
(23)
(24)
а
к
Граничные условия:
2 = 0: V=0,
7 = 1: К„'= 0, = 0.
В дальнейшем штрихи у возмущений опускаются.
Ограничимся рассмотрением плоских возмущений:
V = (и,0^ ). (25)
Применяя к уравнению (23) операцию го^, получаем
д ( ди дw , дt I дг дх )
+ Re
^ д f du dw і d2V0 —' dV0 ^
1 "7 -------------w
— dz
Vo — I--------------------------| + w -
dx V dz dx ) dz
= -Da-l — f* I-Da-1 -if— I u +
к Vdz dx ) dz V к ,
,i du dw
+ДІ----------------
dz dx
Из уравнения (24) имеем
du dw —
— +----+—w = 0.
dx dz —
(26)
(27)
Вводя завихренность а = rotyV , получаем из уравнения (26)
да „ f да і
— + Re I Vo — + w Vo------------------------------Vo w I =
dt V dx —
= -Da-1 — а - Da-1 f — і u + Да .
— + Re I ikV0Q + WV--—- V0" W і =
dt V —
- -Da-1 — Q- Da-l | u + (q""- k2q) .
W " - k 2W + | — W і =-ik Q —
или, вводя новую переменную Q = -ikQ :
(29)
(3O)
— + ik Re I V0Q - V0" W + — VO W і =
dt V —
-Da-1 —Q-Da-1 — I W" + — W і +
—
+ (Q"- k 2q),
W"-k2W + f — W і =Q .
—
(31)
(32)
(33)
Граничные условия:
2 = 0: W = W" = 0, 2 = 1: W = W " = 0.
В качестве функции ф(г), определяющей зависимость пористости от вертикальной координаты, использовалось выражение
ф(г) = Ф0 +(1 - Ф0 ) Л(Ьг) / Ш(Ь). (34)
В этом выражении ф0 - значение пористости при г = 0, параметр Ь определяет толщину переходного слоя между слоями жидкости и насыщенной пористой среды. Основные расчеты проведены при ф0 = 0.1, Ь = 5. Вид функции ф(г) для
ф0 = 0.1, Ь = 5 показан на рис. 4.
(2В)
Будем рассматривать периодические по x возмущения:
w = W(z, t) exp (ikx), а = Q(z, t) exp (ikx) .
Подставляя возмущения в таком виде в (27)-
(28), получим уравнение для амплитуд W и Q в виде
Рис. 4. Вид функции ф(г)
Задача (31)-(33), определяющая эволюцию малых возмущений основного состояния во времени, решалась численно методом конечных разностей. Уравнение (32) решалось на каждом шаге по времени методом скалярной прогонки. Значения завихренности скорости на твердой границе получались по формуле Тома
Ц, = -^1. (35)
Н
Основные вычисления проводились с использованием равномерной сетки с числом узлов, равным 100. Расчеты для каждого набора параметров
к
к
к
к
к
проводились до тех пор, пока временная эволюция возмущений не начинала определяться только наиболее опасной модой, т.е. установления постоянного значения инкремента возмущений, который вычислялся на каждом шаге по времени по формуле
Я = -
1
*н+\ *п-\
ІП
/п+1
/п-1
(36)
Г раница линейной устойчивости стационарного течения относительно данной моды определяется условием X = 0.
5. Результаты расчетов
На рис. 5 приведены нейтральные кривые устойчивости стационарного течения для различных значений отношения толщины пористого слоя к полной толщине Л . Как видно, при Л > 0.5 с увеличением толщины пористого слоя Л, т.е. с уменьшением толщины жидкого слоя, порог неустойчивости повышается. Это связано с тем, что неустойчивость определяется прежде всего возмущениями, развивающимися в жидком слое. При Л < 0.5 уменьшение толщины пористого слоя приводит к повышению порога устойчивости.
Рис. 6. Зависимость минимального числа Рейнольдса Яет1П от отношения толщины пористого слоя к полной толщине й
Рис. 7. Зависимость критического волнового числа ктіп от отношения толщины пористого слоя к полной толщине й
На рис. 8, а-д представлены изолинии действительной части амплитуды возмущений W для различных значений й .
Рис. 5. Нейтральные кривые
Зависимость минимального критического числа Рейнольдса Яетт от Л приведена на рис. 6.
Значение волнового числа наиболее опасных возмущений &тт монотонно растет с увеличением толщины пористого слоя (см. рис. 7).
а
б
01 23456789 10
01 23456789 10
д
Рис. 8. Изолинии действительной части Ж для Л = 0.2 (а), 0.35 (б), 0.5 (в), 0.7 (г),
0.9 (д)
6. Заключение
Решена задача о стационарном течении жидкости над насыщенной пористой средой. Показано, что профиль течения имеет две точки перегиба. Наличие точек перегиба приводит к неустойчивости стационарного течения. Получены нейтральные кривые устойчивости и зависимости минимального числа Рейнольдса Яет1П и волнового числа наиболее опасных возмущений от отношения Л толщины пористого слоя к полной толщине. Показано, что зависимость Кетт (Л) немонотонна: при Л > 0.5 увеличение Л (уменьшение толщины
жидкого слоя) приводит к повышению порога неустойчивости, что связано с тем, что при таких d возмущения локализованы в жидком слое; при d < 0.5 увеличение d приводит к понижению порога неустойчивости.
Список литературы
1. Nield D. A., Bejan А. Convection in Porous Media. 2006. P. 654.
2. Chang M. H., Chen F., Straughan B. Instability of Poiseuille flow in a fluid overlying a porous layer //J. Fluid Mech. 2006. Vol. 564. P. 287-303.
3. Hill A. A., Straughan B. Poiseuille flow in a fluid overlying a porous medium // Ibid. 2008. Vol. 603. P. 137-149.
4. Hill A. A., Straughan B. Poiseuille flow in a fluid overlying a highly porous material // Advances in Water Resources. 2009. Vol. 32. P. 1609-1614.
5. Hill A. A., Straughan B. Global stability for thermal convection in a fluid overlying a highly porous material // Proc. Roy. Soc. London A. 2009. Vol. 465. P. 207-217.
6. Ghisalberti M., Nepf H. M. Mixing layers and coherent structures in vegetated aquatic flows // J. Geophys. Research. 2002. Vol. 107. P. 11.
7. White B. L., Nepf H. M. Shear instability and coherent structures in shallow flow adjacent to a porous layer // J. Fluid Mech. 2007. Vol. 593. P. 1 -32.
8. Auriault J. L. On the domain of validity of Brinkman’s equation // Transport Porous Media. 2009. Vol. 79, N. 2. P. 215-223.
Stability of steady flow above a saturated porous medium
D. T. Baydina, T. P. Lyubimova
Perm State University, Bukirev St., 15, 614990, Perm
In this research stability of steady fluid flow above porous media is studied concerning a problem of fluid stability above a bottom covered with vegetation. The system a only-fluid layer - porous media layer, saturated with the same fluid is described using a single approach (using equations of motion in porous media), and two variables (depending on vertical coordinate) permeability and porosity are introduced to model the presence of two different layers. A surface restricting the system at the bottom is assumed to be solid, upper surface - free and warp-free. Neutral curves of the stability and dependences of minimum Reynolds number R(k) and wave number of the most dangerous perturbations against ratio of porous layer thickness to entire thickness was obtained. It was also shown that dependence Re m (d) is nonmonotonous.
Keywords: stability, porous medium.