Предварительное численное моделирование распространения звука в жидкостях c газообразными включениями

Распространение звука в жидкостях, содержащих газовые пузырьки, широко применяется во многих областях. Например, пузырьки выполняют роль контрастных агентов для улучшения качества визуализации в ультразвуковой диагностике в медицине, в промышленности и химических технологиях для визуализации трехфазной среды «нефть - газ - вода», кавитационных технологиях [1–5]. Ценность эффектов жидкостно-пузырьковой системы в том, что наличие пузырьков резко изменяет акустические свойства среды: скорость звука может варьироваться в пределах до нескольких сотен метров в секунду, возникают сильное рассеяние и затухание, локальные амплитуды давления значительно возрастают вблизи пузырьков резонансного размера и т. д. На эффекте акустической кавитации основано направление сонохимического синтеза – получение наноструктурированных материалов при комнатной температуре и давлении с коротким временем реакции (металлы, полупроводники, графен, композиты с оксидом графена, стабильные наноразмерные структуры из биомолекул) [6–8].

Физическая проблема адекватного описания таких процессов заключается в том, что в ультразвуковом диапазоне частот длины волн становятся сопоставимы с характерными размерами газовых включений [9, 10]. Это приводит к выраженным эффектам дисперсии и интенсивного резонансного рассеяния [11, 12]. Жидкая среда с пузырьками теряет свойства акустической однородности, а локальные волновые явления на границах раздела фаз (включая микроскопические изменения плотности, вязкую диссипацию и образование акустического пограничного слоя) начинают оказывать определяющее влияние на общее распределение акустического давления [13–15]. Классические аналитические модели зачастую ограничены жесткими допущениями (например, условием термодинамического равновесия или пренебрежением пространственной интерференцией). Это создает сложности для их прямого использования при расчетах в реальных условиях [16, 17].

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

Для достижения поставленной цели проведено сопоставление классических моделей, рассмотрены основные режимы взаимодействия ультразвука с газовыми включениями и выполнено качественное численное моделирование распределения акустического давления. Такой аналитический базис необходим для построения адекватных граничных условий и последующего полномасштабного численного моделирования ультразвуковых полей методом конечных элементов в программном комплексе COMSOL Multiphysics [11, 18].

Математические модели ультразвуковых полей в гетерогенных средах

Одной из фундаментальных теорий для расчета скорости звука в жидкостях с газообразными включениями является классическая модель Вуда, в которой реализована идея механического и термодинамического равновесия между газовой и жидкой фазами. Модель предполагает, что скорость звука в неоднородной смеси определяется усредненной плотностью и сжимаемостью ее индивидуальных компонентов [16, 19]. Для разбавленных смесей, когда объемная доля пузырьков мала и теплоемкость несущей жидкости близка к единице, эта модель обеспечивает высокую точность прогнозирования [10, 16]. На относительно низких частотах ультразвука, пока газовые включения «успевают» деформироваться и подстраиваться под фазы сжатия и разрежения проходящей волны без существенного запаздывания, модель Вуда остается надежным аналитическим инструментом [11, 16, 20, 21].

Однако у модели Вуда есть существенные физические ограничения. Во-первых, будучи равновесной, она не учитывает частотную дисперсию и диссипацию энергии, что делает ее малоприменимой для анализа ультразвуковых волн высокой частоты [10, 19, 21, 22–26]. Во-вторых, модель начинает терять адекватность, когда повышается концентрация газовой фазы или в случаях, когда пузырьки содержат твердые или вязкие примеси и начинают вести себя подобно недеформируемым сферам [16, 22]. Также в ней полностью игнорируется теплопередача между газом и окружающей жидкостью, которая при определенных режимах вносит определяющий вклад в затухание волны [27, 28].

Принципиально иной подход к описанию волновых явлений в двухфазных средах предложил Е. Кеннард [17]. В рамках его модели вводится важное уточнение: газовые пузырьки рассматриваются как динамические осцилляторы, совершающие радиальные пульсации под действием переменного акустического давления [17, 27]. Благодаря учету собственной резонансной частоты пульсирующего пузырька (определяется его радиусом, давлением и упругостью газа) модель Кеннарда позволяет описывать зависимость скорости звука от частоты волны, а также резкое возрастание коэффициента затухания вблизи области резонанса [17, 27]. Это приближает ее к реальным физическим условиям распространения ультразвука на средних и высоких частотах [9, 17, 26].

