НЕЛИНЕЙНЫЕ РЕЖИМЫ СТАЦИОНАРНОЙ КОНВЕКЦИИ МАРАНГОНИ В ТОНКОЙ ПЛЕНКЕ ЖИДКОСТИ НА НАГРЕТОЙ ПОДЛОЖКЕ Текст научной статьи по специальности «Физика»

ВЕСТНИК ПЕРМСКОГО УНИВЕРСИТЕТА

2022

• ФИЗИКА •

Вып. 1

УДК 532.5.013.4 PACS 47.55.nb

Нелинейные режимы стационарной конвекции Марангони в тонкой пленке

А. Е. СамойловаЭ.В. Пермякова

1 Пермский государственный национальный исследовательский университет, Пермь, Россия

2 Институт механики сплошных сред ПФИЦ УрО РАН, Пермь, Россия

■ samoilova@psu.ru

1 evelina.v.permyakova@gmail. com@gmail. т

В работе исследуется поведение тонкой пленки жидкости, помещенной на нагретую твердую подложку низкой теплопроводности. В рамках длинноволнового приближения проведено прямое численное моделирование стационарной конвекции Марангони. Используется длинноволновая модель, допускающая деформации свободной поверхности пленки и предписывающая слабую теплоотдачу с поверхности в окружающую среду. Проведено сопоставление результатов численного моделирования с известными предсказаниями слабонелинейного анализа. Ожидаемое совпадение наблюдается в области небольшой надкритичности: при мягком возбуждении конвекции значение амплитуды движения хорошо описывается корневым законом. При жестком возбуждении конвекции показано существование гистерезиса в зависимости амплитуды движения от надкритичности. Профиль тонкой пленки вблизи порога имеет форму близкую к синусоидальной, с возбуждением двух вихрей противоположной закрутки под подъемом толщины. При отдалении от порога возникновения конвекции наблюдается возбуждение стационарных режимов с локальным максимумом в центре провала толщины, т.е. с появлением двух дополнительных вихрей. Величина локального максимума растет с надкритичностью. В ситуации мягкого возбуждения при больших надкритичностях устанавливается режим с удвоенным волновым числом, при этом в минимумах толщины образуются ультратонкие слои. В ситуации жесткого возбуждения при большой надкритичности система эволюционирует к разрыву. При этом в момент разрыва серия локальных максимумов толщины напоминает пальцеобразную неустойчивость. Показано, что при определенных параметрах система эволюционирует от установившегося режима с одной величиной локального максимума к режиму с большей такой величиной через нестационарный режим в виде бегущих волн.

Ключевые слова: конвекция Марангони; тонкая пленка; численный эксперимент

Поступила в редакцию 21.09.21; после рецензии 04.03.22; принята к опубликованию 06.03.22

Nonlinear regimes of the steady Marangoni convection in a thin film heated from below

A.E. Samoilova1, E. V. Permyakova1,2

1 Perm State University, Perm, Russia

2 Institute of Continuos Media Mechanics UB RAS, Perm, Russia f samoilova@psu.ru

* evelina.v.permyakova@gmail.com @gmail.ru

We investigate the dynamics of a thin liquid film that is placed atop a heated substrate of a very low thermal conductivity. The direct numerical simulation of the stationary long-wave Marangoni

© Самойлова А.Е., Пермякова Э.В., 2022

жидкости на нагретой подложке

распространяется на условиях лицензии

Creative Commons Attribution 4.0 International (CC BY 4.0).

instability is performed within the system of coupled partial differential equations. These equations were previously derived within the lubrication approximation; they describe the evolution of film thickness and fluid temperature. We compare our results with previous results of the weakly-nonlinear analysis. A good qualitative agreement is observed for values of the Marangoni number near the instability threshold. In the case of supercritical excitation, our results for the amplitudes are described by the square root dependence on the supercriticality. In the case of subcritical excitation, we found the hysteresis. In the vicinity of the instability threshold the film interface is nearly sinusoidal, with two vortices under elevated region. As the supercriticality increases, the nonlinear stationary regimes arise, with a local elevation at the global minimum and with two additional vortices, correspondingly. In a case of subcritical excitation, relatively high supercriticality results in that these regimes evolve into film rupture via the emergence of secondary humps. For the supercritical excitation, stationary regime with doubled wavenumber occurs. We also revealed the transition through the traveling wave between two nonlinear regimes with different high of local humps.

