Introduction to the generalized theory of non-equilibrium Cahn-Hilliard phase transitions (Thermodynamic problems in continuum mechanics)

Full Text

Введение [Возникновение конвективных течений и их развитие от регулярных форм с последующим переходом к нерегулярным - турбулентным течениям привлекают внимание тем, что они являются ответственными за эффективность многих технологических процессов тепломассопереноса [1, 2]. Такие технологические процессы являются базовыми в химической, нефтехимической, энергетической, металлургической и других отраслях промышленности. Конвективные течения возникают в жидкостях и газах в гравитационном поле при наличии пространственной неоднородности плотности, создаваемой неоднородностью температуры и концентрации компонентов, возникающих в ходе, например, химических реакций или других причин [3-5]. С увеличением разности температур покоящаяся жидкость теряет устойчивость, что затем приводит к возникновению конвективного течения. Дальнейшее увеличение разности температур приводит к неустойчивости первичного конвективного течения, а гидродинамический кризис приводит к кризису теплопередачи. Вторичное конвективное течение, формирующееся при этом далее также может стать неустойчивым, так что в результате последовательности таких кризисов может развиться турбулентное течение [1, 2]. Этим вопросам посвящена обширная литература, в частности монографии [6-13]. Постоянный интерес к подобного рода задачам, в частности, к задаче о неустойчивости Рэлея-Бенара, поддерживается из-за новых практических приложений. Так, в [13, 14] отмечается, что за пределами «чистой» гидродинамики формирование структур, близких к пространственно-периодическим, наблюдается при росте кристаллов, распространении фронтов затвердевания, электрохимических неустойчивостях нематических жидких кристаллов, химических реакционно-диффузионных процессах, автокаталитических реакциях, изгибе тонких пластин и оболочек и многих других. Кроме того, считается, что гидродинамический подход, базирующийся на модели Буссинеска, неприменим в условиях малой гравитации [15-18]. Развитие многих приложений требует более тонкого теоретического анализа, в частности, установления критериев различения термогравитационной неустойчивости Рэлея- Бенара [13, 14] и термокапиллярной неустойчивости Бенара-Марангони [19]. Так, в [1], анализируя экспериментальные результаты, авторы пришли к заключению, что в тонких слоях вертикальный перепад температуры мал и поэтому рэлеевская неустойчивость не реализуется, а реализуется неустойчивость Марангони. В слоях больших толщин вклад неустойчивости Марангони 438 Введение в обобщенную теорию неравновесных фазовых переходов. . . мал, а неустойчивость развивается по механизму Рэлея. Конкуренция двух механизмов наблюдается на масштабе промежуточных толщин.]1 [Интерес к задачам о совместном движении двух и более жидких сред, контактирующих вдоль некоторой поверхности, и их устойчивости обусловлен возможностью использования их в технических приложениях. На возникновение и характер конвективных течений в таких системах оказывают влияние различные факторы. В зависимости от размеров системы, соотношений параметров контактирующих сред может преобладать как свободноконвективный механизм, так и термокапиллярный (см. работы [1, 21, 22] и библиографию к ним). В условиях постоянной силы гравитации термокапиллярное движение, как правило, подавляется гравитационной конвекцией. Однако в случаях микротечений жидкостей важную роль играет эффект Марангони. В условиях пониженной гравитации (на борту космической станции) термокапиллярная конвекция является одной из основных форм движения жидкости. При решении прикладных задач важно знать, какой из механизмов возбуждения конвекции является преобладающим и какое влияние он оказывает на характеристики конвективных движений, возникающих при совместном действии сил гравитации и сил Марангони. В этом случае возрастает влияние таких факторов, как магнитное поле, вибрации, вращение, сложная геометрия системы и др.]2 [В то же время практические приложения демонстрируют значительное усложнение «классических» задач. Так, например, для динамики сезонной циркуляции в водоемах представляет интерес задача Рэлея-Бенара при аномальной зависимости плотности воды от температуры [24]. Развитие нанотехнологий приводит к постановке задач о поведении тонких пленок, когда подогрев снизу приводит к их разрыву и формированию чередующихся «сухих» и «несухих» областей или капель, соединенных ультратонкими пленками [25- 28]. Кроме того, сформулирован ряд задач, в которых данный тип задач о неустойчивости слоя жидкости встроен в более общую задачу тепломассопереса, например, при испарении легкой фракции из многокомпонентной смеси, в результате чего происходит изменение температуры и концентрации вблизи поверхности [5, 29]. К подобного рода задачам в усложненной постановке может быть отнесена задача об устойчивости слоя бинарного электролита между двумя ионообменными мембранами при прохождении электрического тока [30]. Усложнение практических задач, в которые в качестве отдельного блока могут входить задачи о неустойчивости жидкости и развитии конвекции, а также развитие методов моделирования и выполнения вычислительных экспериментов с использованием высокоэффективных численных методов и быстродействующих компьютеров требует разработки обобщенного подхода к моделированию и интерпретации результатов вычислительных и натурных экспериментов. Попытки разработки подобного обобщенного подхода предпринимались ранее, например, путем отнесения большого числа кооперативных явлений, демонстрирующих неустойчивости, к фазовым переходам [31]. В свою очередь, фазовые переходы в различных системах могут быть рассмотрены с единых позиций в рамках термодинамики (термодинамики нерав1 2 Цитировано по работе [20, с. 23-24] Цитировано по работе [23, c. 92] 439 Л у к а ш е в Е. А., Р а д к е в и ч Е. В., Я к о в л е в Н. Н., В а с и л ь е в а О. А. новесных фазовых переходов) [32-34]. Подход, развиваемый в работах [32-34], базируется на том представлении, что любая система может быть охарактеризована как в устойчивом, так и в неустойчивом состояниях термодинамическим потенциалом Гинзбурга- Ландау, имеющим форму свободной энергии Гиббса-Гельмгольца. В данной работе ставится задача распространить этот подход на задачу Рэлея- Бенара. В этом случае система в виде слоя жидкости, подогреваемого снизу, рассматривается первоначально как неустойчивая в термодинамическом отношении, т.е. обладающая избыточной энергией - потенциальной энергией в поле тяготения, которая может трансформироваться в гидродинамическую составляющую конвективного течения. Следует отметить, что подобного типа «термодинамическая» постановка задачи Рэлея-Бенара была сформулирована в [33]. Однако развиваемый в [32-34] подход имеет принципиальные отличия от «энтропийного» подхода, представленного в [33], поскольку в качестве главной характеристики поведения системы используется свободная энергия (энтропия входит одним из членов в уравнение Гиббса-Гельмгольца и является важной, но не единственной термодинамической величиной, определяющей поведение и динамику системы). Кроме того, развиваемый в [30-32] подход включает в себя термодинамический формализм теории неравновесных процессов Кана-Хилларда [36-44], развитый затем в работах [45, 46].]3 Задача Рэлея-Бенара представляется следующим образом. [Пусть имеется горизонтальный слой жидкости бесконечной протяженности. Снизу его подогревают, благодаря чему поддерживается температурный градиент. Выраженный в подходящих безразмерных единицах, этот градиент называется числом Рэлея Ra. Пока число Рэлея не слишком велико, жидкость остается спокойной, а теплота переносится за счет теплопроводности. Однако если Ra превосходит некоторое определенное значение, в жидкости внезапно возникает конвективное движение. <. . . > конвективные структуры весьма регулярны и могут образовывать либо цилиндрические, либо гексагональные конфигурации. Шестиугольники <. . . > представляют собой вид сверху конвективных ячеек. Жидкость поднимается в центре ячейки и опускается у ее границ или наоборот. Задача <. . . > состоит в объяснении механизма этого внезапного перехода типа «беспорядок-порядок» и в предсказании формы и устойчивости ячеек. В более точной теории следует включить в рассмотрение флуктуации. С этой проблемой тесно связана так называемая проблема Тейлора. Пусть между длинными неподвижными внешним цилиндром и концентрическим ему вращающимся внутренним цилиндром находится слой жидкости. Если скорость вращения внутреннего цилиндра, выраженная в подходящих безразмерных единицах (число Тейлора), достаточно мала, течение жидкости происходит вдоль круговых линий тока (течение Куэтта). Но если число Тейлора превосходит критическое значение, то вдоль аксиального направления появляются пространственно-периодические вихри - вихри Тейлора. <. . . > Поле скоростей, давление и температура подчиняются определенным нелинейным уравнениям гидродинамики, которые можно привести к виду с явной зависимостью от числа Рэлея Ra, задаваемого извне. При малых 3 440 Цитировано по работе [20, с. 24-25] Введение в обобщенную теорию неравновесных фазовых переходов. . . значениях Ra мы находим решение, положив компоненты скорости равными нулю. Устойчивость этого решения доказывается путем линеаризации всех уравнений относительно стационарных значений скоростей U , давления p, температуры T , где мы получаем затухающие волны. Если, однако, число Рэлея превосходит определенное критическое значение Racr , решения становятся неустойчивыми. Решения, которые становятся неустойчивыми, определяют набор мод. Реальное поле (u, T ) разлагается по этим модам с неизвестными амплитудами. Для амплитуд мод мы получаем нелинейные уравнения <. . . >. Включая в рассмотрение тепловые флуктуации, мы приходим к задаче, в которой фигурируют детерминированные силы и флуктуирующие силы.]4 Важно учесть, что их совместное действие определяет область перехода, где Ra ≈ Racr . 1. Исследование неустойчивости Рэлея-Бенара методами теории неравновесных фазовых переходов в форме Кана-Хилларда 1.1. Постановка задачи о неустойчивости Рэлея-Бенара Неустойчивость слоя жидкости, подогреваемой снизу (неустойчивость Рэлея-Бенара), обусловлена разностью плотностей верхних и нижних слоев. С термодинамической точки зрения эта неустойчивость обусловлена разницей в энергетическом состоянии частиц жидкости, образующих слой в гравитационном поле. Для оценки тенденций к «расслоению» жидкости с образованием структур макро- и мезомасштаба - конвективных ячеек, вихрей, турбулентных молей будем рассматривать только гравитационную составляющую полной энтальпии. Эта составляющая - потенциальная энергия в поле сил тяготения, которая для частицы жидкости (единичного объема) равна gx, где x - координата по высоте слоя, - плотность, g - ускорение свободного падения. Разницу плотностей верхних и нижних слоев определяет разница температур. Зависимость плотности жидкости от температуры (уравнение состояния при постоянном давлении P = const) имеет вид = ∗ 1 + βT (T - T∗ ) , (1) где ∗ , T∗ - параметры состояния отсчета; βT - коэффициент термического расширения (для оценки в некотором интервале температур он принимается постоянным, например для воды βT = 0.15 · 10-3 K-1 , для масла - 0.7 · 10-3 K-1 ). Для фиксации системы отсчета выберем верхний слой жидкости, температура которого поддерживается постоянной, с параметрами ∗ , T∗ , x∗ ; для нижнего слоя x = x0 = 0, T = T0 ; управляющим параметром будет, соответственно, T0 . Введение безразмерных плотности = / ∗ и температуры T = T /T∗ дает = 4 1 1 + β(T - 1) , (2) Цитировано по [31, c. 287] 441 Л у к а ш е в Е. А., Р а д к е в и ч Е. В., Я к о в л е в Н. Н., В а с и л ь е в а О. А. где β = βT T∗ - безразмерный коэффициент термического расширения. Избыточную энтальпию, характеризующую энергетическую неоднородность слоя жидкости, оценим как разницу относительно энтальпии слоев, близких к однородному состоянию (средний слой жидкости). Избыточную энтропию неоднородной системы представим в форме энтропии «идеальной смеси» частиц жидкости, отличающихся энергией в гравитационном поле. Эти две термодинамические функции позволяют представить свободную энергию Гиббса (или Гельмгольца в зависимости от задачи) с помощью уравнения Гиббса-Гельмгольца. Другая форма этого представления - потенциал Гинзбурга-Ландау, который, как правило, записывается в виде полинома четвертой степени. Графики потенциала Гинзбурга-Ландау для первого и второго случаев имеют один и тот же характерный вид. Однако их принципиальное отличие заключается в том, что изначально для построения зависимости свободной энергии Гиббса используются физические (термодинамические) представления, которые при формализации дают «математический результат» - явный вид зависимости. Потенциал Гинзбурга-Ландау изначально строится формальным разложением в ряд с удержанием членов разложения вплоть до четвертой степени. В таком виде у этого представления появляется возможность описывать явления неустойчивости (если потенциал имеет один минимум, то система термодинамически устойчива; если же потенциал имеет два минимума и один максимум между ними, тогда система является термодинамически неустойчивой, поэтому в ней могут реализоваться процессы типа неравновесных фазовых переходов). После формальной математической процедуры необходима физическая интерпретация коэффициентов этого разложения в ряд (их знаков и численных значений). Однако, как правило, при рассмотрении какой-либо задачи основное внимание уделяется математической стороне, а физическая интерпретация опускается. В связи с этим ниже используется представление потенциала Гинзбурга- Ландау в форме потенциала Гиббса-Гельмгольца. В этом потенциале энтальпия и энтропия могут быть представлены в форме, которая описывает не только идеальные системы, но и системы с учетом взаимодействия и какойлибо другой статистики для представления энтропии. Однако для анализа возможностей подхода на первом шаге следует использовать максимально «идеализированное» описание. Такая постановка задачи, как уже отмечалось, выделяет из полной энтальпии гравитационную составляющую, которая представляет собой потенциальную энергию частицы жидкости в гравитационном поле. Следует отметить, что такая постановка задачи хорошо известна, например, в форме задачи определения седиментационно-диффузионного распределения газов в атмосфере по высоте. Равновесное распределение газов обусловлено двумя противоборствующими тенденциями: стремлением к минимуму энергии и максимуму энтропии, реализуемыми как седиментация и диффузия соответственно. В этой задаче свободная энергия Гельмгольца в виде седиментационно-диффузионного потенциала определяет распределение газов по высоте. Приведенное отступление является обоснованием использования термодинамического подхода к задаче Рэлея-Бенара. Этот подход первоначально 442 Введение в обобщенную теорию неравновесных фазовых переходов. . . предполагает оценку потенциальной возможности реализации термодинамической неустойчивости с помощью критерия в форме потенциала Гинзбурга- Ландау. После чего этот подход предполагает оценку возможности реализации «диффузионного расслоения» исходной системы (макронеоднородной, но микрооднородной) с формированием структур мезо- и макроскопического масштаба с помощью теории неравновесных фазовых переходов с использованием аппарата, берущего начало в теории спинодального распада Кана- Хилларда. Термодинамический подход к задаче Рэлея-Бенара, используюший формализм Кана-Хилларда [36-44], не предполагает отрицания классического «механического» подхода. Здесь «термодинамическое решение» задачи в рамках теории неравновесных фазовых переходов использовано для исследования и иллюстрации возможностей подхода с целью последующего обобщения на другие примеры эволюции неустойчивых систем. Для некоторых систем из этого ряда подход механики сплошных сред, как отмечается в [47], сталкивается с определенными трудностями. 1.2. Конструирование потенциала Гинзбурга-Ландау для задачи Рэлея-Бенара в форме свободной энергии Гиббса-Гельмгольца Известно [48], что распределение температуры по длине тела (стержень, слой жидкости), границы которого поддерживаются при постоянных разных температурах, является линейным. Поэтому для распределения температуры в слое жидкости в задаче Рэлея-Бенара будем полагать T (x) = T∗ - (T0 - T∗ ) x - x∗ , x∗ (3) где T∗ - температура верхнего слоя (минимальная), T0 - температура нижнего слоя (максимальная), являющаяся управляющим параметром, x∗ - высота слоя жидкости, x - текущая координата по высоте слоя (на нижней «горячей» стороне x = x0 = 0). Средний слой жидкости является нейтральным в седиментационном отношении; верхний слой, имеющий минимальную температуру и, соответственно, большую плотность, склонен к положительной седиментации (осаждению), а нижний, имеющий большую температуру и меньшую плотность, - к отрицательной седиментации (всплыванию). По отношению к среднему слою оба этих слоя (верхний и нижний) обладают избыточной энтальпией (гравитационной составляющей): ( ( - 0- ∗ - xm ) > 0, m )g(x0 - xm ) > 0, m )g(x∗ здесь здесь - 0- ∗ m m > 0, x∗ - xm > 0; < 0, x0 - xm < 0, xm - плотность и координата среднего слоя. Для каждого локального слоя жидкости гравитационная составляющая энтальпии может быть записана в виде m, Hg = ( - m )g(x - xm ). (4) 443 Л у к а ш е в Е. А., Р а д к е в и ч Е. В., Я к о в л е в Н. Н., В а с и л ь е в а О. А. Используя зависимости (1) и (2), энтальпию можно представить как функцию температуры. Процесс формирования неоднородного состояния слоя жидкости идет при подводе теплоты, поэтому энтальпия этого (эндотермического) процесса положительна в соответствии с термодинамической системой знаков. В целом процесс формирования неоднородной жидкости и последующего ее «расслоения» на конвективные ячейки в рамках теории неравновесных фазовых переходов будет соответствовать процессу расслоения на фазы в системах с нижней критической точкой расслоения (изменения энтальпии и энтропии в ходе процесса положительны). Выполним некоторые выкладки: 1 Tm = (T0 + T∗ ) = T∗ Tm , 2 Tm = 1 T0 1 + 1 = (T0 + 1), 2 T∗ 2 где T0 = T0 /T∗ , T = T /T∗ , 1 < T < T0 ; - m = 1 ∗ 1 + β(T - 1) - 1 = 1 + β(Tm - 1) = 1 ∗ 1 + β(T - 1) - 2 2 + β(T0 - 1) . Выразим координату x через T . Из (3) следует T = 1 + (T0 - 1)(1 - x), отсюда x=1- x - xm = x∗ (x - xm ) = x∗ T -1 ; T0 - 1 1 1 T -1 x - (1 + x0 ) = x∗ 1 - 2 . 2 2 T0 - 1 Введем обезразмеривающий параметр для гравитационной составляющей энтальпии H∗ : H∗ = ( ∗ - m )g(x∗ = 1 4 - xm ) = ∗ gx∗ (1 - ∗ gx∗ (1 0 )(1 - m )(1 - x0 ) = 1 4 - xm ) = ∗ gx∗ 1- 2 2 + β(T0 - 1) . (5) В итоге для безразмерной энтальпии h(T ) согласно (4) и (5) получаем h(T ) = Hg ( - = H∗ (1 - = - xm ) = )g(1 - xm ) m m )g(x 1 1+β(T -1) - 1- 2 2+β(T0 -1) 2 2+β(T0 -1) = 444 1- 1- T -1 - 12 (1 T0 -1 1 2 (1 + x0 ) 1 1+β(T -1) - 1 + x0 ) = 2 1 2+β(T0 -1) 2 - 2+β(T0 -1) - 2(T -1) T0 -1 . Введение в обобщенную теорию неравновесных фазовых переходов. . . Как уже отмечалось, энергетическую неоднородность системы можно оценить как неоднородность «идеальной» смеси энергетически различных частиц, используя выражение для энтропии идеальной смеси. Поскольку каждый локальный слой жидкости характеризуется энтальпией, для энтропии можем записать s(T ) = -[h ln h + (1 - h) ln(1 - h)], а потенциал Гиббса-Гельмгольца в форме свободной энергии Гиббса-Гельмгольца будет иметь вид g(T ) = h(T ) - γ(T )s(T ), (6) где γ(T ) - аналог температуры для гравитационной составляющей энтальпии. Эту величину можно представить в виде γ(T ) = h(T )/α(T ), где α(T ) - энергоемкость энергетически неоднородной системы (аналог дифференциальной теплоемкости, которая снижается при увеличении температуры). Поскольку с увеличением температуры стремление системы к конвективному перемешиванию возрастает, т. е. вклад «надмолекулярной» энтропийной составляющей (которая может быть названа конфигурационной) растет, тогда параметр µ(T ) = 1/α(T ) должен расти с повышением температуры. При таком представление (6) получает вид g(T ) = h(T ) - µ(T )h(T )s(T ), (7) где µ(T ) - возрастающая функция. Здесь прослеживаются следующие аналогии и тенденции: - ∆H = cp T (cp - изобарная теплоемкость), T = ∆H/cp , т.е. γ(T ) - аналог T , а α(T ) - аналог cp ; - с ростом T дифференциальная теплоемкость ∆cp (∂cp /∂T ) уменьшается, т. е. рост T влечет уменьшение α(T ) и, соответственно, влечет рост величин 1/α(T ) и µ(T ). Для согласования размерностей и введения в последующем обоснованных численных значений характерных параметров следует провести оценку порядка некоторых величин. Для оценки Hg будем полагать ∗ ∼ 1 · 103 кг/м3 , g ∼ 10 м/с2 , x∗ ∼ 1 · 10-2 м, тогда величина ∗ gx∗ будет иметь порядок 1 Дж/м3 или 1 · 10-3 Дж/кг3 . Для согласования порядка величин h(T ), s(T ) и g(T ) будем полагать, что частица жидкости механики сплошной среды (для случая жидкости) имеет линейный размер d ∼ 0.1 ÷ 1 мкм (1 · 10-7 ÷ 1 · 10-5 м). Тогда при d ≈ 0.2 мкм для параметра RT∗ (множитель при энтропийном члене) получим RT∗ ∼ 1 Дж/м3 и условная молекулярная масса частицы жидкости механики сплошной среды составит величину порядка 1 · 109 г/моль (молекулярная масса моля Прандтля). В этом случае величины h(T ), s(T ) и g(T ) будут сопоставимы. Для качественного анализа соответствия полученных выражений поведению системы в задаче Рэлея-Бенара введем новую безразмерную переменную, представляющую изменение температуры. Поскольку интервал измене445 Л у к а ш е в Е. А., Р а д к е в и ч Е. В., Я к о в л е в Н. Н., В а с и л ь е в а О. А. ния ранее введенной безразмерной температуры T от 1 до T0 (1 < T < T0 ), введем переменную θ = (T - 1)/(T0 - 1) с интервалом изменений 0 < θ < 1. Перепишем ранее полученные выражения как функции от θ : = - m = 1 1 + β(T - 1) 1 + β(T0 - 1)θ T -1 x=1- = 1 - θ; T0 - 1 1 ∗ 1 = 1 + β(T - 1) - 2 2 + β(T0 - 1) = ; = 1 ∗ 1 + β(T0 - 1)θ - 2 2 + β(T0 - 1) ; 1 T -1 1 = x∗ (1 - 2θ); x - xm = x∗ (x - xm ) = x∗ 1 - 2 2 2 T0 - 1 h(θ) = 1 1+β(T0 -1)θ - 1- 2 2+β(T0 -1) (1 - 2θ) 2 2+β(T0 -1) ; s(θ) = -[h(θ) ln h(θ) + (1 - h(θ)) ln(1 - h(θ))]; g(θ) = h(θ) - µ(θ)h(θ)s(θ). (8) В задаче Рэлея-Бенара при малом градиенте температуры на толщине слоя жидкости система остается устойчивой по отношению к конвекции, т. е. при малых градиентах температуры с передачей тепловой энергии от источника к стоку справляется механизм теплопроводности. При больших градиентах температуры механизм теплопроводности не может обеспечить большую мощность передачи энергии и система «переключается» на другой - конвективный механизм. Дополнительно при малых градиентах температуры устойчивость системы можно представить введением стабилизирующего фактора - «молекулярной» составляющей энтропии; в этом случае устойчивость реализуется за счет механизма теплопроводности, поэтому энтропию неоднородной в термическом отношении системы как «идеальной» смеси можно представить в виде ϕ(θ) = -[θ ln θ + (1 - θ) ln(1 - θ)]. Стремление системы к термически однородному состоянию обеспечивается теплопроводностью, и чем больше эта тенденция, тем меньше тенденция к реализации конвективного перемешивания, т. е. эти две тенденции имеют противоположную направленность. Поэтому вклад «молекулярной» составляющей энтропии в свободную энергию может быть представлен так: ψ(θ) = η(θ)θϕ(θ), где параметр η(θ) является аналогом параметра µ(θ) в (7). 446 (9) Введение в обобщенную теорию неравновесных фазовых переходов. . . Таким образом, потенциал Гинзбурга-Ландау в форме Гиббса-Гельмгольца с учетом двух тенденций для реализации перемешивания (молекулярного и конвективного) может быть представлен в виде f (θ) = h(θ) - µ(θ)h(θ)s(θ) + ψ(θ). (10) Кроме этого, вместо зависимости ψ(θ) (9) для исследования качественного поведения может быть предложена зависимость ω(θ) = η(θ)ϕ(θ), которая для потенциала Гинзбурга-Ландау дает представление, аналогичное (10): z(θ) = h(θ) - µ(θ)h(θ)s(θ) + ω(θ). (11) Иллюстрации зависимостей g(θ), f (θ) и z(θ) (см. формулы (8), (10) и (11)) приведены на рис. 1. Для качественного анализа тенденций в поведении системы параметры µ = µ(θ) и η = η(θ) приняты постоянными (µ = 2.0, η = 0.2) с учетом того, что с ростом температуры эти параметры растут. Зависимость g(θ) (см.рис. 1, слева) позволяет предположить, что может реализоваться следующий сценарий развития неустойчивости подогреваемого снизу слоя жидкости. В областях, где g(θ) < 0, слои жидкости неустойчивы в отношении конвективного перемешивания. Это области, лежащие ниже и выше средних слоев. В средних слоях g(θ) > 0, и здесь возмущения должны затухать, поэтому в случае возникновения течения верхних слоев вниз, а нижних вверх, из-за нарастания в них возмущений средние слои, оставаясь устойчивыми, будут представлять собою условно «монолитные образования». Эти «монолитные образования» могут двигаться как единое целое, т. е. могут только увлекаться соседними потоками без значительных относительных перемещений внутри себя. Такого рода перемещения в силу затухания в них возмущений должны иметь резко спадающую интенсивность от периферии к центру. Если проследить за ходом зависимости g(θ) от левого максимума вниз 22 d g(θ) < 0) до левого минимума, то можно отметить, что перегиб кривой ( dθ лежит в области, где g(θ) < 0. Это означает, что на участке, где (справа от точки перегиба) слои жидкости находятся в лабильной области, т. е. Рис. 1. Зависимости g(θ), f (θ) и z(θ) при β = 0.2; T0 = 2.0 (µ = 2.0, η = 0.2) [Figure 1. Plots of g(θ), f (θ) and z(θ), when β = 0.2; T0 = 2.0 (µ = 2.0, η = 0.2); see Eqs. (8), (10), and (11)] 447 Л у к а ш е в Е. А., Р а д к е в и ч Е. В., Я к о в л е в Н. Н., В а с и л ь е в а О. А. неустойчивы по отношению к бесконечно малым возмущениям. Поэтому можно предположить, что эти слои являются «генераторами» возмущений. Слои, d2 где dθ 2 g(θ) > 0, находятся в области метастабильности, поэтому являются неустойчивыми к возмущениям, превышающим некоторый критический уровень (аналог зародыша Гиббса в теории неравновесных фазовых передов), которые будут нарастать, формируя устойчивое конвективное движение. Можно также предположить, что для этих слоев внешними по отношению к ним будут возмущения, генерируемые в лабильной области. Последние в своей области могут быстро наращивать свою интенсивность до достаточно высокого уровня, который может обеспечить развитие неустойчивости и в метастабильной области. Следует отметить одну особенность данной системы: лабильные области (нижние и верхние) непосредственно граничат не только с метастабильными, но и с устойчивыми областями (средние слои). Такая ситуация нехарактерна для неравновесных фазовых переходов, где лабильные области граничат только с метастабильными областями. Эта особенность обусловлена «геометрической» спецификой задачи, в которой неустойчивыми оказываются верхние и нижние (по отношению к центральным) слои. Введение поправки на возможность частичного термического снижения неоднородности - представление потенциала Гинзбурга-Ландау в виде f (θ) (см. формулу (10), рис. 1, в центре) или z(θ) (см. формулу (11), рис. 1, справа) на качественном поведении системы принципиально не сказывается. Система, как и следовало ожидать, становится в некоторой степени более устойчивой по отношению к конвекции. Если ситуацию, представленную на рис. 1, избрать в качестве примера сравнения, то можно отметить, что снижение градиента температуры в целом по слою жидкости может привести к тому, что система оказывается устойчивой по отношению к конвекции: g(θ) > 0 для всех слоев жидкости (рис. 2). Увеличение градиента температуры в целом по слою жидкости (рис. 3) приводит к обратному - области, где слои устойчивы по отношению к конвекции, сужаются, а области, где слои являются неустойчивыми, расширяются. Следует отметить также несимметричность зависимости (8), которая обусловлена, прежде всего, зависимостью плотности от температуры. Эта несимметричность фиксирует тот факт, что верхние слои у «холодной» поверхности Рис. 2. Зависимости g(θ), f (θ) и z(θ) при β = 0.2; T0 = 1.2 (µ = 1.2, η = 0.2) [Figure 2. Plots of g(θ), f (θ) and z(θ), when β = 0.2; T0 = 1.2 (µ = 1.2, η = 0.2); see Eqs. (8), (10), and (11)] 448 Введение в обобщенную теорию неравновесных фазовых переходов. . . Рис. 3. Зависимости g(θ), f (θ) и z(θ) при β = 0.2; T0 = 3.0 (µ = 3.0, η = 0.2) [Figure 3. Plots of g(θ), f (θ) and z(θ), when β = 0.2; T0 = 3.0 (µ = 3.0, η = 0.2); see Eqs. (8), (10), and (11)] являются более устойчивыми, чем нижние слои у «горячей». Описанный выше сценарий развития конвективной неустойчивости является предполагаемым, поскольку базируется только на термодинамическом представлении системы, т. е. на ее потенциальных (внутренних) возможностях. В то же время термодинамическое рассмотрение констатирует, что существуют потенциальные возможности к расслаиванию первоначально «конвективно» однородной системы на конвективные ячейки (ячейки Бенара). В целом же развитие конвективной неустойчивости определяется еще и «кинетическими» свойствами системы, т. е. характером развития возмущений во времени. В рамках подхода к описанию неравновесных фазовых переходов, развиваемого в [32, 51], эти свойства обусловлены интенсивностью процессов, условно называемых «отрицательной» диффузией. Ранее было показано, что такой подход может быть распространен на ламинарно-турбулентный переход [32]. Здесь данный подход проходит проверку на хорошо известной задаче о неустойчивости Рэлея-Бенара. Оказалось, что для этой задачи потенциал Гинзбурга-Ландау имеет те же характерные черты, что и для задачи о ламинарно-турбулентном переходе. Основные различия заключаются в различиях «геометрии» систем. Основные тенденции изменения функций g(θ), f (θ) и z(θ) при изменении параметров µ(θ) и η(θ) сводятся к следующим. 1. Увеличение µ(θ) при η(θ) = const приводит к переходу от состояния устойчивости к состоянию неустойчивости. Этот переход соответствует бифуркационному значению параметра µ(θ), например, для g(θ) при η(θ) = 0.2, T0 = 2.0, β = 0.2 это значение равно µ∗ср = 1.45, т. е. этот параметр является аналогом числа Рэлея в классической задаче о неустойчивости Рэлея-Бенара. 2. Увеличение η(θ) при µ(θ) = const приводит к «подавлению неустойчивости» - зона устойчивых к расслаиванию слоев расширяется - более высокое значение параметра η(θ) дает более высокое значение µ∗ср . Таким образом, параметр η(θ) отслеживает тенденцию вклада молекулярной составляющей, ответственной за поддержание термической однородности системы. В этом отношении параметр η(θ), так же как и параметр µ(θ), отслеживает те тенденции, которые в классической постановке задачи 449 Л у к а ш е в Е. А., Р а д к е в и ч Е. В., Я к о в л е в Н. Н., В а с и л ь е в а О. А. представлены числом Рэлея5 Ra = gβ∆T L2 , νλ а именно, параметр η(θ) определяется вязкостью и теплопроводностью. Эти тенденции представимы как η(θ) ∼ νλ, т. е. высокое значение λ способствует выравниванию термической неоднородности и тем самым препятствует наступлению конвективной неустойчивости. Высокое значение ν также препятствует неустойчивости, но за счет другого механизма. Если рост λ снижает градиент температуры по слою жидкости и таким путем стабилизирует его, то рост ν не влияет на сам градиент температуры, но повышает то его значение, которое необходимо для перехода к состоянию неустойчивости по отношению к конвекции. Таким образом, получаем, что классический критерий Рэлея разбивается на два критерия: µ(θ) ∼ Ra1 и η(θ) ∼ 1/Ra2 , которые отслеживают состояния неоднородной жидкости на макроуровне (мезоуровне) и микроуровне соответственно. 2. Анализ термодинамических соотношений для задачи Рэлея-Бенара Как отмечалось выше, глобальную неоднородность системы можно характеризовать как неоднородность распределения энтальпии по высоте слоя жидкости. Это уравнение с учетом гравитационной составляющей может быть представлено в виде [52] (γ - 1) (dW - T dS) - P dV - 2 g dx = P∗ ( - ∗ ) - (γ - 1) T S - = d (γ - 1) W + ∗ + (γ - 1)W + P∗ d + (γ - 1)Sd( T ) + ∗ P 2 g(x - x∗ ) + d + 2 g(x - x∗ ) d . (12) Здесь - (γ - 1)W + P∗ d + (γ - 1)S d( T ) + ∗ = (-γ - 1)W - P∗ ∗ P d + 2 g(x - x∗ ) d = + (γ - 1)T S + S dT P + + 2 g(x - x∗ ) d ; d x - координата по высоте слоя жидкости, - плотность, g - ускорение свободного падения, W - внутренняя энергия, g dx - гравитационная часть внутренней энергии. Величина (γ - 1)T dS + g dx может интерпретироваться как добавка к внутренней энергии (энтальпии), обусловленная неоднородностью системы (по плотности и температуре). Сумма четырех слагаемых 5 Здесь g - ускорение свободного падения, β - коэффициент термического расширения жидкости, ∆T - разница температур между «горячей» и «холодной» стенками, L - характерный линейный размер (высота нагреваемого слоя жидкости), ν - кинематическая вязкость жидкости, λ - температуропроводность жидкости. 450 Введение в обобщенную теорию неравновесных фазовых переходов. . . в (12) - полный дифференциал для локально равновесного состояния (для конкретного состояния системы общая энергия E = const). Дополнительная «ёмкость» энергии неоднородной системы может быть охарактеризована величиной вида g (γ - 1) dS + dx (13) T (аналог теплоемкости однородной системы). Можно полагать, что эта величина будет критерием перехода от режима молекулярной теплопроводности к режиму конвективной теплопроводности, поскольку «ёмкость» жидкости для энергии как молекулярной системы может быть превышена. Задача Рэлея-Бенара является классической задачей, на которой проверяются и оттачиваются различные приемы решения. Кроме того, различные подходы, опробованные на этой задаче, позволяют лучше отразить различные физические стороны этого явления - неустойчивость Рэлея-Бенара и другие, которые могут быть отнесены к неравновесным фазовым переходам. Так, например, в [53, с. 60] отмечается одна характерная черта: «стационарное состояние соответствует минимальному производству энтропии <. . . >. Однако стоит переступить порог неустойчивости, как ситуация изменяется. Вихри Бенара производят больше энтропии, чем это следует из теоремы о минимуме производства энтропии». Таким образом, можно говорить, что величина в (13) характеризует ситуацию, когда энтропии производится больше, чем в случае молекулярной теплопроводности. Поскольку энтропия, согласно Больцману, является «структурно чувствительной» функцией, именно величина в (13) характеризует структурную перестройку системы или неравновесный фазовый переход. В дальнейшем будем считать, что зависимость плотности жидкости от температуры (уравнение состояния при постоянном давлении P = const) имеет вид 1 ∗ . = = 1 + βT (T - T∗ ) 1 + β(T - 1) Для локальных состояний введем безразмерные относительные величины εB = g ∗ x∗ /W∗ , W = W/W∗ , S = S/S∗ ; P = P/P∗ , W∗ = S∗ T∗ , β = βT /T∗ , где εB - безразмерная локальная гравитационная составляющая энтальпии. Тогда в силу (2) уравнение (12) можно переписать в виде 1 W∗ (γ - 1) dW - (γ - 1) T dS - P dV - 2 g dx = ∗ = d (γ - 1) W + - 1 - (γ - 1) T S - εB + 1 W∗ ∗ x-1 + (1 + β(T - 1))2 P∗ dT P -(γ - 1)W - + (γ - 1)T S + (γ - 1) S + + 2 g(x - x∗ ) d . d ∗ Здесь 451 Л у к а ш е в Е. А., Р а д к е в и ч Е. В., Я к о в л е в Н. Н., В а с и л ь е в а О. А. 1 W∗ -(γ - 1)W - ∗ = P∗ + (γ - 1)T S + S ∗ dT P + + 2 g(x - x∗ ) d = d dT 1 -(γ - 1)W - 1 + (γ - 1)T S 1 + W∗ T d dT = T d T d dT =- 1 + β(T - 1) ⇒ βT + 1+ P + 2 εB (x - 1) d ; dT β-1 . = d T βT Потребуем, чтобы -(γ - 1)W - 1 + (γ - 1) β-1 P + 2 εB (x - 1) = 0 β (это требование равносильно реализации локального равновесия), отсюда S+ β-1 S + 2 2 εB (1 - x), β что соответствует уравнению состояния адиабатического процесса в поле силы тяжести. Это соотношение определяет многообразие локального равновесия в фазовом пространстве, на котором приращение энтальпии P = (γ - 1) W + - (γ - 1) 1 (γ - 1) (dW - T dS) - P dV - W∗ 2 g dx = = d (γ - 1) W - (γ - 1) T S + - 1 + εB 1-x (1 + β(T - 1))2 (14) - полный дифференциал. Здесь (W - T S) T =T∗ =1,S=S∗ =1,W =W∗ =1 1-x (1 + β(T - 1))2 T =T∗ =1,x=x∗ =1 = 0, = 0. Таким образом, система уравнений состояния = 1 , P = (γ - 1) W + - (γ - 1) β-1 S + 2 2 εB (1 - x) 1 + β(T - 1) β определяет в фазовом пространстве многообразие локального равновесия, на котором левая часть (14) есть полный дифференциал. 452 Введение в обобщенную теорию неравновесных фазовых переходов. . . 3. Гравитационно-конвективный потенциал для термодинамических систем, в которых реализуется неустойчивость Рэлея-Бенара По аналогии с [51] введем относительные величины: энтальпию, энтропию и свободную энергию Гиббса. Для этого предварительно введем аналог параметра турбулизации [51] (параметр порядка [31]), а именно параметр конвектизации W - TS ξ2 = TS для процесса с избыточной энергией W - T S > 0. Относительную величину, представляющую собой суммарный вклад гравитационной и конвективной составляющих в энтальпию с избыточной энергией, определим в виде (γ - 1) (W - T S) + - 1 + εB h(ξ) = 1-x (1+β(T -1))2 TS = (γ - 1)ξ 2 + νB ; безразмерную величину - 1 + εB νB = 1-x (1+β(T -1))2 TS назовем числом Бенара. Здесь - 1 + εB 1-x (1 + β(T - 1))2 = (1 + β(T - 1)) - (1 + β(T - 1))2 + εB (1 - x) (1 + β(T - 1))2 . Положим Y = 1 + β(T - 1) и найдем корни уравнения Y 2 - Y - εB (1 - x) = 0: Y± = 1 1± 2 1 + 4εB (1 - x) . 1. Следовательно, νB > 0, если 1 1- 2 1 + 4εB (1 - x) < 1 + β(T - 1) < 1 1+ 2 1 + 4εB (1 - x) или max 0; 1 - 1 2β 1+ 1 + 4εB (1 - x) <1+ <T < 1 2β 1 + 4εB (1 - x) - 1 . (15) 453 Л у к а ш е в Е. А., Р а д к е в и ч Е. В., Я к о в л е в Н. Н., В а с и л ь е в а О. А. 2. В то же время νB < 0, если T >1+ 1 1 + 4εB (1 - x) - 1 2β или 0<T <1- 1 2β 1+ 1 + 4εB (1 - x) , при условии, что 1 1- 2β 1+ 1 + 4εB (1 - x) > 0. Энтропию и свободную энергию Гиббса определим следующим образом: s(ξ) = - (γ - 1)ξ 2 + νB ln (γ - 1)ξ 2 + νB - - 1 - νB - (γ - 1)ξ 2 ln 1 - νB - (γ - 1)ξ 2 ; g(ξ) = h(ξ) - αh(ξ)s(ξ). Здесь - 1 - νB νB < ξ2 < . (γ - 1) (γ - 1) Из эксперимента известна полостчатая структура турбулентного потока, в котором перемежаются зоны с большой и малой турбулизацией (см. [37, 38]). Такой же структурой обладает и конвективная неустойчивость Рэлея-Бенара. Для симметричного относительно ξ = 0 потенциала g(ξ) накачка внутренней энергии, вводимая оператором Кана-Хилларда, определяется опорной касательной, т. е. отрицательными минимумами (ξ ± , g(ξ ± )) потенциала g(ξ), где ξ - = -ξ + . Это определяет одинаковые по сути амплитуды |ξ(x, t)| в соседних областях полостчатой структуры решения, что нарушает чередование зон с большой и малой конвектизацией. Более того, в первом приближении границы областей перемежаемости стоят. Эти факты требует коррекции зависимости свободной энергии. В случае, когда выполнено условие (15), имеем νB > 0, структура относительной энтальпии, энтропии и свободной энергии Гиббса такая же, как в [51] для ламинарно-турбулентного перехода. Коррекцию зависимости свободной энергии проводим по аналогии с [51]: g(ξ) =          g(ξ, νB ), νB + ξ 2 1 + α(1-h(ξ))2 (1-νB )2 0<ξ< νB ln(νB )+ +(1 - h(ξ)) ln(1 - h(ξ)) (16) , - где lim g (ξ) = lim g (ξ). ξ→0- 454 1-νB (γ-1) ; ξ→0+ 1-νB (γ-1) < ξ < 0, Введение в обобщенную теорию неравновесных фазовых переходов. . . В случае, когда   - (γ - 1)ξ 2 + β ln (γ - 1)ξ 2 + β -    -(1 - (γ - 1)ξ 2 - β) ln 1 - (γ - 1)ξ 2 - β , ξ 0; s(ξ) = -β ln (γ - 1)ξ 2 + β -   2 )6   ln 1 - (γ - 1)ξ 2 - β) , ξ < 0, - (1-β-(γ-1)ξ 5 (1-β) а h(ξ) = (γ - 1)ξ 2 , получаем g(ξ) =        g(ξ, νB ), νB + ξ2 1+ α(1-h(ξ))2 (1-νB )2 0<ξ< √ 1 - νB ; νB ln(νB )+ +(1 - h(ξ)) ln(1 - h(ξ)) √ , - 1 - νB < ξ < 0, (17) где lim g (ξ) = lim g (ξ). ξ→0- ξ→0+ Дальнейшие расчеты проведены в случае νB > 0. Случай νB < 0 требует более сложной коррекции свободной энергии Гиббса, и мы его рассмотрим позднее. Приведем зависимость графиков g(ξ) для фиксированных νB ∈ {0.1; 0.01} при α ∈ {0.1; . . . ; 0.8}. Эти зависимости проиллюстрированы на рис. 4. Они даны в безразмерных переменных. Рис. 4. Потенциал в форме свободной энергии g для νB = 0.1 (сверху) и νB = 0.01 (снизу) при α = 0.1 (линия 1); 0.2 (линия 2); 0.3 (линия 3); 0.4 (линия 4); 0.5 (линия 5); 0.6 (линия 6); 0.7 (линия 7); 0.8 (линия 8) [Figure 4. Plots of the potential in the form of the free energy g for νB = 0.1 (top) and νB = 0.01 (bottom), when α = 0.1 (line 1); 0.2 (line 2); 0.3 (line 3); 0.4 (line 4); 0.5 (line 5); 0.6 (line 6); 0.7 (line 7); 0.8 (line 8)] 455 Л у к а ш е в Е. А., Р а д к е в и ч Е. В., Я к о в л е в Н. Н., В а с и л ь е в а О. А. По точкам (νB , α) ∈ {(0.5, 0.7); (0.2, 1.1); (0.1, 1.5); (0.05, 1.7); (0.01, 1.9)} строим функцию α = α(νB ). Подставляя построенное α(νB ) в (17), получим (1) (1) потенциал gB (ξ, νB ). Графики результирующего потенциала gB приведены на рис. 5. (1) Рис. 5. Графики результирующего потенциала gB для νB = 0.01 (линия 1); 0.06 (линия 2); 0.11 (линия 3); 0.16 (линия 4); 0.21 (линии 5 и 7); 0.26 (линии 6 и 8); 0.31 (линия 9); 0.36 (линия 10); 0.41 (линия 11); 0.46 (линия 12) (1) [Figure 5. Plots of the resultant potential gB for νB = 0.01 (линия 1); 0.06 (line 2); 0.11 (line 3); 0.16 (line 4); 0.21 (lines 5 and 7); 0.26 (lines 6 and 8); 0.31 (line 9); 0.36 (line 10); 0.41 (line 11); 0.46 (line 12) ] 4. К построению модели конвективной неустойчивости Рэлея-Бенара как неравновесного фазового перехода Приведенное выше построение иллюстрирует тот факт, что термодинамический метод позволяет установить потенциальные возможности эволюции системы - в данном случае стремление к расслаиванию на устойчивые и неустойчивые слои жидкости по отношению к конвекции. Следует отметить, что использованный термодинамический анализ позволяет только установить тенденции развития процесса перехода к конвекции, т. е. потенциальную возможность такой реализации. В то же время выбор траектории и скорости эволюции системы определяются кинетикой процесса развития неустойчивости. Таким образом, термодинамический анализ позволяет заключить, что в слое жидкости, подогреваемом снизу, за счет неоднородного распределения плотности в гравитационном поле возникают области устойчивых и неустойчивых состояний. Неустойчивые состояния являются аналогами метастабильных и лабильных состояний в теории неравновесных фазовых переходов. Это не означает, что в областях устойчивости не будут наблюдаться признаки конвекции при ее развитом состоянии; просто в этих областях диффузия возмущений будет их размывать, тогда как в областях неустойчивости, особенно в лабильной области, процесс «отрицательной» (кановской) диффузии будет их концентрировать. Можно предположить, что области неустойчивости (в отношении конвекции) первоначально гомогенной системы являются источниками возмущений, а области устойчивости - стоками. Все эти предположения могут быть проверены путем численного эксперимента, т. е. путем расчетов с использованием математической модели, приведенной ниже, включающей, кроме уравнений гидродинамики (первоначальной в самой простой форме), математическую модель неравновесного фазового перехода, сформулированную в терминах теории спинодального распада 456 Введение в обобщенную теорию неравновесных фазовых переходов. . . Кана-Хилларда. Эта математическая модель конвективной неустойчивости Рэлея-Бенара как аналога неравновесного фазового перехода может усложняться начиная с исходной термодинамической модели. Так же как в [51] при моделировании конвективной неустойчивости Рэлея-Бенара, сформулируем три феноменологические гипотезы: этот процесс является (I) процессом эволюции системы с избыточной энергией и представляет собой (II) неравновесный фазовый переход, механизм которого - (III) диффузионное расслаивание. Процесс диффузионного расслаивания на потенциально устойчивые и неустойчивые в конвективном отношении слои жидкости далее развивается с формированием конвективного течения и может быть представлен в рамках формализма механики сплошной среды (в простейшем иллюстративном варианте - это система уравнений Эйлера). В этом случае стартуем с системы уравнений d ∗ T - divU = 0, dt β d U + ∇P = ge1 , dt (18) d W + P divU = 0, dt γ d S + νT γ-1 e-κS divU = 0 dt при справедливости уравнений состояния = ∗ 1 + β T - T∗ , P = (γ - 1) W + P∗ ∗ - βT - T∗ T∗ S + 2 2 g x∗ - x , βT где U - поле скоростей, e1 = (1, 0, 0), W - внутренняя энергия, S - макроэнтропия, d F = ∂t F + (U · ∂x )F. dt Из первого и четвертого уравнений (18) имеем 2γ-1 d ∗ (T S) + - S + νT γ-1 e-κS divU = 0. dt β Положим Q = T S. Тогда 2γ-1 κ d ∗ Q+ - Q + νT γ-1 e- T Q divU = 0. dt Tβ Отсюда d ∗ T - divU = 0, dt β d U + ∇P = ge1 , dt 2γ-1 κ d ∗ Q+ - Q + νT γ-1 e- T Q divU = 0, dt Tβ 2γ-1 d ∗ (W - T S) + P + S - νT γ-1 e-κS divU = 0. dt β 457 Л у к а ш е в Е. А., Р а д к е в и ч Е. В., Я к о в л е в Н. Н., В а с и л ь е в а О. А. В безразмерных переменных x∗ = U∗ t∗ , W∗ = T∗ S∗ , P∗ = 2 ∗ U∗ , S∗ T∗ = U∗2 , γ γ-1 β = βT /T∗ , κ = κS∗ , ν = νeκ T∗ / ∗ S∗ , t = t/t∗ , x = x/x∗ , U = U /U∗ , P = P/( ∗ U∗2 ), = / ∗ , W = W/W∗ имеем d 1 T- 1 + β(T - 1) divx U = 0, dt β d U + 1 + β(T - 1) ∇x P = ge1 , dt d Q + 1 + β(T - 1) × dt 2γ-1 1 × - Q + ν(T ) γ-1 exp -κ Tβ d (W - Q) + 1 + β(T - 1) × dt 2γ-1 1 Q - ν(T ) γ-1 exp κ × P+ βT Q T (19) -1 Q T divx U = 0, -1 divx U = 0. Здесь divx , ∇x - операторы в безразмерных координатах x. Система уравнений состояния примет вид = 1 1 + β(T - 1) , P = (γ - 1) W + - (γ - 1) β-1 S + 2 2 εB (1 - x). β Теперь введем параметр конвекцизации ξ2 = W -Q Q при условии избыточности энергий W -Q > 0. Тогда в переменных {ξ, U , T , Q} последнее уравнение (19) примет вид d (W - Q) + 1 + β(T - 1) × dt 2γ-1 1 Q × P- Q + ν(T ) γ-1 exp -κ -1 T βT = 2ξ Q d 1 ξ+ -ξ 2 1 + β(T - 1) dt ξQ + 1 + β(T - 1) 458 P+ 1 βT - 1 Tβ 2γ-1 divx U = 2γ-1 Q + ν(T ) γ-1 exp -κ Q - ν(T ) γ-1 exp -κ Q T -1 Q T -1 divx U = 0. + Введение в обобщенную теорию неравновесных фазовых переходов. . . Теперь введем вязкость и накачку внутренней энергии. Известно, что вязкость жидкостей уменьшается с ростом температуры и может быть представлена формулой Френкеля-Андраде Evisc , (20) RT где R - универсальная газовая постоянная, Evisc - энергия активации вязкости. По усредненным значениям подобранных констант C и Evisc рассчитываются значения вязкости при различных температурах согласно (20). В безразмерных переменных формула Френкеля-Андраде будет иметь следующий вид: η = C exp Evisc 1 η = ε exp -1 η∗ R T , где Evisc /(RT ) = Evisc /(RT ), T /T∗ = T , ε = η∗ L∗ /(U∗ ∗ ), ε = 1/Re∗ , Re∗ - число Рейнольдса возмущаемой однородной среды. В дальнейшем мы рассмотрим только динамику жидкости: 1 d T- 1 + β(T - 1) divx U = 0, dt β d U + 1 + β(T - 1) ∇x P = dt = ε 1 + β(T - 1) exp d Q + 1 + β(T - 1) dt - 1 Evisc 1 -1 R T 2γ-1 Tβ Q + ν(T ) γ-1 exp -κ Q T -1 = ε 1 + β(T - 1) exp d 1 ξ+ -ξ 2 1 + β(T - 1) dt ξQ + 1 + β(T - 1) = 1 P+ n ∂ xj 2 j=1 1 βT D T - 1 Tβ 2γ-1 2γ-1 ∂xj Q (1 + β(T - 1))2 Q T divx U = Evisc 1 -1 R T Q + ν(T ) γ-1 exp -κ Q - ν(T ) γ-1 exp -κ ∆ U + εB e 1 , -1 Q T -1 ∆Q, + divx U = gξ (ξ) - ε2 (divx A ∇x ξ) , D = D/(U∗ L∗ ). В поле силы тяжести поверхностное натяжение неоднородно   a(T ) 0 0 1 0 , a > 1. A= 0 0 0 1 459 Л у к а ш е в Е. А., Р а д к е в и ч Е. В., Я к о в л е в Н. Н., В а с и л ь е в а О. А. Зависимость параметра a(T ) от T будет выбранна при расчете двумерной модели. Возможный вариант такой зависимости приведен в [55]. 5. Численный эксперимент для одномерной модели и анализ результатов Одномерная модель в отличие от двумерной и трехмерной не может продемонстрировать развитие неустойчивости Рэлея-Бенара с зарождением структур вихревого типа. Однако именно на одномерной модели проще всего в вычислительном эксперименте проследить отклик системы на возмущения, которые в конечном счете и являются причиной возникновения и последующего развития вихревых структур конвективного течения. В этой статье мы ограничимся численным экспериментом одномерной модели для случая положительной постоянной Бенара (νB > 0, т. е. в интервале температур условия (15)).6 Итак, в одномерном случае в безразмерных переменных модель (далее - система) примет вид 1 d T- 1 + β(T - 1) ∂x U = 0, dt β d Evisc 1 -1 U + 1 + β(T - 1) ∂x P = ε 1 + β(T - 1) exp R T dt d Q + 1 + β(T - 1) dt - 1 Q 2γ-1 Q + ν(T ) γ-1 exp -κ Tβ -1 T + 1 + β(T - 1) P+ - 1 βT = 1 Tβ ∂x U = Evisc 1 -1 R T = ε 1 + β(T - 1) exp d 1 ξ+ -ξ 2 1 + β(T - 1) dt ξQ Q 2γ-1 Q + ν(T ) γ-1 exp -κ 2γ-1 Q - ν(T ) γ-1 exp -κ 1 ∂ 2 x D d = ∂t + U ∂ x , dt T ∂x = Q T T -1 -1 Q 1 + β(T - 1) 2 gξ (ξ) 1 1 + β(T - 1) для случая, когда постоянная Бенара νB = 6 460 1 Q - 1 + εB (1 - x) 1 + β(T - 1) 2 ∂x2 U + εB , > 0, Случай νB < 0 будет рассмотрен в следующих публикациях авторов. ∂x2 Q, + ∂x U = - ε2 ∂x2 ξ , Введение в обобщенную теорию неравновесных фазовых переходов. . . где k - постоянная Больцмана, γ - адиабатическая постоянная, εB - безразмерное ускорение силы тяжести. Для x ∈ [0, 1] рассмотрим самовозбуждение однородного (в отношении конвекции) состояния. Управление температурой на нижней поверхности x=0 дается формулами T = T00 + VT t, x=0 Q00 1, T00 Q 1, x=0 = Q00 VT 1; (21) 1; начальное значение ξ0 = ξ00 = 0.5 выбирается из зоны лабильности потенциала для 1 (1 - x0 ) 0 νB = 0 - 1 + εB 2 > 0. 1 + β(T - 1) 0 Q0 На рис. 6-8 приводятся графики параметра конвектизации ξ, температуры T и скорости U для самовозбуждения однородного (в отношении конвекции) состояния, полученные из решения рассматриваемой системы в определенные временные срезы. На рис. 6 видим развитие возмущения, первоначально заданного ξ t=0 ≡ ξ0 на отрезке x ∈ [0, 1]. Возрастание температуры на нижней поверхности x = 0 (в силу условия (21)) порождает конвективную неустойчивость на верхней поверхности x = 1, что указывает на нелокальность возмущения. Последнее говорит о невозможности применения к этой задаче классической теории возмущений. Ранее отмечалось, что термодинамический анализ позволяет заключить, что верхние слои жидкости являются в большей степени неустойчивыми по отношению к конвекции, чем нижние слои, с оговоркой, что «кинетические» факторы - теплопроводность и вязкость - могут существенным образом изменить ситуацию. Таким образом, данная система может дать большое число различных сценариев благодаря различиям термодинамического и кинетических вкладов на разных стадиях развития процесса конвекции. Выполненные предварительные расчеты это подтверждают. На нижней поверхности на начальной стадии конвективная неустойчивость не наблюдается, поскольку здесь тенденция к всплыванию, реализация которой необходима для развития конвективной неустойчивости, тормозится теплопроводностью. В последующие моменты времени на границе прогретой области возникает осцилляция, т. е. зарождается конвективная неустойчивость за счет того, что механизм теплопроводности не обеспечивает выравнивания температуры и, как следствие, выравнивания плотности. Таким образом, здесь конвективное течение является следствием возникновения силы Архимеда, направленной вверх, т. е. «отрицательной» силы гравитации. На рис. 7, 8 приведено развитие возмущений температуры T и скорости U , которые первоначально заданы на отрезке x ∈ [0, 1]: T t=0 ≡ T0 , U t=0 ≡ 0. Развитие возмущений на верхней поверхности x = 1 тормозится вязкостью. На рис. 9 приведено развитие возмущений производства энтропии St = ∂t (Q/T ) при первоначально заданных T t=0 ≡ T0 , U t=0 ≡ 0 на отрезке x ∈ [0, 1]. 461 Л у к а ш е в Е. А., Р а д к е в и ч Е. В., Я к о в л е в Н. Н., В а с и л ь е в а О. А. Рис. 6. Графики параметра конвектизации ξ в определенные временные срезы при ξ00 = 0.5, U00 = 1, Q00 = 4, T00 = 3, VT = 2 [Figure 6. Plots of the parameter of convection ξ at certain time intervals, when ξ00 = 0.5, U00 = 1, Q00 = 4, T00 = 3, VT = 2] Рис. 7. Графики температуры T в определенные временные срезы при ξ00 = 0.5, U00 = 1, Q00 = 4, T00 = 3, VT = 20 [Figure 7. Plots of the temperature T at certain time intervals, when ξ00 = 0.5, U00 = 1, Q00 = 4, T00 = 3, VT = 20] Рис. 8. Графики скорости U в определенные временные срезы при ξ00 = 0.5, U00 = 1, Q00 = 4, T00 = 3, VT = 20 [Figure 8. Plots of the velocity U at certain time intervals, when ξ00 = 0.5, U00 = 1, Q00 = 4, T00 = 3, VT = 20] 462 Введение в обобщенную теорию неравновесных фазовых переходов. . . Рис. 9. Графики производства энтропии St в определенные временные срезы при ξ00 = 0, U00 = 0, Q00 = 4, T00 = 3, VT = 20 [Figure 9. Plots of entropy production St at certain time intervals, when ξ00 = 0, U00 = 0, Q00 = 4, T00 = 3, VT = 20] Для пояснения потенциальной возможности перехода рассматриваемой системы в целом от термически неоднородного, но конвективно устойчивого состояния к конвективно неустойчивому вернемся к правому рис. 5 и проведем на нем вертикаль, отвечающую начальному значению параметра конвектизации ξ = ξ0 . Как видим, с уменьшением модуля константы Бенара νB > 0 при увеличении температуры T в краевом условии вертикаль из устойчивой зоны переходит в зону лабильности, т. е пересекает участок графика свободной энергии Гиббса между точками перегиба. Вход в зону лабильности приводит к срыву локального равновесия, вызванного появлением осцилляций конвективной неустойчивости (нелокального возмущения). Отметим, что возникновение и развитие конвективной неустойчивости Рэлея-Бенара существенно отличается от ламинарно-турбулентного перехода и перехода к развитой турбулентности [51]. Для ламинарно-турбулентного перехода построенная в [51] модель дает следующий сценарий. В случае отрицательного наклона опорной касательной к графику g и при начальных значениях ξ00 в ядре лабильной области характерен режим быстрого перехода к турбулентности, который при своем развитии имеет общие черты с процессом в модели «хищник-жертва», т. е. реализуется «колеблемость» типа перемежаемости [54]. При t → ∞ поток стабилизируется и переходит к течению, в котором ξ(x, t) → ξ∞ x - U∞ t , ε где ξ∞ (x) - строго монотонно убывающая функция (типа кинка [57]) с гра+ > 0, ξ - < 0 (аналог решения уравничными значениями ξ x=0 = ξ∞ = ξ∞ x=1 нения Колмогорова-Петровского-Пискунова). В случае положительного на463 Л у к а ш е в Е. А., Р а д к е в и ч Е. В., Я к о в л е в Н. Н., В а с и л ь е в а О. А. клона опорной касательной к графику g самовозбуждение турбулентности развивается на порядок быстрее, чем в предыдущем случае отрицательного наклона опорной касательной, и при этом процесс стремится к однородному распределению параметра турбулизации ξ(x, t) → ξ∞ . Таким путем поток приходит к состоянию развитой турбулентности (ξ ≡ const), т. е. к состоянию с избыточной энергией ξ∞ > ξ00 . Этот сценарий может быть проиллюстрирован описанием экспериментальных результатов, которое приведено Л. Г. Лойцянским в [54]. «Тщательное исследование потока в трубе при рейнольдсовых числах, близких к критическим, показало, что в одном и том же фиксированном сечении трубы и при том же значении рейнольдсова числа Re = U ср d/ν может происходить чередование ламинарных и турбулентных режимов. Это явление получило наименование перемежаемости <. . . >. Причина перемежаемости режимов течения заключается в том, что турбулентность <. . . > образуется вначале в дискретных областях потока в виде «облачков» или «пятен» <. . . >, в случае трубы заполняющих поперечное сечение трубы «пробками», которые могут достигать протяженности вдоль трубы порядка нескольких десятков диаметров трубы, причем эта протяженность зависит от рейнольдсова числа потока. Основной количественной характеристикой явления перемежаемости служит доля времени существования турбулентного режима в данном сечении трубы. Эту безразмерную величину, равную нулю, если течение все время ламинарное, и единице, если течение сохраняет турбулентную форму, называют «коэффициентом перемежаемости» <. . . > γ. Величина эта зависит как от рейнольдсова числа потока, так и от расстояния x от входа в трубу. <. . . > Коэффициент перемежаемости резко возрастает в области критического значения числа Рейнольдса, причем в ближних к входу сечениях позже, чем в дальних» [54, с. 525, 526]. «<. . . > передние границы пробок при закритических режимах движутся быстрее задних, вследствие чего «пробки» растягиваются, заполняя при своем движении все большие и большие объемы трубы. Вместе с тем передний край одной «пробки» догоняет задний край смежной «пробки». Все это приводит к тому, что при закритических значения Re в удалении от входа в трубу устанавливается сплошное турбулентное движение. <. . . > Скорость передней стенки «пробки» сначала меньше, но с ростом рейнольдсова числа становится больше, чем скорость турбулентного потока на оси трубы, а скорость задней стенки значительно меньше этой скорости. Можно еще заметить, что до критического значения числа Рейнольдса скорость передней грани «пробки», наоборот, меньше скорости задней грани; это приводит к сокращению длин образующихся «пробок» и их исчезновению в ламинарном потоке» [54, с. 527]. Это качественное описание экспериментальных результатов [54] полностью соответствует результатам, полученным в [51]. В первом случае фронт турбулизации движется по потоку с некоторой скоростью, составляющей собственно скорость нарастания фракции «турбулентной» жидкости. То же происходит и на задней стенке, поэтому в лабораторной системе координат задняя стенка движется против потока, но в то же время сносится потоком. Здесь фазовый переход осуществляется от ламинарного к турбулентному течению. 464 Введение в обобщенную теорию неравновесных фазовых переходов. . . При докритических режимах ситуация обратная: происходит «растворение» зародышей «турбулентной» жидкости в нетурбулентной; здесь имеет место направление фазового перехода от турбулентного режима к ламинарному - реламинаризация. Для конвективной неустойчивости Рэлея-Бенара ситуация иная, поскольку могут существовать две области лабильности, т. е. могут реализоваться два встречных процесса развития конвекции, стартующих с нижней (x = 0) и с верхней (x = 1) поверхностей, при невозмущенной центральной области. На нижней поверхности, где гравитация «отрицательна» и мала, на начальной стадии развития возмущения существенное влияние оказывает теплопроводность. Конвективные возмущения сначала возникают на верхней границе, когда прогрев всего слоя жидкости при достаточной его высоте становится существенным, а теплопроводность не может устранить градиент температуры и, соответственно, плотности по высоте слоя. Такое поведение может быть установлено в рамках одномерной модели. В многомерной модели (двумерной, трехмерной) развитие возмущений по вертикали от верхней и нижней поверхностей при их горизонтальном смещении друг относительно друга должно приводить к возникновению вращательных движений, т. е., в конечном счете, к формированию ячеек Бенара. Возникновение этих движений в верхнем слое жидкости может привести к формированию термической неоднородности верхней поверхности, что, в свою очередь, должно привести к чередующейся неоднородности поверхностного натяжения в силу ее зависимости от температуры. Таким образом, может реализоваться неустойчивость Бенара-Марангони, усиливающая эффект первоначальной неустойчивости Рэлея-Бенара. Но исследование этого эффекта не входит в цели и задачи настоящей статьи. Конкурирующие интересы. Мы не имеем конкурирующих интересов. Авторский вклад и ответственность. Все авторы принимали участие в разработке концепции статьи и в написании рукописи. Авторы несут полную ответственность за предоставление окончательной рукописи в печать. Окончательная версия рукописи была одобрена всеми авторами. Финансирование. Исследование выполнялось без финансирования.

About the authors

Evgeniy A Lukashev

Joint-stock company Turaevo Machine-Building Design Bureau “SOYUZ”

Email: elukashov@yandex.ru
Dr. Techn. Sci., Professor; Chief Specialist 10, Turaevo st., Lytkarino, 140080, Russian Federation

Evgeniy V Radkevich

M. V. Lomonosov Moscow State University

Email: evrad07@gmail.com
Dr. Phys. & Math. Sci., Professor; Professor; Dept. of Differential Equiations Faculty of Mechanics and Mathematics, Vorob’evy gory, Moscow, 119899, Russian Federation

Nikolay N Yakovlev

Joint-stock company Turaevo Machine-Building Design Bureau “SOYUZ”

Email: amntksoyuz@mail.ru
Cand. Phys. & Math. Sci.; General Manager 10, Turaevo st., Lytkarino, 140080, Russian Federation

Ol’ga A Vasil’eva

Moscow State University of Civil Engineering

Email: asiljeva.ovas@yandex.ru
Cand. Phys. & Math. Sci., Associate Professor; Associate Professor; Dept. of Applied Mathematics 26, Yaroslavskoe shosse st., Moscow, 129337, Russian Federation

References

Supplementary files