Тем не менее и модель Кеннарда имеет ряд допущений, ограничивающих ее универсальность. Затухание колебаний в рамках данного подхода принимается постоянным, хотя в действительности диссипативные процессы носят выраженный нелинейный характер [18, 27, 28]. Модель Кеннарда не рассматривает поступательное движение пузырьков как целого под действием сил акустического излучения (сил Бьеркнеса), а также пренебрегает теплопроводностью газа и детальным расчетом вязкостных потерь на межфазной границе, которые заметно изменяют характеристики волнового поля [28 - 31]. Для преодоления этих ограничений требуется учет нелинейных уравнений динамики (таких как уравнение Рэлея – Плессета) и тепловых эффектов Просперетти [18, 28, 32]. Дальнейшее развитие этот подход получил в модели Коммандера - Просперетти [33], учитывающей термическую и вязкую диссипацию через комплексный коэффициент сжимаемости смеси, что снимает часть ограничений модели Кеннарда, отмеченных выше.

Материалы и методы исследования

Для проведения численного анализа реальных ультразвуковых полей в трехмерных геометриях феноменологические эффекты и классические подходы собираются в дифференциально-волновые модели, реализуемые в программном комплексе COMSOL Multiphysics [11, 18, 34]. Математическим аппаратом в данном случае служит линейное волновое уравнение Гельмгольца, адаптированное для многофазных систем [1, 18, 35, 36]:

V' (-к^р))--^ = Qm ,                  (1)

где p - акустическое давление, р t - эффективная (средняя) плотность гетерогенной среды, являющаяся функцией плотности жидкой и газовой фаз, а Q m – источник массы (или звука) [1]. Параметр K eq представляет собой эквивалентное волновое число смеси, которое рассчитывается с учетом фазовой скорости звука в несущей жидкости, частоты ультразвукового возбуждения и объемной концентрации пузырьков [1, 11, 37].

Решение уравнения (1) методом конечных элементов в COMSOL Multiphysics дает возможность уйти от пространственной однородности и рассчитывать интерференционные картины, а также распределение акустического давления на границах раздела сред с учетом геометрии излучателей и свойств стенок рабочих емкостей [11, 18, 38–44].

Анализ режимов взаимодействия ультразвука с газовыми включениями

В результате сопоставления теоретических моделей и анализа условий прохождения упругих волн через неоднородные среды мы остановились на физическом параметре, определяющем характер изменения акустического поля – отношении длины волны ультразвука к эффективному размеру (на данном этапе к диаметру) пузырька ( λ/d ). В зависимости от величины этого отношения и наличия резонансных явлений, взаимодействие ультразвука с пузырьками газа в жидкости можно разделить на характерные режимы [9, 10, 16].

• 1.    Коротковолновый (геометрический) режим (λ/d≲1)

• 2.    Дифракционный режим (λ/d≈1…10)

• 3.    Длинноволновый (рэлеевский) режим вне резонанса (λ/d≫1, например, λ/d>10)

• 4.    Резонансный режим

В этом режиме длина волны много меньше размера газовых включений. Ультразвук подчиняется законам геометрической акустики (лучевое приближение). Волновой фронт не огибает пузырьки, а отражается и преломляется на границе раздела фаз, из-за чего пузырьки работают как акустические линзы и формируют за собой выраженные зоны «акустической тени».

Этот режим наблюдается, когда длина волны сопоставима с размером пузырька [11, 12]. Волновая природа излучения играет ключевую роль: происходит сложная дифракция (огибание препятствия), а пузырьки начинают работать как мощные вторичные точечные источники. Рассеянные сигналы интерферируют с первичной волной и отражениями от стенок, формируя сложную картину чередования зон повышенного и пониженного давления [1, 9, 18]. В численной модели это проявляется в виде интенсивного направленного рассеяния и локального усиления давления на границах фаз [1, 11].

Когда длина волны существенно превосходит диаметр пузырьков и частота далека от резонанса, включения практически не искажают волновой фронт. Интенсивность рассеяния стремится к минимуму, пузырьки остаются акустически прозрачными [9, 27, 45]. Для моделирования таких сред корректно использовать равновесную модель Вуда, рассматривая смесь как квазиоднородную среду с усредненными параметрами [16, 19, 46].

Специфика газовых включений в жидкости заключается в наличии выраженного низкочастотного резонанса [47]. При совпадении частоты ультразвука с собственной частотой пульсаций пузырька (которая геометрически обычно соответствует λ / d ≫1) происходит качественная перестройка поля [9, 48]. Колеблющаяся граница фаз начинает активно поглощать и переизлучать энергию, сечение рассеяния возрастает на порядки [9, 17]. В этом режиме классические равновесные модели перестают работать, и требуется переход к динамическим моделям (Рэлея - Плессета, Кеннарда), учитывающим нелинейные радиальные пульсации газа [1, 17, 28, 49–52].

Результаты исследования и их обсуждение

На рисунке 1 приведены результаты качественного численного моделирования распределения акустического давления в жидкости с газовыми включениями. Расчет был выполнен для трех характерных частот: 10, 15 и 20 кГц.

На частоте 10 кГц реализуется длинноволновый режим ( λ/d ≫1). В этих условиях волновой фронт был искажен незначительно, а газовые включения оставались преимущественно акустически прозрачными.

