CC BY
78
ПРИКЛАДНАЯ МАТЕМАТИКА И МЕТОДЫ МАТЕМАТИЧЕСКОГО МОДЕЛИРОВАНИЯ
УДК 532.529
Д. Ч. Ким
ДИССИПАТИВНАЯ МОДЕЛЬ НЕЛИНЕЙНЫХ ВОЛН В ЖИДКОСТИ С ПУЗЫРЬКАМИ ГАЗА
Разработана двухуровневая диссипативная математическая модель жидкости с распределенными пузырьками газа. Проведено численно-аналитическое исследование нелинейных волн в построенной модели. Существование диссипативной уединенной волны возможно лишь при наличии двух пар конкурирующих факторов, одной из которых является баланс между слабой нелинейностью и дисперсией линейных волн, а другой — баланс между притоком и оттоком теплоты в газовую фазу. В длинноволновом приближении получены уравнения Кортевега-де Фриза. Разработаны компактные разностные схемы четвертого порядка точности для нелинейного уравнения теплопроводности и обобщенного неоднородного волнового уравнения для долговременного моделирования волновых процессов.
Постоянный интерес исследователей к волновым задачам жидкости с пузырьками газа обусловлен широким распространением таких сред в природе и их разнообразным применением в современной промышленности.
Жидкость с газовыми пузырьками является сложной нелинейной дисперсной диссипативной распределенной динамической системой. Многим нелинейным эффектам из различных областей механики и физики спустя некоторое время обнаруживались аналоги и в этой системе [1]. Например, проводя аналогию между волнами сжатия в рассматриваемой среде и поверхностными волнами на мелкой воде, авторы работы [2] получили уравнение Кортевега-де Фриза (КдФ). В работе [3] рассмотрено акустическое фазовое эхо, которое аналогично циклотронному резонансу в плазме, находящейся в магнитном поле. При этом специфические свойства, присущие пузырьковой среде, оставались не выявленными. К ним можно отнести виброплотность и вибровязкость, исследованные в работе [4], и новый тип "солитона" с осциллирующей структурой, обнаруженный авторами работы [5]. Недостаточно изученным является также теплообмен между пузырьками и жидкостью, особенно его влияние на структуру и динамику нелинейных волн.
В работах [6-9] показано, что затухание колебаний стенок пузырька определяется теплообменом между фазами, а не вязкостью несущей жидкости. В работе [6] использована конечно-разностная схема второго порядка точности по пространственным переменным, а в [7] реализован ряд спектральных методов для решения таких задач. В работе [8] решена система уравнений, описывающая нестационарное одномерное движение газожидкостной среды (включая законы сохранения массы, импульса, энергии смеси, осцилляций теплопроводных пузырьков с учетом относительного движения фаз) конечно-разностными схемами первого и второго порядка точности. Авторы работы [9] для описания структуры ударных волн в двухфазной смеси использовали решение уравнения нестационарной теплопроводности внутри пузырька методом прямых, заключающимся в том, что уравнение в частных производных сводится к системе обыкновенных дифференциальных уравнений простой заменой производных по координате конечно-разностными аналогами первого и второго порядка точности. В оригинальной работе [10] авторы отказались от традиционной модели тепловой диссипации в пользу использования гиперболического уравнения теплопроводности. Это, по-видимому, оправдано для бы-стропротекающих процессов.
Следует сказать, что краевая задача динамики двухфазной газожидкостной смеси, корректно учитывающая теплообмен между фазами даже с упрощающими предположениями и допущениями, является сложной и требует новых подходов для численного решения и, в частности, повышения точности численных схем. Для настоящей работы нами разработаны вычислительные алгоритмы, основанные на компактных разностных схемах четвертого порядка точности по пространственной координате как для нелинейного уравнения теплопроводности внутри газового пузырька [11], так и для обобщенного неоднородного волнового уравнения типа Лайтхилла [12].
Диссипация кинетической энергии радиальных колебаний пузырька связана с необратимостью процессов в газе и жидкости. Полагают, что при сжатии, когда температура газа становится выше температуры жидкости, газ отдает в жидкость больше теплоты, чем получает от жидкости в процессе своего расширения, если его температура оказывается ниже температуры жидкости [12]. Диссипация тепловой энергии упрощенно учитывается эффективной вязкостью в уравнении Рэлея-Лэмба [14-15] или в уравнении Кортевега-де Фриза-Бюргерса (КдФБ) [16-17]. Феноменологический учет диссипации в модельном уравнении КдФБ соответствует затуханию, вызванному сдвиговыми напряжениями, которое увеличивается пропорционально
квадрату волнового числа и не является адекватным физическим механизмом нестационарного теплообмена пузырька с окружающей жидкостью. Хорошо известно, что нелинейные эффекты, проявляющиеся при распространении акустических волн, во многом зависят от частотного закона поглощения в среде [18]. Поэтому практика подмены истинного межфазного теплообмена квадратичным законом затухания представляется неоправданной.
Другим общепринятым приближенным подходом является учет диссипации соответствующими эффективными безразмерными числами Нуссельта [13]. В этом случае интенсивность нестационарного теплообмена пузырька с жидкостью описывается фиксированным числом Нуссельта, определяемым из приближения тонкого теплового пограничного слоя. При этом универсальных соотношений в виде алгебраических формул для чисел Нуссельта также нет, и в зависимости от свойств компонентов среды, а также степени нестационарности рассматриваемого процесса формулы для чисел Нуссельта получают из дополнительных условий или из опыта. Следует признать, что этот подход существенно упрощает вычислительную задачу и в то же время позволяет учесть в грубом приближении межфазный теплообмен, но затрудняет получение надежных количественных результатов, не дает полной картины неоднородных нестационарных тепловых полей внутри пузырька, а также другие тонкости изучаемого явления.
Таким образом, используемые в работах [13-17] допущения о характере диссипации колебаний пузырьков нельзя считать адекватными в волновых задачах. В настоящей работе приведены математическая модель жидкости пузырьковой структуры с теплопроводными пузырьками газа для описания нелинейных волн, результаты компьютерного моделирования и аналитического исследования образования стационарных уединенных волн при распространении акустических возмущений.
Математическая модель жидкости с распределенными пузырьками. В работе [5] предложена двухуровневая волновая динамическая модель жидкости, содержащей равномерно распределенные пузырьки газа. Жидкая фаза описывается на макроуровне обобщенным неоднородным волновым уравнением Лайтхилла, в котором газовая фаза представлена в виде точечных монопольных источников стока-истока жидкости. Газовые пузырьки описываются на микроуровне уравнением Рэлея-Лэмба или его модификациями. Принято допущение, что давление в пузырьке однородно по радиусу и зависит только от времени, так как рассматриваются случаи, когда скорость стенок пузырька значительно меньше скорости звука в газе. Предполагается, что пузырьки могут только пульсировать. Их осцилляциями (т.е. смещением
центра пузырьков относительно жидкости) пренебрегают — они не дают вклада в сжимаемость среды, а дипольное излучение малоэффективно по сравнению с монопольным. Газ считается идеальным с постоянной теплоемкостью, нейтральным по отношению к окружающей его жидкости. Локальные параметры состояния для жидкости или газа определяются на основе уравнений состояний для каждой фазы. Действующее давление и температура таковы, что не достигают области фазового перехода газ-жидкость. В этих условиях массообмен между дисперсионной и дисперсной фазами отсутствует.
В соответствии с изложенным, в пузырьковой жидкости принципиально важен корректный учет диссипации тепловой энергии. Для этого в модели [5] политропическое поведение газа пузырьков заменим уравнением нестационарной теплопроводности. Окончательно система уравнений, описывающая одномерное движение среды, имеет вид
д2р 2 д2р _ д f д dt2 - СдХ = dt Vöt n
1 +
pR3
1 - ^0
(1)
1 -R rR+2 r2
i - i» =
i + RV-
c p
Pb (t) - Ps ( t + R ) - Po
+
R dPB (t)
cp
dt
(2)
(ÖT + u^) - dpg = V(kVT), 0 <r<R, t> 0; (3) Y - 1 T \дЬ дг) dt v h ' ' w
^ = "i(Y - 1)кдТ dt R\(Y 1)Кдг
- YPgR
r=R
P - Po = c2(p - po); 4
Ф = ö nnR
3
(4)
(5)
(6)
где
1 Л дТ 1 dPg
u =- (y — ---r——
YP \ дг 3 dt
Pb = Pg (t) - R -
2a 4pR
R
Здесь и далее: t — время; х — координата вдоль распространения волны; с — скорость звука в жидкости без пузырьков; р — плотность чистой жидкости; р = р/р0 — безразмерная плотность жидкости; р — объемное газосодержание; R — радиус пузырька, Л = Я^)/Яо — безразмерный радиус пузырька; Р — давление жидкости; Рд — давление газа; Рв — давление на внешней поверхности пузырька; Р^ — амплитуда звукового возмущения; Т — температура газа; и — радиальная
компонента скорости газа; г — радиальная координата; к — коэффициент теплопроводности; 7 — показатель политропы; ц — коэффициент динамической вязкости; а — коэффициент поверхностного натяжения. Точки над переменными означают производную по времени, черта сверху — безразмерную величину, нижний индекс "0" у параметров указывает на невозмущенное значение.
Распространение плоских волн в двухфазной жидкости пузырьковой структуры описывается обобщенным неоднородным волновым уравнением Лайтхилла (1). В этом уравнении сделана поправка к величине газосодержания из-за сжимаемости несущей фазы. В качестве дополнительного уравнения состояния, связывающего давление жидкости, газа и радиус пузырька, используется уравнение Рэлея-Лэмба (2) (формулировка Келлера) для описания динамических и дисперсионных характеристик пузырька. Уравнение теплопроводности (3) описывает закон нестационарного теплообмена между фазами. Жидкость можно считать термостатом, задающим постоянную температуру Т^ стенки пузырьков, а пузырьки — нелинейными точечными осцилляторами, взаимодействующими с термостатом через поверхность раздела фаз. Это позволяет ограничиться решением только внутренней задачи теплопроводности для газовой фазы. Обыкновенное дифференциальное уравнение (4) для давления газа получается при принятых в модели допущениях и граничных условиях. Формула (5) является акустическим уравнением состояния жидкости. Объемное газосодержание монодисперсных пузырьков определяется формулой (6).
Сделаем следующие замечания об отличительных особенностях приведенной модели. Полагаем справедливым основное допущение механики сплошной среды — расстояния, на которых параметры течения смеси меняются существенно (вне поверхности раздела фаз), много больше размера пузырьков и расстояния между ними. Это позволяет описать газожидкостную смесь как совокупность двух квазиконтинуумов, заполняющих один и тот же объем. Каждый компонент становится сплошной средой во всем объеме, занимаемым смесью, и для нее можно применить законы сохранения парциальной массы, импульса и энергии [19]. Уравнения сохранения виртуальной сплошной среды записаны в терминах плотности, давления, скорости и т.д. для каждого компонента смеси, что впервые было сделано в работе [20], и принято односкоростное приближение (равновесие по скоростям, по терминологии работы [19]). Уравнения гидродинамики линеаризованы, а нелинейные и дисперсионные эффекты вносятся только дискретной фазой. В окончательной модели не используются параметры фиктивной гомогенной среды (плотность смеси, давление смеси
и т.п.). Это позволяет избежать ряда несоответствий, присущих квазигомогенным моделям. И наконец, межфазный теплообмен учитывается прямым решением уравнения теплопроводности с помощью специального высокоточного численного метода.
Граничные и начальные условия. Монохроматическая плоская волна распространяется в положительном направлении оси х. Точное решение нелинейного волнового уравнения (1) не может быть получено, поэтому оно интегрируется численно как эволюционная краевая задача. Начальные условия для него задаются в виде
Р (х, 0) = Р„, ™ = 0. (7)
Двухточечные граничные условия задаются слева (х = 0) в виде синусоидального источника, а справа (х = Ь) в виде безотражательного условия Зоммерфельда:
ч л л • / ч ,'dP(x,t) дР(x,t) Р(0,t) = ро + ps sm(Wt), ( + c-pdx-tl
= 0, (8)
х=Ь
где ш — частота внешнего источника. Второе условие необходимо ставить на границах ограниченной расчетной области для численного решения задачи, относящейся к неограниченным областям.
При пульсациях пузырьков возникает задача сопряженного нестационарного теплообмена, когда на движущейся границе жидкость-газ должны быть выполнены соответствующие граничные условия. Для уравнения (3) граничные условия диктуются конечностью температуры в центре и необходимостью ее равенства температуре жидкости на стенке (условие Дирихле):
Нш г2 к(Т)= 0, Т (Я,*) = Ть. (9)
Второе условие на межфазной границе является хорошим приближением. Оно следует из большого различия в значениях теплофизиче-ских параметров жидкости и газа. Если приравнять тепловые потоки на границе фаз, то можно оценить
■ (крСр )с 11/2
Ts - TL Tg — Ts
= 4-10
-3
(крСр )ь_
где Т^ — температура пристенного слоя жидкости при учете прогрева; Тс — температура газа в пузырьке. Оценка показывает, что температура границы Т? отличается от температуры жидкости Т^ на пренебрежимо малую величину.
Начальное условие для поля температуры газовой фазы имеет вид
Т(г, 0) = Ть, 0 < г < Я. (10)
Для уравнений (2), (4) и (6) ставятся естественные начальные условия:
-R(x, 0)
R(x, 0) = R0,
dt
Pg (x, 0) = Po + Rf, R0
<(x, 0) = <po.
= 0,
(11) (12) (13)
Отметим, что точность аппроксимации граничных условий не ниже точности аппрокцимации численных схем.
Алгоритмы численного интегрирования задачи. Рассмотрим основные разностные схемы численного интегрирования нелинейной системы (1)-(6). Для этого в уравнении теплопроводности пузырька введем тепловой потенциал Ф и безразмерную радиальную координатУ У-
т
Ф = / к(в)М, y =
Tl
R(t)'
(14)
Центру пузырька соответствует у = 0, а подвижной стенке у = 1. Для численных расчетов удобно ввести новую функцию связанную с тепловым потенциалом формулой $ = уФ. Для нелинейного уравнения теплопроводности предлагается компактная трехслойная разностная схема на девятиточечном шаблоне с тремя узлами с индексом п +1 на верхнем временном слое tn+l, тремя узлами с индексом п на промежуточном слое Ьп и тремя узлами с индексом п — 1 на нижнем временном слое Используя асимптотическое разложение, получаем схему повышенного порядка точности (0(т2 + к4)) [11]
МП + 1 + ßn-l h2 iV Ж
h2 ду фП
Жп = ГЛ0--—--— АОт1 - фп -— , (15)
t s ° 2 12 0 zn 12 0 Z n
в которой
(Ф)П = У
Y - 1
д Ф
YPgR2 \ dy
дФ ' дУ
дФ • yR дФ .
— - DPg + У-— , (16)
у=11 -у R -у
жп =
rßn+l _ Ж'п-1
2Т :
Zn =
D(pn ,т n)
n\2 '
(Rn)
(17)
где к — пространственный шаг по у; т — временной шаг; Л0 — обычная аппроксимация двойного дифференцирования по координате у; фп — "источниковый" член в уравнении (3) после всех преобразований.
r
Трехслойная компактная разностная схема повышенного порядка точности (0(т2 + h4)) для неоднородного нелинейного волнового уравнения (1) была реализована в работе [12]:
/ Pn + 1 _ 2 Pn + Pn— 1 \
A -^-J = ЛХРП + / (18)
A = E + 12 ЛХ, / = /П + 12 Л/,
где E — единичная матрица; ЛХ — центральная аппроксимация двойного дифференцирования по x; /s описывает монопольные источники Лайтхилла в правой части уравнения (1). Этот член содержит производную давления (или плотности) по времени в каждом узле пространственной сетки. Здесь используется формула дифференцирования назад четвертого порядка точности. Волновое уравнение интегрируется как эволюционная краевая задача.
Для контроля корректности вычислений теплообмена между газовой полостью и окружающей средой используется закон сохранения массы m пузырька на каждом временном шаге
1
m = 3PR3/ Т(Ууdy. (19)
о
Во всех приведенных результатах расчетов относительная ошибка в определении массы не превышала 1 % для тепловой задачи одиночного пузырька. Для второй задачи о распространении длинноволновых возмущений относительная ошибка на порядок (и более) была меньше, чем в предыдущей задаче.
Уравнения (2) и (4) интегрировали явным методом Рунге-Кутты четвертого порядка с автоматическим выбором шага. В целом численная схема имеет четвертый порядок по пространству и второй по времени 0(т2 + h4)). Пакет программ написан на языке Фортран-77 с применением технологии параллельного программирования MPI. Расчеты проводились на многопроцессорном вычислительном комплексе МВС-1000М и впервые получены результаты с использованием исследуемой математической модели.
Обсуждение результатов, Общие для всех расчетов параметры системы вода-воздух: P0 = 100 кПа, TL = 293 K, 7 = 7/5, ро = 103 кг/м3, а = 0,072 Н/м, р = 10—3 кг/(м-е), c = 1500 м/с, шо = 1 уЗтРОТРо — собственная частота пузырька Миннаэрта. Зависимость коэффициента теплопроводности воздуха от температуры задана в виде к(Т) = AT + + B, где A = 5,52810—5 Дж/(м-с-К2), B = 1,165 • 10—2 Дж/(м-с-К).
Для тестирования программы повторены расчеты, аналогичные работе [6], с параметрами PA = 60 кПа, R0 = 50 мкм и ш/ш0 = 0,78.
- . 1 \2 1
11 I 2 .
—-*io\ XL« 5 6 7
V
Ш _р_
Тъ р 100 -,2 р
9
1 /•" ч
s4
■ J 1........."
О 1.047 2.094 3.141 4.188 5.235 6.282 0 1.047 2.094 3.141 4.188 5.235 6.282
а б
Рис. 1. Профили безразмерных радиуса, температуры и давления пузырька в поле гармонической волны:
а — радиус (1) и температура газа (2) в центре пузырька как функции безразмерного времени; Ь — давление газа (1) и безразмерное давление гармонической волны (2)
Получено хорошее, не только качественное, но и количественное, согласие в пределах 20 %.
Тепловая динамика пульсирующего пузырька. Первая рассмотренная в данной работе задача — изучение распределения температуры в одиночном теплопроводящем пузырьке на стадиях сжатия и разрежения под действием вынуждающего гармонического звукового поля.
Безразмерные временные профили температуры в центре (тонкая линия) и радиуса (жирная линия) пузырька за один полный период стационарных вынужденных колебаний приведены на рис. 1, а (показан один период колебаний). Здесь на кривой 2 точками и цифрами показаны выбранные моменты безразмерного времени, а соответствующие им температурные профили показаны на рис. 2. Сравнение кривой нормированного радиуса с аналогичной кривой на рис. 3 из работы [6] показывает хорошее качественное и количественное согласие. Кривая температуры центра пузырька в работе [6] не приводится, но температурные профили для точек А, В и С, указанные на графике радиуса
Рис. 2. Радиальные профили температуры газа соответствующие фазе расширения (а) и фазе сжатия (б) пузырька. Номера кривых соответствуют моментам безразмерного времени, указанным оцифрованными точками на кривой 2 рис. 1, а
Рис. 3. Образование отрицательных солитонов радиуса пузырьков (сплошные линии) и положительных солитонов давления (штриховые линии)в газожидкостной смеси без диссипации из начального синусоидального возмущения
(рис.3) приведены в работе [6] на рис. 12. Эти точки соответствуют нашим точкам ¿д, ¿1, ¿6 и соответственно температурным профилям 9, 1 и 6, показанным на рис. 2. Нетрудно убедиться, что они также хорошо согласуются. Хорошее совпадение получено при сравнении кривой давления (кривая 1) на рис. 1, Ь с аналогичной кривой давления, приведенной на рис.4 работы [6]. Характерно, что учет межфазного теплообмена приводит к резкому уменьшению пикового давления в пузырьке по сравнению с политропическим законом поведения газа. На рис. 1, Ь кривой 2 показана безразмерная амплитуда гармонической вынуждающей силы. Аналогичные сравнения для другого радиуса пузырька Я0 = 10 мкм также показали хорошее совпадение результатов.
Отметим, что на значительном промежутке одного периода пульсации пузырька температура в центре ниже температуры на стенке (между временными точками ¿4 и ¿8). Минимальная температура в центре пузырька (точка ¿5) оказалась равной Тт1П(г = 0) = 0,68Ть, а соответствующее этому моменту времени давление — равным Рт1П = 0,78Ро.
Эффект теплопроводности, который приводит к существенно немонотонным температурным распределениям внутри пузырька проиллюстрирован на рис. 2. Стадии расширения соответствует рис. 2, а. Температурные распределения 1-2 (пунктирные кривые, левая вертикальная ось) — монотонные, а распределения 3-4 (сплошные линии, правая ось) — немонотонные. Кривая 5 соответствует моменту времени ¿5, при котором достигается минимум температуры в центре пузырька. Он предшествует моменту времени ¿6, при котором достигается максимум радиуса. Очевидно, что эти точки не совпадают, так как газ получает теплоту из жидкости. Это объясняется тем, что при до-расширении пузырька до своего максимального значения Я за время ¿6 — ¿5 газ получает теплоты из жидкости больше, чем отдает, когда охлаждается из-за расширения. Этот факт был отмечен и в работе [6].
Далее рассмотрим вновь обнаруженные эффекты, не отмеченные в работах [6-7]. На рис.2,б показаны температурные распределения на стадии сжатия. Кривые 6, 7 и 9 (сплошные линии) относятся к правой вертикальной температурной оси, а кривые 10-11 (пунктирные линии) относятся к левой температурной оси. Для профиля 7 максимум температуры достигается примерно на расстоянии Я/4 от границы жидкости. Этот максимум по мере сжатия плавно перемещается к центру пузырька. Интересно отметить, что при сжатии пузырька момент времени, при котором достигается максимальная температура в центре, также предшествует моменту времени, при котором достигается минимальный радиус. Этот эффект трудно заметить на графиках на рис. 1, а, но он четко выражен на рис. 6 и связан с тем, что при сжатии до минимального радиуса пузырек отдает теплоты жидкости больше, чем получает, когда прогревается за счет сжатия. По-видимому, выявление этих эффектов связано с использованием более точной численной схемы.
Аналитическая модель нелинейных волн в приближении КдФ. Во второй задаче рассмотрено распространение акустических волн в газожидкостной смеси с содержанием газа = 10-4, радиусом пузырьков Д0 = 0,5 мм, амплитудой звукового возмущения = 5 кПа с частотой ш = ш0/20.
Прежде чем переходить к изучению диссипативной системы, мы исследовали гамильтонову подсистему, содержащуюся в полной системе (1)-(6). В современной физике нелинейных явлений в гамиль-тоновых системах хорошо известно, что образование стационарных уединенных волн возможно лишь при уравновешивании слабой нелинейности и дисперсии. Существует несколько физических ситуаций, приводящих к такому равновесию, и имеются соответствующие им нелинейные уравнения. В исследуемой ниже гамильтоновой модели жидкости с пузырьками газа слабая нелинейность ЯЯХ, уравновешиваемая слабой дисперсией Яжжж, приводит к уравнению КдФ.
Хорошо известное уравнение КдФ для смеси жидкости с однородно распределенными пузырьками газа было получено Вейнгарденом в предположении, что колебания пузырька линейные и изотермические, а амплитуда гидродинамических возмущений конечна [2]. Следует отметить, что радиальное движение пузырька является изначально нелинейным, и на основании работы [21] можно утверждать, что и на низких частотах ш ^ ш0 нелинейность среды определяется также нелинейностью пузырьков, которая превосходит гидродинамическую нелинейность на два порядка.
Ранее, в работе [5] было обнаружено, что начальное возмущение в жидкости с нелинейными, но адиабатическими пузырьками газа,
распадается на последовательность волн, близких к уединенным. На рис. 3 показан пример пространственной эволюции гармонического (на входе) возмущения. Волновое поле относится к моменту времени t = 44 мс. Верхняя координата X относится к профилю радиуса, изображенному тонкой сплошной кривой, а нижняя абсцисса — к профилю радиуса R, изображенному жирной кривой, и профилю давления, изображенному пунктирной кривой. Результаты моделирования показывают, что синусоидальное возмущение вначале ведет себя подобно обычной нелинейной волне в газовой динамике: крутизна фронтов увеличивается и формируется N-образная волна. На расстоянии 13 м начинает сказываться дисперсия, обусловленная наличием пузырьков. На тонкой кривой профиля радиуса видны зародыши солитонов, возникающие в течение отрицательного полупериода. Дальнейшая эволюция этих зародышей приводит к рассыпанию синусоидальной волны на совокупность солитонов, которые показаны цифрами 1.. .4 на профиле радиуса (жирная линия), а также на профиле давления (штриховая линия); вершины их лежат на одной прямой. Сравнение с результатами работы [22] указывает на полное качественное совпадение этих процессов и на существование своего уравнения КдФ в рассматриваемой системе. Действительно, если пренебречь всеми видами диссипации, поверхностным натяжением и учесть, что концентрация пузырьков в единице объема n = const, ^ ^ 1, R/с ^ 1, то в этом частном случае полная система (1)-(6) выглядит так:
1 dP dV
рс2 dt dx
+ — = 4nR2nR; (20)
дК дР
"Ж + & =0 (21)
3 • 1
ДД + -И2 = -(Рд - Р); (22)
2 р
рд = Р° (IГ' (23)
где V — гидродинамическая скорость частиц жидкости. Из этой га-мильтоновой системы [23] с использованием метода многомасштабных разложений получено уравнение КдФ для возмущения радиуса пузырьков, совершающих нелинейные колебания [24].
Уравнение КдФ можно вывести иначе, используя формализм нелинейной волновой динамики. Оно естественно возникает как главный член аппроксимации во всех консервативных волновых системах со слабой нелинейностью и слабой дисперсией. Разложим искомые
функции в степенные ряды по малому безразмерному параметру возмущения е, который можно положить равным максимальному из отношений
е _ R - До Р - Po Р - Ро До Ро Ро
Тогда формальные разложения примут вид
V _ eV' + ... , R _ До + еД' + ... , Р _ Ро + еР' + .... (24)
Подставим разложение (24) в подсистему (20)-(23) и сохраним члены до второй степени е включительно. В результате для волн, бегущих в одну сторону вдоль характеристики x - cit _ const, в новых переменных т _ t и n _ x — cit после ряда преобразований получим искомое уравнение для микропеременной R' в виде (штрихи опущены)
3c c3
Дт - е—1(y +1)RRn + Дтп _ 0, (25)
2 До 2^о
где ci — скорость звука в низкочастотном приближении, определяемая из формулы 1/ci _ 1/c2 + ^ор/(7Ро). Коэффициенты при нелинейном и дисперсионном членах уравнения (25) определяются только дискретной газовой фазой — пузырьками. Знак минус при нелинейном члене полученного уравнения соответствует знаку амплитуды волны возмущения радиусов R (сплошные линии на рис. 3). Отрицательные консервативные солитоны получаются при условии баланса между нелинейным искажением профиля волны и его дисперсионным расплы-ванием. В нашем случае оба этих эффекта в газожидкостной суспензии связаны только с наличием пузырьков. Важно отметить, что в настоящей работе уравнение КдФ (25) выведено с использованием нового сильного механизма нелинейности, который кардинально отличается от механизма нелинейности, использованного в работе Вейнгардена [2].
Экспериментально проще измерять давление в жидкой фазе, являющееся переменной макроуровня, для которой также можно вывести уравнение КдФ. Между давлением газа и радиусом полости нетрудно установить связь: Р'д/Рдо _ -37Д'/До. В длинноволновом приближении динамическая нелинейность пузырька, определяемая вторым нелинейным членом слева в уравнении его радиального движения (22), существенно меньше нелинейности уравнения состояния газа (23), как показано в [21]. В силу медленных пульсаций микроинерцией пузырька также можно пренебречь, тогда из уравнения (22) следует Р'д _ Р'. Это приближение полностью подтвердилось в численных расчетах. Теперь уравнение (25) можно переписать для макропеременной Р'
(штрихи опущены) в виде уравнения КдФ:
С1 (7 +1) Рр + _е| 2 7 Ро П 2^,
Рт + ^^ -Pn + ^ Pnnn = 0, (26)
которое описывает солитоны положительной полярности возмущения давления жидкого компонента и показаны на рис. 3 (штриховая кривая).
Теперь не представляет труда сформулировать уравнение КдФ для остальных переменных смеси. Принимая во внимание связь между возмущениями скорости и давления, получаем
К + ККп + Кт = 0. (27)
Особого внимания заслуживает тот факт, что из гамильтоновой системы (20)-(23) [23] выведено уравнение КдФ, которое, как известно, также является гамильтоновой системой. Уравнение КдФ точно интегрируется с помощью метода обратной задачи рассеяния (МОЗР). Однако МОЗР не применим в случаях, когда в динамике системы важную роль играют диссипативные или межфазные процессы. В таких ситуациях приходится прибегать к численным моделям и методам.
Численная модель диссипативных уединенных волн в жидкости с пузырьками газа. Интерес представляют реальные среды с диссипацией, в которых имеет место приток и отток энергии, а в результате их баланса возможно существование диссипативных солитонов. В связи с этим интересно проанализировать принципиальные особенности возникновения стационарных возбуждений в нелинейной диспергирующей среде с потерями. В настоящей работе на примере жидкости пузырьковой структуры продемонстрированы качественные отличия нелинейных волновых процессов в средах с диссипацией и без нее.
Эволюция звуковой волны в смеси жидкости с теплопроводными пузырьками показана на рис. 4. Видны особенности процесса формирования уединенной волны. Вначале за одно колебание пузырька количество теплоты, отданное нагретым газом в жидкость, больше количества теплоты, поглощенного охлажденным газом из жидкости, за счет чего и происходит рассеяние тепловой энергии и демпфирование колебаний пузырька газа. В результате, как видно из графиков, амплитуда первых волн сильно убывает, но этот процесс не приводит к полному затуханию, когда имеются нелинейные и дисперсионные эффекты. Начиная с некоторой амплитуды, затухание первого гребня волны прекращается. В процессе эволюции он объединяется со вторым гребнем (как видно на штрихпунктирной кривой рис. 4, а, которая соответствует моменту времени £ = 50,6 мс), а тот — в свою очередь с
Рис. 4. Первый, второй (а) и третий (Ь) этапы эволюции исходной звуковой волны в нелинейной, диспергируюшей и диссипативной газожидкостной системе
третьим гребнем, как показано на сплошной кривой в момент времени £ = 68,5 мс. В результате к моменту времени £ = 110 мс образуется предельная конфигурация — стационарная уединенная волна.
Уединенная волна радиусов показана на рис. 4, б в момент времени £ = 212,4 мс (тонкая сплошная кривая, вторая правая вертикальная ось). У пузырьков, попавших в область, занятую уединенной волной, суммарные приток и отток теплоты уравниваются. Это возможно потому, что колебания пузырька нелинейные. При этом фаза сжатия пузырька, как известно, короче фазы расширения, в результате чего общее количество теплоты, отданное пузырьком в фазе сжатия, может быть равно общему количеству теплоты, полученного им в фазе расширения. В этом случае не происходит рассеивания внутренней энергии газа и пузырек находится в тепловом равновесии с окружающей жидкостью. Следовательно, амплитуда уединенной волны стабилизируется за счет фоновой теплоты жидкого термостата.
Действительно, на рис. 4, б показан результат взаимодействия динамических и тепловых процессов при осцилляциях газовых пузырьков. Пунктирная кривая относится к градиенту температуры, жирная
сплошная линия — к радиальной скорости стенки (левая вертикальная ось, точка означает производную по времени) в момент времени £ = 212,4 мс. Тепловой поток ( связан с градиентом температуры на стенке соотношением ( = —к(Т)(дТ/дг)|г=д. Как следует из графиков, у пульсирующего пузырька знак градиента температуры на стенке тесно коррелирован со знаком скорости стенки, т.е. сдвиг фазы между ними практически отсутствует. Оценки показывают, что характерное время тепловой релаксации больше периода сжатия-расширения пузырька. При расширении газ в окрестности стенки охлаждается, так как теплопроводность не успевает компенсировать это охлаждение пристенных слоев газа. В результате градиент температуры положителен и теплота втекает из жидкости в газ. И наоборот, при сжатии газ нагревается, скорость отрицательна, градиент температуры отрицателен и газ отдает теплоту в жидкость — это эффект движущейся границы [24].
В эксперименте измерение сигналов датчиками проводят в выбранной точке. Рассмотрим осциллограмму волнового процесса, снятую в пространственной точке, например X = 142,04м и показанную на рис.5. Пока проходит уединенная волна давления (кривая 2), радиусы пузырьков (кривая 1) меньше начального значения. Как следует из рис. 5, знак градиента температуры зависит от фазы колебания пузырька, т.е. сжатия или расширения. На переднем фронте импульса давления пузырек сжимается, производная температуры на поверхности стенки (кривая 3) у него отрицательна и теплота отдается газом в жидкость, а на спаде импульса давления, наоборот, — пузырек начинает расширяться, производная температуры на стенке становится
Рис. 5. Временные профили в точке X = 142,04 м безразмерного радиуса Я (1), давления Р (2) и безразмерного градиента температуры (дТ/ду)\у=\ на границе
фаз (3)
положительной и теплота возвращается жидкостью в газ. При этом средний за цикл тепловой поток равен нулю. Для подтверждения этого была вычислена тепловая энергия, отводимая и подводимая через сферическую поверхность газовой полости, которая определяется интегралом
г г
= / / д(г)|г=л^ = -4^/ Л2 0 £ 0
dt.
r=R
Кривые Л(£) и (дТ/дг)|г=д, входящие в этот интеграл, показаны на рис.5. Как и следовало ожидать, интеграл практически равен нулю. Если принять во внимание, что жидкость считается термостатом, то можно доказать, что этот численный результат не противоречит второму началу термодинамики. В результате достигнутого баланса между притоком и оттоком теплоты в газовую фазу и равенства нелинейного кручения и дисперсионного расплывания профиля волны возникает квазистационарная уединенная волна, распространяющаяся без затухания с постоянной скоростью, как сквозь прозрачную среду. Очевидно, что возможность термодинамической обратимости межфазного теплообмена связана с принятым в данной работе приближением об изотермичности жидкости. Данный эффект обнаружен впервые, благодаря более точной численной схеме, и он важен для анализа волновых процессов в жидкости с пузырьками.
Обратим также внимание на асимметрию формы диссипативной уединенной волны на рис. 5 в отличие от симметричных солитонов уравнения Кортевега-де Фриза-Бюргерса или Буссинеска. Эту асимметрию профиля "солитона" можно наблюдать во многих экспериментах. Например, на рис. 2 работы [25] проведено сравнение экспериментальных профилей с расчетными профилями уравнения Буссинеска. Солитон образуется быстрее, если вместо звуковых волн рассмотреть распространение импульсов, как это имеет место в [25]. Поэтому эксперимент относится к более короткому временному интервалу, но механизм возникновения асимметрии один и тот же.
Интересные и полезные результаты можно извлечь из рассмотрения пространственного температурного профиля пузырьков, входящих в уединенную волну, который представлен на рис. 6 в момент времени £ = 194,8 мс. Из него следует, что в пространственной области, занятой уединенной волной, все пузырьки находятся в сжатом состоянии Л < 1 (кривая 2). Несмотря на это, кривая температуры в центре пузырьков Т0 (кривая 1) разбивается на две примерно равные подобласти или полуволны. Такое волнообразное поведение (с отрицательной полуволной) не может соответствовать физически допустимому, если
120 140 160 m 180
Рис. 6. Пространственные профили температуры T0 (1), радиуса пузырьков R (2) и безразмерного градиента температуры (дТ/ду)\y=i (3)
принять политропический закон поведения газа, так как в этом случае любому состоянию с меньшим радиусом должна соответствовать большая температура T0 > TL, а в изотермическом случае T0 = TL. Волнообразное поведение T0 возможно только при немонотонном температурном профиле внутри пузырьков. Обратим внимание на то, что процесс очень медленный, но поведение газа не изотермическое, его температура не постоянна и не равна TL вопреки интуитивно принятому представлению для низкочастотного приближения во всех литературных источниках.
Соответственно, в первой подобласти (т.е. на переднем фронте волны сжатия) температура центров выше температуры окружающей жидкости, а во второй подобласти (на заднем фронте волны сжатия) температура центров ниже температуры окружающей жидкости. В первой полуволне у пузырьков градиент температуры на стенке (кривая 3 на рис. 6) отрицателен и все пузырьки этого гребня отдают теплоту в жидкость, а во второй полуволне градиент температуры положителен и все пузырьки, попавшие в этот второй гребень, получают теплоту из жидкости.
Было проведено также моделирование с другими частотами и амплитудами возбуждающего акустического возмущения. Качественно картина эволюции акустической волны осталась прежней. В результате выявлено, что уменьшение только одного параметра частоты в два раза (ш = ш0/40) при неизменных прочих параметрах задачи привело к росту амплитуды уединенной волны почти вдвое. Данный результат можно объяснить с позиций второго начала термодинамики. Чем медленнее происходит цикл сжатия и расширения пузырька, тем он ближе к обратимому процессу, а значит, тем меньше потери энергии на первом и втором этапах формирования уединенной волны.
Прежде чем перейти к заключению, отметим, что данная модель обобщена на случай волн конечной амплитуды, т.е. в уравнениях гидродинамики учтены члены до второго порядка малости включительно. Численное моделирование показало, что в длинноволновом приближении результаты, полученные по исходной полулинейной (1)-(6) и обобщенной полностью нелинейной моделям, практически совпадают [26, 27]. Это еще раз подтверждает доминирующую роль нелинейности пузырьков по сравнению с гидродинамической нелинейностью, как было впервые установлено теоретическими оценками в работе [21].
Заключение. Проведенный анализ показал, что солитоны в консервативной пузырьковой жидкости возникают из-за баланса между эффектами возрастания крутизны волн, связанными с нелинейностью пузырьков, и расползанием волн вследствие дисперсии, связанным также с пузырьками газа.
Проведено обобщение двухуровневой волновой динамической модели жидкости с распределенными пузырьками на случай теплообмена между фазами. Разработаны компактные разностные схемы четвертого порядка точности для нелинейного уравнения теплопроводности и неоднородного нелинейного волнового уравнения.
Численное моделирование на основе уравнений (1)-(6) динамики искажения профиля гармонической на входе волны позволяет выделить три характерных этапа. На первом этапе существенны диссипа-тивные эффекты, поэтому гармоническая волна действительно сильно затухает по амплитуде, но до определенного уровня. На втором этапе доминируют нелинейные и дисперсионные эффекты, а диссипация несущественна. На этом этапе в результате конкуренции между механизмом, связанным с нелинейным увеличением крутизны и механизмом дисперсионного расплывания волны, формируется стационарная уединенная волна. На третьем этапе уединенная волна, после того как она приняла свою окончательную форму, фактически распространяется без затухания и с постоянной скоростью. На третьем этапе жидкая и газовая фазы обмениваются равными порциями энергии, если жидкость считать термостатом. Выполненные исследования позволяют сделать вывод, что газожидкостная смесь с пузырьками является прозрачной для уединенных волн. Установлена асимметрия формы уединенной волны в диссипативной пузырьковой среде.
Общепринятый в настоящее время феноменологический учет диссипации в виде некоторых специально подобранных коэффициентов "эффективной" вязкости или "эффективных" чисел Нуссельта в общем случае нельзя признать адекватным механизмом учета межфазного теплообмена.
Автор выражает глубокую сердечную благодарность Ю.И. Шокину
за научные консультации и поддержку данной работы, Г.А. Михайлову
за постоянное внимание и всестороннюю поддержку и Ю.Н. Морокову
за плодотворное обсуждение компактных разностных схем.
СПИСОК ЛИТЕРАТУРЫ
1. Наугольных К. А., Островский Л. А. Нелинейные волновые процессы в акустике. - М.: Наука, 1990.
2. WijngaardenvanL. On the equations of motion for mixtures of liquid and gas bubbles // J. Fluid Mechanics. - 1968. - T. 38. - P. 465-474.
3. Котельников И. А. Эффект эха в жидкости с пузырьками газа // Известия вузов. Радиофизика. - 1983. - T. 26. - №. 10. - C. 1227-1234.
4. Fedotovsky V. S., Vereschagina A. V., Derbenev A. V. Experimental research of resonance sound dispersion in bubbly media // Proc. Int. Conf. "Fluxes and structures in fluids", Moscow, Russia, June 20-23, 2005. Ed. by Yu. D. Chashechkin and V. G. Baydulov. Selected Papers, pp. 119-123. The Conference is dedicated to 250th Anniversary of the M.V. Lomonosov Moscow State University.
5. K i m D. C. A nonlinear wave dynamical model for two-phase flows and its numerical solutions // Computer Physics Communication. - 2002. - T. 147. - P. 526529.
6. Prosperetti A., Crum L. A., Commander K. W. Nonlinear bubble dynamics // J. Acoust. Soc. Am. - 1988. - V. 83. - № 2. - P. 502-514.
7. Kamach V., Prosperetti A. Numerical integration methods in gas-bubble dynamics // JASA. - 1989. - V. 85. - № 4. - P. 1538-1548.
8. K a m e d a M., M a t s u m o t o Y. Shock waves in a liquid containing small gas bubbles // Phys. Fluids. - 1996. - V. 8. - № 2. - P. 322-335.
9. А й д а г у л о в P. P., Х а б е е в Н. С., Ш а г а п о в В. И. Структура ударной волны в жидкости c пузырьками газа с учетом нестационарного межфазного теплообмена // Прикл. мех. и техн. физика. - 1977. - № 3. - C. 67-74.
10. R o m e n s k i E. I., T o r o E. F. Compressible two-phase flows: two-pressure models and numerical methods // Computational Fluid Dynamics Journal. - 2004. -Vol. 13. - P. 403-416.
11. K i m D. C. Application of fourth-order compact finite difference schemes for two-phase bubbly flow problems // Proc. 4th International Conference on Computational Heat and Mass Transfer. Paris-Cachan, France, May 17-20, 2005. ISBN: 2-74300801-6. CD ROM, 562.pdf.
12. Kim D. C. Hydraulic transients in two-phase bubbly flows // Proc. of the 8th Korea-Russia Intern. Sym. on Science and Technology, Russia, Tomsk, June 26-July 3, 2004, vol. 2, pp. 225-229. ISBN: 0-7803-8383-4; Library of Congress USA, Card Number: 2004101655.
13. Нигматулин Р. И. Динамика многофазных сред. - М.: Наука, 1987.
14. B r e n n e n C. Cavitation and bubble dynamics. Oxford University Press, 1995.
15. А л е к с е е в В. Н., Р ы б а к С. А. Распространение стационарных звуковых волн в пузырьковых средах // Акустический журнал. - 1995. - T. 41. - № 5. -C. 690-698.
16. Кузнецов В. В., Н а к о р я к о в В. Е., П о к у с а е в Б. Г., Шрей-б е р И. Р. Жидкость с пузырьками газа как пример среды Кортевега-де Фриза-Бюргерса // Письма в ЖЭТФ. - 1976. - T. 23. - Вып. 4. - C. 194-197.
17. Накоряков В. Е., П о к у с а е в Б. Г., Шрейбер И. Р. Волновая динамика газо- и парожидкостных сред. - М.: Энергоатомиздат, 1990.
18. А в е р ь я н о в М. В., Б а с о в М. С., Х о х л о в а И. Р. Стационарные и квазистационарные волны в диссипативных системах четного порядка // Акустический журнал. - 2005. - T. 51. - C. 581-588.
19. Куропатенко В. Ф. Модель многокомпонентной среды // ДАН РАН. -2005. - Т. 403. - С. 761-763.
20. C r i g h t o n D. G., F o w c s W i 11 i a m s J. E. Sound generation by turbulent two-phase flow // J. Fluid Mech. - 1969. - V. 36. - P. 585-603.
21. Заболотская Е. А. Генерация второй гармоники звуковой волны в жидкости с равномерно распределенными воздушными пузырьками // Акустический журнал. - 1975. - T. 21. - C. 934-937.
22. Z a b u s k y N. J., K г u s k a 1 M. D. Interaction of "solitons" in a collisionless plasma and the recurrence of initial states // Phys. Rev. Lett. - 1965. - V. 15. -P. 240-243.
23. C a f l i s c h R. E., M i k s i s M. J., P a p a n i c o l a u G. C., L u Ting. Effective equations for wave propagation in bubbly liquids // J. Fluid Mechanics. -1985. -V. 153. - P. 259-273.
24. Ким Д. Ч. Возможность существования стационарных уединенных волн в диссипативной газожидкостной среде // Письма в ЖТФ. - 2006. - T. 32. - Вып. 5. -C. 38-46.
25. Д о н ц о в В. Е. Распространение волн давления в газожидкостной среде кластерной структуры // Прикладная механика и техническая физика. - 2005. -T. 46. -№ 3. - С. 50-60.
26. K i m D. C. Fully nonlinear wave dynamical model for bubbly liquids containing thermal conductive bubbles // 6th International Conference on Multiphase Flow, ICMF 2007, Leipzig, Germany, July 9-13, 2007, Paper No S6-Fri-B-66.
27. К и м Д. Ч. Физическая природа акустических солитонов в жидкости с распределенными пузырьками газа // Доклады РАН. - 2008. - Т. 418. - № 5. -C. 619-623.
Статья поступила в редакцию 17.01.2007
Дин Чер Ким родился в 1946 г. окончил Новосибирский государственный технический университет в 1970 г. Канд. физ-мат. наук, доцент кафедры физики Сибирского государственного университета телекоммуникаций и информатики (СибГУТИ). Автор более 40 научных работ в области нелинейных волновых процессов в газожидкостных системах и волоконно-оптических линиях связи.
D.Ch. Kim (b. 1946) graduated from the Novosibirsk State Technical University in 1970. Ph. D. (Phys.-Math.), assoc. professor of physics department of the Siberian State University of Telecommunications and Information Technology. Author of more than 40 publications in the field of nonlinear wave processes in gas-liquid systems and fiber-optical communication links.
