On one class of analytic solutions of the stationary axisymmetric convection Bénard–Maragoni viscous incompreeible fluid

Full Text

Введение. Не вызывает сомнений, что конвективное течение жидкости, связанное с неоднородным нагревом, является самым распространенным во Вселенной [1]. Систематические исследования конвективных течений вязкой несжимаемой жидкости датируются тысяча восемьсот восемьдесят восьмым годом, когда Томсон наблюдал образование сотовых структур в мыльной воде и мясном бульоне [2]. Это эмпирическое наблюдение легло в основу экспериментальной гидродинамики и регулярно используется для подтверждения или опровержения новых теоретических концепций, позволяющих количественно и качественно оценить основные характеристики течений. В работах Бенара [3–5] отражено продолжение экспериментальных исследований конвективных течений и описана структура течения, представляющая собой пчелиные соты. При интерпретации опыта Бенаром в качестве механизма образования конвекции при подогревании плоского горизонтального слоя указывалась связь вязкости жидкости с термокапиллярными силами, возникающими на свободной поверхности [3–5], а не подъемная сила, как предполагали другие исследватели [1, 2]. Несмотря на то, что изучение конвекции жидкости началось с плоских течений, не меньший интерес вызывают осесимметричные течения. Осесимметричная конвекция Бенара— Марангони очень важна при изучении ряда технологических процессов, например, при описании выращивания кристаллов методом Чохральского [6], в Сергей Николаевич Аристов (д.ф.-м.н.), главный научный сотрудник, лаб. гидродинамической устойчивости. Просвиряков Евгений Юрьевич (к.ф.-м.н.), докторант, доцент, каф. основ конструирования. 110 Об одном классе аналитических решений стационарной осесимметричной конвекции . . . астрофизике [7], геофизической гидродинамике [7–10], метеорологии и климатологии [7, 9, 10]. Точное решение уравнений, описывающих термокапиллярную конвекцию, впервые было найдено Бирихом [11]. Развитие и обобщение решения Бириха, интерпретацию физического смысла механизмов термокапиллярной конвекции можно проследить по работам [11–16]. Отметим, что до сих пор исчерпывающе не решена задача об описании осесимметричных термокапиллярных течений вблизи экстремальных значений температуры, в которых производная, вычисленная по радиальной координате, обращается в нуль. Построенные решения [11–16] не позволяют изучить анонсируемый класс течений, поэтому целью настоящей работы является построение точного решения системы Обербека—Буссинеска термокапиллярной конвекции, описывающего экстремальное температурное поле при числе Грасгофа, стремящемся к нулю. Найденное решение при необременительных ограничениях может быть использовано в качестве начального приближения для поиска устойчивых и неустойчивых решений осесимметричной конвекции Бенара—Марангони. 1. Постановка задачи. Запишем систему уравнений в цилиндрической системе координат, описывающих термокапиллярную конвекцию при установившемся осесимметричном движении вязкой несжимаемой жидкости (система уравнений Обербека—Буссинеска [17, 18]): 2 Vϕ ∂Vz ∂P Vr ∂Vr + Vz − =− + ν ∆Vr − 2 , Vr ∂r ∂z r ∂r r ∂Vϕ Vr Vϕ Vϕ ∂Vϕ + Vz + = ν ∆Vϕ − 2 , Vr ∂r ∂z r r ∂Vz ∂Vz ∂P Vr + Vz =− + gβT + ν∆Vz , ∂r ∂z ∂z ∂T 1 ∂ (rVr ) ∂Vz ∂T + Vz = χ∆T, + = 0. Vr ∂r ∂z r ∂r ∂z (1) В системе уравнений (1) введены следующие обозначения: Vr , Vϕ , Vz — радиальная, азимутальная и вертикальная скорости соответственно; P — отклонение давления от гидростатического, делённое на постоянную среднюю плотность жидкости ρ; T — отклонение от средней температуры; ν, χ — коэффициенты кинематической вязкости и температуропроводности жидкости ∂ ∂ ∂2 соответственно; ∆ = 1 ∂r r ∂r + ∂z 2 — оператор Лапласа, определённый в r цилиндрической системе координат. При достаточно малом нагреве решения системы (1) могут быть разложены в ряды Тейлора с различным числом удерживаемых членов для анализа течений [19, 20]. Решение системы (1) будем искать в следующем классе [19]: Vr = ru, Vϕ = rv, Vz = w, T = T0 + r2 r2 T11 , P = P0 + P11 . 2 2 (2) Класс решений (2) обобщает решения, полученные Бирихом [11], Карманом [17], Сидоровым [20], Линем [21]. Подставим класс (2) в (1), получим полиномиальные выражения, которые в силу линейной независимости степенного базиса формируют систему уравнений, распадающуюся на две под111 С. Н. А р и с т о в, Е. Ю. П р о с в и р я к о в системы: ∂u ∂2u − v 2 = −P11 + ν 2 , ∂z ∂z ∂T11 ∂ 2 T11 2uT11 + w =χ ; ∂z ∂z 2 ∂w ∂P11 = 0, = gβT11 , ∂z ∂z ∂v ∂2v 2uv + w = ν 2, ∂z ∂z u2 + w 2u + (3) ∂P0 ∂2w ∂T0 ∂ 2 T0 ∂w =− + ν 2 + gβT0 , w = χ T11 + . (4) ∂z ∂z ∂z ∂z ∂z 2 Запишем граничные условия, учитывая, что на верхней, свободной, границе задан локализованный источник тепла: w T = θ r 1− 2 R 2 . (5) Здесь θ — характерная разность температур, R — максимальное радиальное значение. Для основной системы (3) восьмого порядка граничные условия имеют следующий вид: при z = −h : v = w = при z = 0 : w = ∂w = 0; ∂z ∂v θ ∂u σ = 0, P11 = 0, T11 = − 2 , = − T11 , ∂z R ∂z η где σ = ∂α/∂T температурный коэффициент поверхностного натяжения α, η — коэффициент динамической вязкости [18]. Для вспомогательной системы третьего порядка граничные условия имеют следующий вид: при z = −h : T0 + T11 R2 θ = 0; при z = 0 : T0 = , P0 = A. 2 2 (6) Здесь A — атмосферное давление [13, 17, 18]. Граничное условие на твёрдой поверхности задается из условия, согласно которому имеет место нулевая температура на радиальной границе. Раскладывая в ряды Тейлора граничные условия, убеждаемся, что класс решений (2) удовлетворяет граничным условиям. 2. Решение редуцированной системы уравнений Обербека—Буссинеска. Для решения краевой задачи (3)–(6) приведем уравнения (3), (4) и граничные условия (5)–(6) к безразмерному виду. В качестве базиса выберем следующие размерные комплексы: h — вертикальный (поперечный) характерный размер; R — радиальный характерный размер; θ — температура; gβθh4 /(νR2 ) — вертикальная скорость; gβθh — давление, делённое на постоянную плотность. Используя выбранный базис, получим краевую задачу, записанную в безразмерной форме для систем (3) и (4): ∂w ∂P11 = 0, = T11 ; ∂z ∂z ∂v ∂2v Gr 2uv + w = 2, ∂z ∂z 2u + 112 ∂u ∂2u − v 2 = −P11 + 2 ; ∂z ∂z 2T ∂T11 1 ∂ 11 Gr 2uT11 + w = ; ∂z Pr ∂z 2 Gr u2 + w (7) Об одном классе аналитических решений стационарной осесимметричной конвекции . . . ∂w ∂T0 ∂P0 ∂2w 1 2 ∂ 2 T0 =− + δ 2 2 + T0 , Grw = δ T11 + ; ∂z ∂z ∂z ∂z Pr ∂z 2 ∂w = 0; при z = −1 : v = w = ∂z ∂v ∂u при z = 0 : w = = 0, P11 = 0, T11 = −1, = −MgT11 ; ∂z ∂z Grδ 2 w при z = −1 : T0 + (8) T11 1 = 0; при z = 0 : T0 = , P0 = 1. 2 2 (9) Здесь Gr = gβθh5 /(ν 2 R2 ) — осесимметричное число Грасгофа, записанное для вертикального и радиального характерных масштабов; Pr = ν/χ и Mg = = σ/(ρh2 gβ) — число Прандтля и Марангони соответственно; δ = h/R — безразмерный геометрический параметр. Отметим, что если осуществлять обезразмеривание с одним характерным масштабом, то есть h = R, то получим классическое число Грасгофа: Gr = gβθh3 /ν 2 [1, 2, 6–18, 20]. При числе Грасгофа, равном нулю (Gr = 0), краевая задача (7)–(9) имеет аналитическое решение: 1 z2 1 z − 1, P11 = 15 Mg − − z, 4 4 2 1 z4 z3 3 1 u = 15 Mg − − + Mgz + Mg − , 4 4! 3! 8 96 1 z5 z4 z2 3 1 w = 15 −2Mg + + 2 − 2Mg + − Mg + , 2 5! 4! 2 4 48 1 z3 z2 15 15 5 δ2 1 − Mg + δ2 + − Mg + Mgδ 2 − z+ , T0 = 15δ 2 4 3! 2! 8 2 2 8 2 1 z4 z3 15 15 5 δ2 z2 − Mg + 3δ 2 + − Mg + Mgδ 2 − + P0 = 45δ 2 4 4! 3! 8 2 2 8 2! 1 + − 2δ 2 Mg z + 1. 2 T11 = 15 Mg − (10) 3. Исследование полученных решений. Проанализируем знаки значений гидродинамических полей (10). Рассмотрим для начала вертикальную скорость Vz = w. Найденное решение принадлежит классу многочленов пятой степени. Отметим, что в силу граничных условий (8), (9) многочлен w можно факторизовать следующим образом [22]: w = z (z + 1)2 g = z(z + 1)2 1 Mg 2 Mg 1 1 3Mg − z + − z+ − . 16 4 2 24 48 4 Поскольку координата z определена на отрезке [−1; 0], то z(z + 1)2 0, следовательно, знак скорости будет определяться знаком квадратичного полинома g, коэффициенты которого непрерывно зависят от числа Марангони. При Mg ∈ [0; 1/36] скорость w принимает только отрицательные значения на области определения координаты (рис. 1), при Mg ∈ (1/36; 1/12) — как положительные, так и отрицательные значения, а при Mg 1/12 — положительные значения (рис. 1). Таким образом, в зависимости от числа Марангони на 113 С. Н. А р и с т о в, Е. Ю. П р о с в и р я к о в свободной (верхней) и твёрдой (нижней) границах скорость может одновременно иметь разные знаки (рис. 1). Полученный эффект смены знака вертикальной скорости на свободной границе при изменении числа Марангони наблюдался в экспериментальных и теоретических исследованиях [1, 11, 13, 18]. Знак радиальной скорости Vr = ru, очевидно, зависит от знака полинома u, который факторизуется следующим образом: 5 3 5Mg 1 2 5Mg − z − + z + 8 32 8 96 5Mg 1 3Mg 1 + + z+ − . 8 96 8 96 На рис. 2 изображено множество значений функции u, что в силу класса (2) и области определения переменной r ∈ [0; 1] с точностью до положительного множителя определяет множество значений радиальной скорости. При Mg ∈ [0; 1/36) ∪ (1/12; +∞) полином s имеет один корень в области определения, а при Mg ∈ [1/36; 1/12] — два корня. Данная локализация корней поясняет геометрию множества значений радиальной скорости. Используя решение (10), приведем (рис. 3) распределение температуры на твёрдой границе (z = −1), зависящей, как нетрудно видеть, только от числа Марангони. Отметим, что на твёрдой границе при Mg ∈ [0; 11/60] жидкость имеет отрицательную температуру, при Mg > 0 — положительную. Если рассматривать z ∈ (−1; 0), то знак и распределение температуры существенно зависят от величины δ, причем существуют такие значения числа Марангони, что параметр δ принимает значения, большие шести. Это обстоятельство подчеркивает важность учета и введения различных масштабных параметров, а не одного [19]. Отметим, что для обобщённого давления можно провести анализ, аналогичный температурному, и уже на твёрдой границе (z = −1) проявляется зависимость распределения знаков приведенного давления от числа Марангони Mg и параметра δ. На рис. 4–8 представлены изолинии обобщённой (модифицированной) функции тока [17], давления, отнесенного к плотности, и температуры при числе Марангони Mg = 0 и δ = 1/2. Для уточнения типа изолиний, изображенных на рис. 4 – 8, приведем алгоритм, позволяющий определить характер изолинии. Не ограничивая общности, рассмотрим алгоритм для температурного поля. Известно, что тип функции нескольких переменных (геометрически — поверхность) можно определить, анализируя знак гауссовой кривизны, который определяется знаками собственных значений ее матрицы Гессе [23–26]. Матрица Гессе, записанная для температуры, имеет следующий вид: u = (z + 1)s = (z + 1) H(T ) = ∂2T ∂r2 ∂2T ∂z∂r ∂2T ∂r∂z ∂2T ∂z 2 = T11 r ∂T11 ∂z r ∂T11 ∂z ∂ 2 T0 ∂z 2 . Изолинии температуры будут принадлежать к параболическому типу тогда, когда матрица Гессе имеет собственные значения одного знака, гиперболического — собственные числа разных знаков [23–26]. Воспользовавшись методами локализации собственных значений [22, 27, 28], получим, что изолинии температуры и других гидродинамических полей могут принадлежать только к гиперболическому типу. 114 Об одном классе аналитических решений стационарной осесимметричной конвекции . . . Рис. 1. Множество значений вертикальной скорости w при различных числах Марангони Mg (закрашенная область соответствует положительным значениям, незакрашенная — отрицательным) Рис. 2. Множество значений функции u при различных числах Марангони Mg (закрашенная область соответствует положительным значениям, незакрашенная — отрицательным) Рис. 3. Множество значений температуры на твёрдой границе при различных числах Марангони Mg (закрашенная область соответствует положительным значениям, незакрашенная — отрицательным) Рис. 4. Изолинии модифицированной функции тока при Mg = 10, δ = 1/2 (изолинии построены для значений, принадлежащих отрезку [0; 1, 4], с шагом, равным 0,1) Рис. 5. Изолинии обобщённого давления при Mg = 0,01, δ = 1/2 (изолинии построены для значений, принадлежащих отрезку [0, 25; 0, 9], с шагом, равным 0,05; для значений, принадлежащих отрезку [0, 9; 1] с шагом, равным 0,01) Рис. 6. Изолинии обобщённого давления при Mg = 10, δ = 1/2 (изолинии построены для значений, принадлежащих отрезку [−30; 40], с шагом, равным 5) 115 С. Н. А р и с т о в, Е. Ю. П р о с в и р я к о в Рис. 7. Изолинии температуры при Mg = = 0,01, δ = 1/2 (изолинии построены для значений, принадлежащих отрезку [−1; 0,4], с шагом, равным 0,05) Рис. 8. Изолинии температуры при Mg = 10, δ = 1/2 (изолинии построены для значений, принадлежащих отрезку [0; 60], с шагом, равным 5) На рис. 5, 6 и рис. 7, 8 показана зависимость направления теплового потока от величины числа Марангони. Пограничное значение числа Марангони вычисляется методами локализации корней степенных решений (10), которые были использованы выше [22, 27, 28]. Заключение. Построено новое точное решение системы уравнений Обербека—Буссинеска, описывающее осесимметричное термокапиллярное течение вязкой несжимаемой жидкости при локальном параболическом нагреве свободной поверхности. Рассмотрен случай, когда число Грасгофа стремится к нулю. Определены области постоянств знаков у множеств значений решений в зависимости от числа Марангони и безразмерного параметра δ. Найдено граничное значение числа Марангони, при котором происходит смена знака вертикальной и радиальной скоростей. Показано, что изолинии гидродинамических и температурного полей могут принадлежать только к гиперболическому типу. Отметим, что полученные решения можно обобщить в следующих случаях. Во-первых, можно получить точные решения для расширенного класса (2), совершая преобразование сдвига для давления и температуры. Во-вторых, найденное решение может служить стартовым (начальным) решением при использовании итерационных методов нахождения точных решений конвекции Бенара—Релея при числах Грасгофа, больших нуля.

About the authors

Sergey N Aristov

Institute of Continuous Media Mechanics UB RAS

Email: asn@icmm.ru
(Dr. Phys. & Math. Sci.), Chief Researcher, Lab. Hydrodynamical Instability 1, Academician Koroleva st., Perm, 614013, Russia

Evgeii Yu Prosviryakov

Kazan National Research Technical University named after A. N. Tupolev

Email: evgen_pros@mail.ru
(Ph. D. Phys. & Math.), Doctoral Candidate, Dept. of Bases of Construction 10, Karl Marx st., Kazan, 420111, Russia

References

Supplementary files