На частотах 15 и 20 кГц ситуация менялась. Эти значения попали в зоны дифракционного и, возможно, резонансного рассеяния. Здесь уже наблюдались выраженные интерференционные эффекты. Модель демонстрировала формирование сложных зон с перепадами давления, а также локальные градиенты на поверхностях пузырьков [53].

Задача получения точных количественных характеристик на данном этапе не ставилась. Первоначальная цель состояла в том, чтобы качественно проверить возможность данной численной модели, подтвердить теоретические выводы о перестройке поля. Тем не менее, опираясь на визуальный анализ полей давления (рис. 1), были сделаны следующие заключения.

Во-первых, при переходе к 15 и 20 кГц градиенты давления вблизи поверхности пузырька явно усиливались. Это объяснялось формированием парциальных стоячих волн на границе раздела фаз. Система переходила в состояние, близкое к резонансному, что хорошо согласовывалось с предсказаниями модели Кеннарда [17, 27].

Во-вторых, интерференционная картина закономерно усложнялась. На более высоких частотах четко проступали «акустические тени» и зоны фокусировки. Данный эффект коррелировал с ожидаемым ростом интенсивности рассеяния и коэффициента затухания [9, 11, 54].

Рисунок 1 - Распределение акустического давления на стенки сосуда в жидкости с пузырьками газа на частотах 10, 15 и 20 кГц

Получение точных количественных данных потребовало отдельного параметрического исследования. Для этого было необходимо разместить виртуальные точечные датчики в расчетной области [34]. Подобный подход дал возможность оценить кратность усиления давления, точно определить собственную резонансную частоту по формуле Миннарта [47] (2) для заданного диаметра пузырька, а также верифицировать квадратичную зависимость затухания [37].

f o =

\

2nR ^

3ГР о Р ,

где f o - собственная резонансная частота колебаний (Гц), R - равновесный радиус пузырька (м), y — показатель адиабаты газа внутри пузырька (для воздуха y~1.4), P o — статическое гидростатическое давление окружающей жидкости (Па), р - плотность окружающей жидкости кг/м³).

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

Практические рекомендации по настройке численных моделей в COMSOL Multiphysics

Проведенный анализ позволлил сформулировать основные требования для настройки численных моделей в среде COMSOL Multiphysics:

• 1.    Выбор математического описания. Для длинноволнового режима λ/d ≫1 допустимо использовать упрощенные модели эффективной среды. Для резонансного и коротковолнового диапазонов расчет должен строиться на основе дифференциального уравнения Гельмгольца с явным заданием геометрии каждого пузырька как упругого или пульсирующего препятствия.

• 2.    Учет диссипации. В коротковолновом режиме модель должна учитывать вязкостные и тепловые потери на границе раздела сред, которые приводят к дополнительному затуханию из-за образования акустического пограничного слоя.

• 3.    Дискретизация расчетной области. Для адекватного воспроизведения интерференционных эффектов размер элементов конечно-элементной сетки должен жестко привязываться к минимальной длине волны в системе. Сгущение сетки обязательно должно быть предусмотрено на границах раздела фаз «жидкость - газ», где присутствуют максимальные градиенты давлений.

• Пример расчетной сетки приведен на рисунке 2.

Рисунок 2 – Расчетная сетка для метода конечных элементов; она имеет локальные сгущения у источника звука (круглая область слева), на границах поверхностей и вокруг пузырьков воздуха.

Сетка предназначена для моделирования распространения акустических волн в жидкостях с газовыми включениями

Заключение

Проведенный анализ показал, что определяющим параметром, задающим режим взаимодействия ультразвука с газовыми включениями, является отношение длины волны к диаметру пузырька ( λ/d ). На его основе выделены три режима – длинноволновый, резонансный и коротковолновый – каждый из которых потребовал своего подхода к описанию: от квазиодно-родной модели Вуда до динамического подхода Кеннарда и полного численного решения уравнения Гельмгольца.

На этой основе сформулированы практические рекомендации по настройке численных моделей в COMSOL Multiphysics: выбор математического описания в зависимости от λ/d , учет вязкотермической диссипации на границе раздела фаз и требования к дискретизации расчетной сетки вблизи пузырьков. Такая настройка модели заранее, до проведения полномасштабных расчетов, позволяет избежать избыточного сгущения сетки там, где это физически не оправдано, и тем самым экономить вычислительные ресурсы.

Систематизированные допущения и рекомендации задают фундамент для последующих исследований. В качестве следующих этапов планируется построение в COMSOL численных моделей нелинейной динамики пульсирующих пузырьков, нелинейных эффектов контролируемой акустической кавитации (в том числе на основе модифицированного уравнения Гельмгольца с распределенной нелинейной диссипацией [55], а также моделирование распределения ультразвуковых полей в жидких средах с нанодобавками.

Работа выполнена в рамках государственного задания СФУ (№ FSRZ-2020-0012) и гранта SECIHTI CBF-2025-G-470.