Keywords: Marangoni convection; thin film; numerical simulation

Received 21.09.21; revised 04.03.22; accepted 06.03.22 doi: 10.17072/1994-3598-2022-1-49-57 1. Введение

В данной работе численно исследуется поведение тонкой пленки жидкости, помещенной на нагретую твердую подложку низкой теплопроводности [1-3]. В такой системе при определенных условиях может возникать стационарная конвекция, главный механизм которой - термокапиллярный эффект, т.е. зависимость поверхностного натяжения жидкости от температуры [4-6]. Характеристики конвективной неустойчивости при этом сильно зависят от свойств системы. Так, при достаточно сильном поверхностном натяжении жидкости, которое препятствует деформациям свободной поверхности, возникает т.н. ячеистая или коротковолновая конвективная неустойчивость [7]. Характерный размер конвективных структур (ячеек) при этом сравним с характерным пространственным масштабом системы - толщиной слоя жидкости. Если же недеформируемая поверхность очень слабо отдает тепло во внешнюю среду, то неустойчивость возникает в виде сильно вытянутых структур [4]. Такой тип неустойчивости называется длинноволновым; он возникает также в ситуации, когда свободная поверхность может деформироваться [8, 9].

Основное внимание в работе уделяется длинноволновой конвективной неустойчивости вследствие того, что рассматриваемая система - тонкая пленка жидкости - обладает подходящими пространственными характеристиками: горизонтальный масштаб в ней намного больше, чем толщина слоя. В этом случае для описания динамики системы удобно использовать не полные гидродинамические уравнения, а специальную длинноволновую модель [1]. В данной работе мы используем длинноволновую модель, полученную в [10]. Особенность этой модели состоит в том, что авторам удалось скомбинировать условия для возникновения

длинноволновой неустойчивости: слабую теплоотдачу со свободной поверхности и конечную деформируемость. В результате были получены новые интересные результаты. В частности, предсказано существования мягкого возбуждения неустойчивости. Прежде полагалось, что развитие неустойчивости Марангони в тонкой пленке всегда приводит к разрыву вследствие жесткого возбуждения конвекции [11, 12]. Авторам работы [10] удалось показать существование стабильных стационарных конвективных структур в тонкой пленке жидкости.

В настоящей работе было проведено прямое численное моделирование динамики тонкой пленки в рамках уравнений в частных производных, описывающих эволюцию толщины пленки и температуры жидкости. Результаты сопоставлялись с предсказаниями слабонелинейного анализа из [10]; также исследована динамика системы вдалеке от порога возникновения неустойчивости.

Настоящая статья организована следующим образом. В главе 2 приводится постановка задачи. Показан вывод длинноволновых уравнений из полных гидродинамических уравнений и граничных условий. Глава 3 посвящена исследованию плоских стационарных структур неустойчивости Марангони. Здесь обсуждается методика вычислений, а также приведены основные результаты. В главе 4 подведены итоги и обозначены планы будущего исследования.

2. Постановка задачи

2.1. Уравнения и граничные условия

В работе рассматривается поведение тонкой пленки жидкости, размещенной на твердой нагретой подложке (см. рис. 1). Предполагается, что подложка обладает значительно более низкой теплопроводностью, чем жидкость. Таким образом,

нижняя граница слоя жидкости считается теплоизолированной, на ней выполняется условие фиксированного теплового потока. Верхняя граница тонкой пленки является свободной; тепловой поток на этой границе по закону Ньютона уходит в окружающее пространство.

2, > н 1 И(х,у,()

V/ Ь ^^ к->! ' У\

0

л:' /; > у;

— + V-УТ = У2Т, Шу V = 0 , дt

которые дополняются граничными условиями [13]:

дТ

z = 0: V = 0, — = -1,

дz

дh дТ

г = h(х,у, t) : — = ^ — и -УН, — = —Б1Т,

дt дп

(2)

Рис. 1. Схематическое изображение тонкой пленки жидкости на твердой нагретой подложке (не в масштабе, предполагается Н << Ь). Горизонтальными пунктирными линиями обозначено равновесное положение поверхности пленки. Вертикальные пунктирные линии ограничивают произвольный объём, в нашей модели пленка представляется бесконечным слоем жидкости

Основным механизмом массопереноса в такой системе является конвекция Марангони, вызванная термокапиллярным эффектом. В этом случае поверхностное натяжение жидкости уменьшается с ростом температуры по линейному закону, характерному для большого круга жидкостей: а = а0 -атТ. Здесь и далее Т обозначает температуру жидкости, отсчитываемую от температуры окружающего газа, которая считается неизменной. Конвективным переносом, вызванным эффектом теплового расширения, в данной системе можно пренебречь ввиду малой толщины пленки.

Сила тяжести считается нормальной, однако предполагается, что свободная поверхность пленки может значительно деформироваться вследствие развития конвекции. Для очень тонких пленок (~0.1-0.01 мм) это условие выполнимо даже для воды, которая обладает весьма высоким поверхностным натяжением.

В качестве примера реальной системы подходит слой силиконового масла с кинематической вязкостью V = 100 сСт толщиной Н = 0.1 мм. В качестве теплоизолирующей подложки в таком случае может выступать плита из органического стекла. Дальнейшие рассуждения, однако, будут проводиться в безразмерном виде.

Конвекция Марангони в данной системе описывается уравнением Навье-Стокса, уравнениями переноса тепла и сохранения массы [6]:

т(0-о-п = — Мат(0-У(т|г_А), п■ о■ п = р — СаК.

Здесь и далее мы примем обозначение для нормальной и касательной компонент скорости жидкости v = (^ м>). Единичный вектор нормали п к свободной поверхности тонкой пленки определяется как п = (к — УН)/. Кривизна свободной поверхности вычисляется согласно соотношению К = - ^ п. Тензор вязких напряжений обозначен с, а пара ортогональных друг другу касательных векторов к искривленной поверхности пленки - т(,).

Задача (1)-(2) содержит следующие безразмерные комбинации физических параметров:

Рг = ■

Оа =

X

И

X

Б1 =

дИ

— + V■Vv \ = -Ур + У2V-Оак, Рг \дt 1

(1)

ат АИ „ ап И

Ма = —-, Са = .

Это числа Прандтля, Галилея, Био, Марангони и параметр капиллярности, соответственно.

Как видно, в данной системе может реализоваться состояние механического равновесия, соответствующее теплопроводному прогреву пленки:

V = 0 , Н = 1, р = Оа(1 — 2) , Т = 1 — г + БГ1. (3)

При превышении порогового значения числа Марангони Мас механическое равновесие сменяется конвективным движением. Однако прежде чем переходить к изучению конвективного движения, произведем некоторое упрощение используемой математической модели.

2.2. Длинноволновое приближение

Действительно, изучение конвекции в рамках задачи (1)-(2) представляется весьма трудоемким. Математическую модель удается значительно упростить, сосредоточившись на изучении длинноволновой конвекции, соответствующей крупномасштабным движениям пленки. Принимая отношение толщины пленки Н к характерному горизонтальному масштабу движений Ь за малый параметр е = Н/Ь, можно провести длинноволновое разложение в рамках уравнений (1)-(2). В результате, в ведущем порядке по малому е тепло- мас-

V

к

соперенос в тонкой пленке можно описать двумя уравнениями в частных производных [10]:

f = V- j, j = VP + Mah- Vf, ôt 3 2

h °© = V-( hV©)- 1(Vh)2+ j-Vf +

V

^h4

h3

—VP + — MaVf |. 8 6 1

(4)

(5)

Здесь Р = ОаА - СУ^ обозначает давление под искривленной поверхностью, / = 0 - h соответствует возмущению температуры на искривленной поверхности. В меру малости параметра е здесь также перемасштабированы число Био и параметр капиллярности: Са = е 2 С, Ы = е 2в.

Уравнения (4)-(5) описывают крупномасштабную эволюцию тонкой пленки с помощью двух функций: толщины пленки Щх,у,() и отклонения температуры жидкости 0(х,у,/) от равновесного значения (3). Состоянию механического равновесия (3) отвечают следующие решения уравнений (4)-(5): h = 1, 0 = 1.

Линейный анализ устойчивости равновесия относительно нормальных возмущений ~ ехр(Х/ + 1кхх+1куу) позволяет определить порог возникновения конвекции. В работе [10] получено уравнение для инкремента X:

X2 + X

p + k 2I 1 +

Ga + Ck2 - Ma

(Ga + Ck2) 3

k2(P + k2) - Mak'

(72 + Ga + Ck2) 144

= 0,

где к2 = кх2+ку2. Для монотонной моды несложно определить границу возникновения неустойчивости, положив X = 0. Нейтральная кривая определяется выражением, полученным в [10]:

Ma(k ) =

48(P + k 2)(Ga + Ck2) k2 (72 + Ga + Ck2)

(6)

Минимум функции Ма(к) отвечает критическим, т.е. наиболее опасным возмущениям. Характерный размер этих возмущений определяется критическим волновым числом:

k 2 =

PGaC + ,/ 72 pGaC {Ga + 72 -PC ) C (72 -PC )

(7)

Пороговое значение Mac несложно получить, подставив выражение (7) в формулу (6). Следует отметить, что конечные значения kc можно получить только при вС < 72. Это связано с ограничениями длинноволнового приближения (см. ниже). При вС > 72 в (7) получаем kc ^ да, что на самом деле соответствует коротковолновым возмущениям.

При этом (6) дает пороговое значение Mac = 48, которое совпадает с результатом, полученным Пирсоном в [4] для случая недеформированной свободной границы.

Сопоставление результатов линейного анализа в рамках полной системы гидродинамических уравнений (1)-(2) и уравнений в приближении тонких пленок (4)-(5) было проведено в работе [14]. Показано очень хорошее количественное совпадение результатов не только для длинноволновых возмущений, но и для ячеистых (с конечным значением волнового числа вплоть до 1). Однако область применения длинноволнового приближения оказывается ограничена значениями, которые могут принимать два безразмерных параметра - число Био Bi и параметр капиллярности Ca. Действительно, при выводе модельных уравнений (4)-(5) используются предположения: Bi ~ O(s 2) и Ca ~ O(s - 2). Первое требование соответствуют слабой теплоотдаче с поверхности пленки, тогда именно длинноволновые возмущения возбуждают конвекцию Марангони. Второе условие отвечает такому поверхностному натяжению, которое, с одной стороны, допускает деформации свободной поверхности, а с другой, не дает развиться коротковолновым возмущениям.

Во всех дальнейших вычислениях принято С = 1. В действительности выбор значения перемасштабированного параметра капиллярности С влияет лишь на значение малого параметра s, но не влияет на результаты вычислений.

3. Плоские структуры стационарной конвекции Марангони

3.1. Численное моделирование

В данном разделе мы ограничиваемся рассмотрением плоской задачи. В работе [10] показано, что это не влияет на порог возникновения конвекции, но, безусловно, ограничивает круг возможных наблюдаемых конвективных режимов.

Будем полагать, что рассматриваемая система однородна вдоль оси OY. Система уравнений (4)-(5) в таком случае упростится:

dh_ hL OP

ôt ôx I 3 Ox

h2 ôf .__

2 ôx

— = — | —— + Ma—P = Gah -C— , (8)

,ô2 h

ôx2

hô©=Afhô©ï im -p/+

ôt ôx l ôx ) 2l ôx )

h ôP , h2 ôf I ôf ô I h ôP h\, ôf ^

--+ Ma--— I —+—I--+ —Ma —

3 ôx 2 ôx ôx ôx 8 ôx 6 ôx

(9)

Данная система решается эволюционным методом с применением метода конечных разностей. В качестве расчетной области выбран промежуток [0, Ь] на оси ОХ, Ь = 2п/ к. Значение волнового

числа k выбирается соответствующим наиболее опасным возмущениям при данных значениях параметров Ga, C, в. На границах ставятся периодические граничные условия: функции h(x,t) и 0(x,t), а также их производные должны сшиваться при x = 0 и x = L.

Численное моделирование производится при помощи явной схемы. Для интегрирования по времени используется метод Рунге-Кутты пятого порядка. Шаг по времени менялся в пределах Дт =(0.01 ^ 0.001)(Дх)4, чтобы удовлетворить сходимость явной схемы. Шаг по пространственной координате определялся как Дх = L/N. Количество узлов сетки N варьировалось от 100 до 200 в зависимости от наблюдаемых режимов конвекции. Для разрешения и сходимости схемы при вычислениях с небольшой надкритичностью, а также в отсутствие сложной нелинейной динамики (см. рис.2, 6) приемлемым оказывается N = 100. При моделировании сложных четырехвихревых режимов, перехода к нестационарным режимам (см. рис. 4, 5) число узлов выбиралось равным 150 или 200. Аппроксимация пространственных производных производилась стандартными центральными разностями. Важным этапом является задание начальных условий. Малое отклонение от механического равновесия задавалось в виде h(x,0) = 1+0.01cos(kcx). Для исследования динамики системы на неустойчивой ветви решения в качестве начальных условий использовались данные вычислений для установившегося режима с конечной деформацией пленки. Алгоритм вычислений был реализован на языке Fortran90, на котором авторами была написана и отлажена программа.

3.2. Результаты при Ga = 1, в = 4

Для расчетов выбраны несколько значений безразмерных параметров. При Ga = 1, в = 4 конвекция возбуждается при Ma >Mac = 5.75. По линейной теории критические возмущения характеризуются волновым числом kc = 1.5. Слабонелинейный анализ предсказывает, что конвекция возбуждается через прямую вилочную бифуркацию при данных значениях параметров (см. рис.4 в работе [10]). Также в [10] исследована конкуренция различных конвективных паттернов, которые могут возникать в данной системе: валов, квадратов, шестиугольников. При Ga = 1, в = 4 вблизи порога возникновения конвекции Mac должны наблюдаться устойчивые валы с характерным размером 2л/ kc.

На рис. 2 изображены профили толщины пленки h(x) и температуры 0(x) после установления стационарной конвекции Марангони. Вблизи порога возникновения конвекции (сплошные линии, Ma = 5.8) профили толщины и температуры хорошо описываются решением линейной задачи h(x) = 1 + Acos(kx), 00(x) = 1 + Bcos(kcx). Такой

режим мы будем называть двухвихревым по количеству вихрей противоположной закрутки, образующихся в жидкости под искривленной поверхностью (на одном пространственном периоде Е). При отдалении от порога (штриховые линии, Ma = 6.5) профили значительно деформируются: в центре расчетной области температура жидкости проседает, а минимум функции h(x) при этом расширяется в полочку.

Рис. 2. Профили толщины пленки h и температуры 0 при Ga = 1, в = 4. Сплошные линии - Ma = 5.8, штриховые линии -Ma = 6.5

Рис. 3. Амплитудная кривая стационарной конвекции при Ga = 1, в = 4. Штриховая линия соответствует корневому закону, полученному из слабонелинейного анализа

Итоговая кривая зависимости hmax(Ma) представлена на рис. 3. На ней квадратами отмечены точки, полученные численно, а штриховая линия соответствует корневому закону hшax = 1 + K(Ma -Mac)1/2, где коэффициент K известен из слабонелинейного анализа [10].

Как видно, уже при Ma > 6.2 наблюдается расхождение численных результатов с предсказаниями слабонелинейного анализа. Действительно, с ростом числа Марангони поведение тонкой пленки становится существенно нелинейным (см. рис. 4, 5).

Как показано выше, при значительной надкри-тичности синусоидальный профиль деформируется: провал в центре расчетной области расширяется (см. рис. 2). Дальнейшее увеличение надкритичности приводит к формированию стационарного режима с половинным пространственным периодом или характеризующимся удвоенным волновым числом 2кс.

Ранее появление такого режима было показано в [10], однако он полагался неустойчивым. Для случая идеально теплопроводной подложки его возникнование связывалось с пальцеобразными структурами при образовании разрыва пленки [11].

Рис. 4. Нелинейная динамика пленки жидкости при Ga = 1, в = 4. Сплошные линии — Ма = 7, штриховые линии — Ма = 8: а — профили установившихся конвективных течений; б — эволюция системы к ним На рис. 4 приведены профили толщины пленки при Ма = 7 (сплошные линии) и Ма = 8 (штриховые линии), а также изображено установление стационарных режимов. Как видно, при Ма = 7 формируется локальный максимум при Ь/2. Ранее в [10] было показано, что при этом двухвихревой стационарный режим сменяется четырехвихревым: в центре области под прогибом пленки образуются два слабых конвективных вихря с противоположной закруткой. При Ма = 8 максимум при Ь/2 ста-

новится глобальным, мощность четырех вихрей -сопоставимой. Стоит отметить, что ранее в [10] четырехвихревой режим полагался неустойчивым. Однако в наших расчетах этот режим сохраняется и при больших временах вычислений. Важно, что при этом минимум толщины пленки составляет ~0.03 от равновесной толщины пленки. Можно полагать, что формируется не разрыв, а ультратонкая пленка жидкости.

На рис. 5 изображена динамика пленки при промежуточном значении числа Марангони, Ma = 7.3. Здесь система находится между двумя притягивающими стационарными режимами: с небольшим подъемом в центре, как при Ma = 7, и развитым четырехвихревым режимом, как при Ma = 8. Как видно из графика hmax(t), система переходит из одного установившегося состояния в другое, при этом в переходном режиме реализуются бегущие структуры.

3.3. Результаты при Ga = 2, в = 2

При Ga = 2, в = 2 линейная теории предсказывает возникновение конвекции при Ma >Mac = 5.05 с характерным пространственным масштабом 2п/ kc, kc = 1.45. При таком наборе параметров система ведет себя принципиально иным образом: динамика вблизи точки бифуркации описывается обратной вилочной бифуркацией. Ранее в работах [11, 12] было показано, что возникновение стационарной неустойчивости Марангони через обратную вилочную бифуркацию приводит к разрыву пленки. Однако в указанных работах подложка обладала противоположными тепловыми свойствами: она полагалась хорошо проводящей тепло. В нашей работе подложка обладает низкой теплопроводностью (по сравнению с жидкостью), что играет важную роль в образовании т.н. сухих пятен, т.е. разрыве пленки.

L/2 х а

L/2 х б

L 0 10 20 30 40 50 60 70 80 t

в

Рис. 5. Динамика тонкой пленки при Ма = 7.3. Профили толщины пленки, соответствующие двум установившимся режимам: а — при t = 10 (сплошная линия) и / = 50 (штриховая линия); б — переходному режиму с бегущими структурами при t = 30, 32, 35 (сплошная, штриховая и штрих-пунктирная линии); в — эволюция максимума и минимума толщины пленки со временем

Рис. 7. Амплитудная кривая стационарной конвекции при Ga = 2, в = 2. Штриховая линия соответствует неустойчивой ветви решения, возникающей через обратную вилочную бифуркацию. Стрелками обозначены переходы между устойчивыми решениями

На рис. 7 изображена амплитудная кривая, т.е. зависимость hmax(Ma). Как видно, при критическом значении числа Марангони Mac неустойчивость равновесия возникает жестким образом: формируется конвективное течение с конечной амплитудой. Важно, что при этом не формируется разрыв слоя, как описано в предыдущих исследованиях для теплоизолированной подложки [11, 12]. На рис. 7 можно также видеть характерный для жесткого возбуждения конвекции гистерезис. Пороговое

значение числа Марангони Ma* = 4.9 оказывается ниже, если в качестве начального условия в вычислениях использовать профиль установившегося течения.

Рис. 6. Профили толщины пленки h и температуры 0 при Ga = 2, в = 2. Сплошные линии - Ma = 5.1, штриховые линии -Ma = 5.6

Характерный вид конвективных структур при небольших надкритичностях представлен на рис. 6. Как и в предыдущем разделе, при небольшой надкритичности профиль толщины пленки близок к синусоидальной форме. А с увеличением Ma в нижней части профиля возникает локальный максимум (здесь он расположен вблизи границ расчетной области, а не в середине, как на рис. 4, 5).

(а)

(б)

Рис. 8. Динамика тонкой пленки при Ma = 5.7. Профили толщины пленки, соответствующие: а - двум установившимся режимам: при t = 15 (сплошная линия) и t = 36 (штриховая линия); б - переходному режиму с бегущими структурами при t = 25, 27, 30 (сплошная, штриховая и штрих-пунктирная линии)

Как и для параметров, рассмотренных в разделе 3.2, динамика тонкой пленки становится существенно нелинейной при отдалении от порогового значения Mac. На рис. 8 показана динамика тонкой пленки при Ma = 5.7. Из положения механического равновесия система переходит в установившийся режим с локальным максимумом. Однако этот режим сменяется режимом с профилем, сдвинутым на L/2. Смена стационарных режимов происходит через образование бегущих конвективных структур. В отличие от случая мягкого возбуждения, рассмотренного в предыдущем разделе, здесь второй стационарный режим приводит в конечном итоге к образованию разрыва пленки.

4. Заключение

В данной работе проведено прямое численное моделирование длинноволновой неустойчивости Марангони. Изучена нелинейная динамика тонкой пленки жидкости на нагретой твердой подложке, плохо проводящей тепло.

Рассмотрены два набора параметров задачи, соответствующих различным толщинам слоя и теплоотдачи во внешнюю среду. При одном наборе параметров проведенный ранее слабонелинейный анализ предсказывает возникновение стационарной конвекции Марангони через прямую вилочную бифуркацию. При втором наборе параметров предсказано жесткое возбуждение конвекции, т.е. через обратную вилочную бифуркацию. В обоих случаях проведено сопоставление результатов численного счета и аналитических данных слабонелинейного анализа. Значение порога возникновения конвекции совпадает с точностью до

ошибки вычислений. Качественное поведение системы неплохо описывается аналитическими кривыми вблизи порога: для мягкого возбуждения амплитудная кривая с хорошей точностью совпадает с аналитической кривой, для случая жесткого возбуждения наблюдается гистерезис. Однако довольно быстро динамика тонкой пленки становится существенно нелинейной. Причем, когда профили толщины пленки ещё вполне описываются линейной теорией при небольших надкритично-стях, распределение температуры сильно отклоняется от предсказаний линейной теории.

Подтверждено существование четырехвихре-вых режимов, т.е. режимов не синусоидального профиля, а с локальным подъемом в центре глобального минимума профиля пленки. Обнаружено, что при определенной надкритичности система со временем переходит из одного установившегося режима в другой с явным удвоением волнового числа, т.е. 2k. Этот переход реализуется через нестационарный режим в виде бегущих волн несинусоидального профиля. В случае мягкого возбуждения нелинейный режим с 2k образует провалы по толщине почти до самой подложки. Однако это не приводит к разрыву, но образуется ультратонкая пленка ~0.03 толщин невозмущенной пленки. В ситуации жесткого возбуждения при небольших надкритичностях разрыв пленки не наблюдается, образуются стабильные стационарные конвективные структуры. С повышением значения числа Марангони реализуется сценарий образования разрыва через пальцеобразную неустойчивость [11].

В заключение отметим, что в теории гидродинамической неустойчивости немало примеров, когда добавление третьего измерения в задачу существенно влияет на результаты. В рассматриваемой системе было предсказано множество устойчивых трехмерных конвективных паттернов разной симметрии (квадраты, гексагоны, ромбы) [10, 15]. Поэтому численное моделирование динамики трехмерной пленки представляется важным этапом будущих исследований.

Авторы благодарят А.А. Алабужева за плодотворные дискуссии и полезные комментарии. Исследование выполнено за счет гранта Российского научного фонда (проект № 20-71-00070).

Список литературы

1. Oron A., Davis S. H., Bankoff S. G. Long-scale évolution of thin liquid films. Rev. Mod. Phys, 1997. V. 69. P. 931-980.

2. Craster R. V., Matar O. K. Dynamics and stability of thin liquid films. Rev. Mod. Phys, 2009. V. 81. P. 1131-1198.

3. Холпанов Л. П., Шкадов В. Я. Гидродинамика и теплообмен с поверхностью раздела. М. : Наука, 1990. 271 с.

4. Pearson J. R. A. On convection cells induced by surface tension. J. Fluid Mech, 1958. V 4. P. 489500.

5. Davis S.H. Thermocapillary instabilities. Annu. Rev. Fluid Mech.1987. V. 19. P. 403-435.

6. Гершуни Г. З., Жуховицкий Е. М. Конвективная устойчивость несжимаемой жидкости. М.: Наука, 1972. 392 с.

7. Benard H. Les tourbillons cellulaires dans une nappe liquid. Rev. Gen. Sci. Pures Appl, 1900. V.11. P. 1261-1271.

8. Scriven L. E., Sternling С. V. On cellular convection driven by surface-tension gradients: effects of mean surface tension and surface viscosity. J. Fluid Mech, 1964. V 19. P. 321-332

9. Smith K. A. On convective instability induced by surface-tension gradient. J. Fluid Mech, 1966. V. 24. P. 401-414

10. Shklyaev S., Alabuzhev A. A., Khenner M. Longwave Marangoni convection in a thin film heated from below. Phys. Rev. E. 2012. V. 85, 016328

11. Oron A. Nonlinear dynamics of three-dimensional long-wave Marangoni instability in thin liquid films, Phys. Fluids. 2000. V. 12. P. 1633-1645

12. VanHook S. J., Schatz M. F., Swift J. B., McCor-mick W. D., and Swinney H. L. Long-wavelength surface-tension-driven Benard convection: Experiment and theory. J. Fluid Mech, 1997. V. 345. P. 45-78.

13. Ландау Л.Д., Лифшиц Е.М. Теоретическая физика. Гидродинамика. M.: Наука, 1986. 736 с.

14. Samoilova A. E., Lobov N. I. On the oscillatory Marangoni instability in a thin film heated from below. Physics of Fluids. 2014. V. 26, 064101

15. Samoilova A. E., Nepomnyashchy A. Marangoni patterns on a rhombic lattice in a thin film heated from below. Physics of Fluids. 2021. V. 33, 014101.

References

1. Oron A., Davis S. H., Bankoff S. G. Long-scale evolution of thin liquid films. Rev. Mod. Phys, 1997, vol. 69, pp. 931-980.

2. Craster R. V., Matar O. K. Dynamics and stability of thin liquid films. Rev. Mod. Phys, 2009, vol. 81, pp. 1131-1198.

3. Holpanov, L. P., Shkadov, V. Ya. Hydrodynamics and heat and mass transfer with interface, Moscow: Nauka, 1990. 271 p. [in Russian]

4. Pearson J. R. A. On convection cells induced by surface tension. J. Fluid Mech, 1958, vol. 4, pp. 489-500.

5. Davis S.H. Thermocapillary instabilities. Annu. Rev. Fluid Mech, 1987, vol. 19, pp. 403-435.

6. Gershuni G. Z., Zhukhovitskii E. M. Convective stability of incompressible fluids. Jerusalem, Israel: Keter Publishing House, 1976, 330 p.

7. Bénard H. Les tourbillons cellulaires dans une nappe liquid. Rev. Gén. Sci. Pures Appl, 1900, vol. 11, pp. 1261-1271.

8. Scriven L. E., Sternling C. V. On cellular convection driven by surface-tension gradients: effects of mean surface tension and surface viscosity. J. Fluid Mech, 1964, vol. 19, pp. 321-332

9. Smith K. A. On convective instability induced by surface-tension gradient. J. Fluid Mech, 1966, vol. 24, pp. 401-414

10. Shklyaev S., Alabuzhev A. A., Khenner M. Longwave Marangoni convection in a thin film heated from below. Phys. Rev. E, 2012, vol. 85, 016328

11. Oron A. Nonlinear dynamics of three-dimensional long-wave Marangoni instability in thin liquid films. Phys. Fluids, 2000, vol. 12, pp. 1633-1645.

12. VanHook S. J., Schatz M. F., Swift J. B., McCor-mick W. D., Swinney H. L. Long-wavelength surface-tension-driven Benard convection: Experiment and theory. J. Fluid Mech, 1997, vol. 345, pp. 45-78.

13. Landau L. D. and Lifshitz E. M. Fluid Mechanics. London: Pergamon Press, 1987. 536 p.

14. Samoilova A. E., Lobov N. I. On the oscillatory Marangoni instability in a thin film heated from below. Phys. Fluids, 2014, vol. 26, 064101

15. Samoilova A. E., Nepomnyashchy A. Marangoni patterns on a rhombic lattice in a thin film heated from below. Phys. Fluids, 2021, vol. 33, 014101.

Просьба ссылаться на эту статью в русскоязычных источниках следующим образом: Самойлова А.Е., Пермякова Э.В. Нелинейные режимы стационарной конвекции Марангони в тонкой пленке жидкости на нагретой подложке // Вестник Пермского университета. Физика. 2022. № 1. С. 49-57. doi: 10.17072/1994-3598-2022-1-49-57

Please cite this article in English as:

Samoilova A.E., Permyakova E.V. Nonlinear regimes of the steady Marangoni convection in a thin film heated from below. Bulletin of Perm University. Physics, 2022, no. 1, pp. 49-57. doi: 10.17072/1994-3598-2022-1-4957

Сведения об авторах

1. Анна Евгеньевна Самойлова, канд. физ.-мат. наук, научный сотрудник кафедры теоретической физики, Пермский государственный национальный исследовательский университет, ул. Букирева, 15, Пермь, 614990

2. Эвелина Владимировна Пермякова, студентка бакалавриата физического факультета, Пермский государственный национальный исследовательский университет, ул. Букирева, 15, Пермь, 614990

Author information

1. Anna E. Samoilova, Candidate of Physical and Mathematical Sciences, Research Fellow at Theoretical Physics department, Perm State University, Bukirev str. 15, 614990, Perm, Russia

2. Evelina V. Permyakova, Bachelor Degree Student of Physical Faculty, Perm State University, Bukirev str. 15, 614990, Perm, Russia