Моделирование колебаний капли жидкости, лежащей на вибрирующем недеформируемом основании

УДК 519.63:532.5

МОДЕЛИРОВАНИЕ КОЛЕБАНИЙ КАПЛИ ЖИДКОСТИ,

ЛЕЖАЩЕЙ НА ВИБРИРУЮЩЕМ НЕДЕФОРМИРУЕМОМ ОСНОВАНИИ

КУЗЬМИН И. М., САРМАКЕЕВА А. С., ЧЕРНОВА А. А.

Институт механики Уральского отделения РАН, 426067, г. Ижевск, ул. Т. Барамзиной, 34

АННОТАЦИЯ. Рассматриваются вопросы математического моделирования методом объема жидкости в ячейке (Volume of Fluid — VOF) колебаний капли жидкости, расположенной на вибрирующем с низкими частотами недеформируемом жестком основании. Исследуются вопросы получения начальной равновесной конфигурации капли, учета как непосредственно контактного угла в тройной точке подложка-жидкость-газ, так и его изменения в процессе колебаний капли. Определяется влияние метода учета поверхностного натяжения и угла смачивания на амплитуду колебаний капли жидкости. Численные результаты сопоставляются с известными экспериментальными данными.

КЛЮЧЕВЫЕ СЛОВА: капля жидкости, волны Фарадея, свободная поверхность, метод объема жидкости (Volume of Fluid - VOF), контактный угол.

ВВЕДЕНИЕ

Собственные и вынужденные колебания жидкости, как ограниченной стенками сосуда, так и свободно лежащей на различных поверхностях (капли жидкости различного объема), широко применяются в технике и медицине, а анализ таких колебаний и их характерных особенностей, соответственно, является актуальной задачей [1].

Отметим, что передача энергии от вибрирующего тела к расположенной на нем жидкости приводит к возбуждению поверхностных капиллярных волн - ряби Фарадея, динамика которой зависит как от параметров приложенной к телу нагрузки, так и от физических свойств жидкости. При этом, режимы вибрационного воздействия, обеспечивающие возникновение того или иного вида поверхностных эффектов в малом объеме жидкости, на сегодняшний день, остаются малоизученными. Известен ряд экспериментальных работ, посвященных вопросам исследования вибрационного распыления жидкости. Так в работе [2] рассматривается колебание капли жидкости, объемом 100 мкл, лежащей на вибрирующей с ультразвуковыми частотами нежесткой подложке. В качестве подложки используется тонкая металлическая диафрагма, также совершающая изгибные колебания. Изучаются процессы атомизации капли и параметры, влияющие на скорость и интенсивность атомизации жидкости. В работе [3] экспериментально исследуется межфазная динамика капли жидкости (в том числе поверхностные эффекты: осесимметричные стоячие и азимутные волны, образование и трансформация устойчивой узловой решетки волн, образование гребней на свободной поверхности), лежащей на вибрирующей с ультразвуковой частотой нежесткой диафрагме. Выявлены и указаны режимы колебаний, в том числе и амплитудно-частотные характеристики, на которых начинается процесс атомизации. В [4] приводятся результаты экспериментального исследования процесса распыления капли жидкости, свисающей с поверхности твердой недеформируемой вибрирующей пластины, в том числе и выявленные резонансные частоты колебаний пластины и минимальные значения амплитуд, при которых начинается процесс атомизации. В работе [5] экспериментально изучены вертикальные колебания капли, расположенной на вибрирующей несмачиваемой подложке с малым гистерезисом контактного угла. Экспериментально получены два типа колебаний свободной поверхности капли и выявлен переходный режим, также проанализировано влияние контактного угла и

гистерезиса на моду колебания жидкости. В экспериментальных работах [6 - 8] исследуются колебания жидкости малого объема (5 мкл), лежащей на недеформируемой вибрирующей с низкими частотами (до 800 Гц) подложке. Приводятся данные о колебаниях свободной поверхности капли в зависимости от положения вибрирующей подложки, также рассматриваются мгновенные профили капли и внутренние капиллярные течения. Однако, экспериментальные данные не позволяют выявить влияние всех параметров на изучаемый процесс, что значительным образом затрудняет создание теории ультразвукового распыления жидкости, учитывающей особенности всех протекающих подпроцессов. При этом, численные работы, направленные на изучение динамики капель, посвящены, в основном, задачам моделирования висячих капель жидкости, задачам взаимодействия двух и более капель между собой, вопросам моделирования падения капель одной жидкости в другую, а также исследованию свободных и вынужденных горизонтальных колебаний капли жидкости. При этом, вопросы математического моделирования процессов, протекающих в отдельной капле жидкости, лежащей на вертикально вибрирующей с низкими частотными и малыми амплитудами недеформируемой подложке, рассматриваются только в работах [2, 9]. Однако, в [9] рассматриваются вопросы численного моделирования динамики капли, свисающей с вибрирующего твердого стержня, методом конечных элементов с применением метода Галеркина. Заметим, что данные о тестировании и валидации используемых методов и алгоритмов в работе отсутствуют. В работе [2] исследуется колебание одиночной полусферической капли жидкости объемом 30 мкл на твердом стержне, при этом, возможное движение контактной линии в постановке не учитывалось. Движение стержня задается синусоидальным законом и, как и силы поверхностного натяжения, включается в уравнение Навье-Стокса. Полученная система уравнений решается проекционным методом маркера в ячейке (MAK), с использованием для дискретизации по пространству метода конечных объемов. Приведенные в [2] данные хорошо согласуются с экспериментальными, однако, в работе [2] приводятся результаты решения только плоской двумерной задачи, без учета влияния контактного угла на исследуемые процессы. Отметим, что экспериментальные исследования показывают, что внутренние течения в капле характеризуются явной трехмерностью, а одним из определяющих параметров является угол смачивания.

Также изучение особенностей процессов, сопровождающих вынужденные колебания капли жидкости, целесообразно проводить, во-первых, на каплях малого объема - до 10 мкл, а, во-вторых, вибрации основания не должны приводить к разбрызгиванию капли жидкости, то есть режим вибрационного воздействия не должен превышать порог резонанса [6 - 8]. Исследованию внутренних течений в малом объеме жидкости - капле и их влиянии на механизм возникновения и порога возбуждения волн посвящена данная работа.

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

ПОСТАНОВКА ЗАДАЧИ

Рассмотрим задачу о движении капли жидкости, вызванном вертикальными вибрациями жесткой подложки. Пусть область П £ М3(рис. 1, а) заполнена двухфазной средой и П = П1 U Г U П2. Для определенности будем считать, что в подобласти П1 находится газ, а в П2 - жидкость; Г является границей раздела фаз. В начальный момент времени t = 0 подложка П0 начинает вибрировать под действием акустических возмущений. Расположенная на подложке П0 капля жидкости П2 также начинает совершать вертикальные колебания.

Экспериментальная установка подробно описана в работе [7], необходимо отметить то, что поверхность подложки тщательно очищалась с последующим нанесением гидрофобного покрытия. Деионизированная капля воды, объемом 5 мкл с помощью пипетки размещалась в центре подложки. При этом, измеренный контактный угол водной капли составлял 115°. В данной работе исследуются только колебания капли под действием низкочастотных вибраций основания — мода 2 (частота 85 Гц).

я*

<->

а) б)

Рис. 1. Равновесная конфигурация капли: а) — расчетная область; б) — геометрические характеристики

МАТЕМАТИЧЕСКИЕ МОДЕЛИ И АЛГОРИТМЫ

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

дрЩ + у-

У •и1 = 0,

дЬ

где индекс 1={ 1,2} - номер среды, р - плотность, и - вектор скорости, р - давление, = ^(Уи^ + и!) - тензор вязких напряжений, ^ - динамическая вязкость, д - вектор ускорения свободного падения.

На границе раздела сред Г1 выполняются условия динамического равновесия

иг = и2,

где е - единичный нормальный вектор, внешний по отношению к П2, К - кривизна поверхности Г1, $ - коэффициент поверхностного натяжения, 0О - угол смачивания.

Геометрические характеристики поверхности раздела фаз Г1 однозначно определяются

как

Уа Уа

,К = У

где а - скалярная функция

1УаУ \Уа\ _( 1,хеп2,

0 1,хеп21 1 а(1,х) = {1,хеп1иг1у

имеющая в данном случае смысл объемной концентрации жидкости.

Тогда функции плотности давления и вязкости определим, соответственно, как р = р2 + (р1 — p2)a, Р = Р2 + Pi — Р2а и р = р2 + (P-i — P2)a и перепишем систему (1), (2) с учетом (4), (3) следующим образом

^ + V

V -U = 0.

Система (5), (6) дополняется уравнением переноса для а

%=U-Va.

Силу поверхностного натяжения, действующую в узком переходном слое, который в пределе представляется бесконечно тонким и совпадает с Г1 можно свести к объемной

G = —<rKVaS(a)(1 — cos0o).

Численное решение системы

^+V -ии = —Vp + V  pVu + G.

Запишем аппроксимацию уравнений (6), (8) в полудискретном виде:

I¡feNMSru;>+1=0

1п,П+1-А-,,П -¡л,™-1

~2FM+UM Gn + 'LfeNM If  и}1?1 — 'L/eNm P Sf

здесь uJ - скорость в центре ячейки с номером М на n-ом временном шаге; VM - объем ячейки; Ff - поток через грань с номером /, NM - множество номеров соседних ячеек; S; - вектор внешней нормали к грани с номером /, по модулю равный площади этой грани.

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

Уравнение (7) решается на каждом шаге по времени после отыскания решения системы (5), (6). Для этого (7) удобно переписать в консервативном виде

j^ + V -

где ur = u2 — u1 - скорость относительного движения фаз на границе раздела Г1. В исходной модели Г1 является бесконечно-тонкой и ur = 0, однако из-за численной диффузии при вычислениях неизбежно появляется переходная область. Следуя методу контрольного объема, запишем аппроксимацию

~Тм + + Zf£<M If

„i, Л qurS Va

где Fl(a) = camax-q,q =—

W a Q S IVal

Скалярный коэффициент c; является параметром управления искусственным сжатием решения в области разрыва для компенсации эффекта численной диффузии.

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

ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ

В рассматриваемой задаче [7] капля деионизированной воды объемом 5 мкл при помощи микропипетки располагается на оси симметрии недеформируемой подложки. Геометрические характеристики покоящейся на неподвижной подложке капли (рис. 1, б) снимались до инициации вибраций подложки. Так, угол смачивания составлял в0 = (115 ± 1)0 [7], а диаметр контакта и высота капли составляли Ок = 2,02 мм и к = 1,52 мм, соответственно.

Поскольку численное моделирование производится методом УОБ, получение численной равновесной конфигурации капли является самостоятельной задачей. В рамках получения равновесной конфигурации рассматривалось две основных постановки: уравновешивание капли под действием силы тяжести и поверхностного натяжения с полусферы (рис. 2, а) и с цилиндра (рис. 2, б). В ходе установления капли измерялись как диаметр контакта Ок, высота капли к, так и величина контактного угла 0О.

с1

<->

а) б)

Рис. 2. Начальная форма жидкости в задаче установления: а) - полусфера; б) - цилиндр

Выявлено, что использование в качестве начальной формы жидкости полусферы позволяет получить равновесную конфигурацию за меньшее число шагов по времени, однако приводит к уменьшению высоты равновесной капли к на 5 - 8 %. При этом, определение местоположения жидкости в начальный момент времени как цилиндрического столба приводит не только к существенному увеличению времени расчета и числа временных шагов (более чем в 1,4 раза), но и к дополнительным осцилляциям поля скорости капли жидкости в процессе уравновешивания. Геометрические параметры полученных равновесных конфигураций капли приведены в таблице.

Таблица

Геометрические характеристики равновесной конфигурации капли

Параметры Начальная форма жидкости Эксперимент [7]

полусфера цилиндр

Диаметр контакта Вк, мм 2,8 2,4 2,02

Высота капли к, мм 1,47 1,508 1,51

Угол смачивания 0О, ° 112,4 114,3 115 ± 1

Рассмотрим методы учета контактного угла в линии трехфазного контакта газ-жидкость-подложка. Угол смачивания в определяется не только физико-химическими свойствами жидкости, но и свойствами материала подложки и не является константой даже при неизменности жидкости и подложки. Поскольку при колебаниях подложки динамическая составляющая контактного угла отлична от нуля, возможно оценить ее влияние на деформацию профиля разреза колеблющейся капли, и, как следствие, на локальную высоту капли. На рис. 3 представлен график изменения высоты профиля разреза капли от времени, полученный в результате моделирования вибрационного процесса с учетом динамического изменения угла смачивания (кривая 1) и с применением только модели статического краевого угла (кривая 2).

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

мп

0,2 0,1 О -0,1 -ол

-0,3

- 1  3

 »   ш я

г х 1, У \ ¿у

/

^_> • -

10

15

£ 1,1 с

Рис. 3. Графики изменения высоты капли, полученные с учетом: 1 — динамического угла смачивания (гистерезис не превышает 1°); 2 — статического контактного угла; 3 — эксперимент [7]

При вибрации подложки с частотой 85 Гц (рис. 4) зафиксирована форма колебаний, представляющая одну кольцевую линию узлов на свободной поверхности, соответствующая второй моде колебаний [7].

Рис. 4. Формы колебаний капли для моды 2: а) - численная; б) - экспериментальная [7]

Рассмотрим график изменения высоты капли при вибрации подложки (рис. 5). Из него видно, что зависимость изменения высоты капли от времени, полученная в результате моделирования, качественно близка к экспериментальной, а зафиксированное расхождение частот колебаний жидкости связано с тем, что экспериментальная кривая (рис. 5, кривая 1) отклонена от классической синусоиды, соответствующей частоте 85 ГЦ.

Н мм

0,2 0,1 о -0,1 -0,2 -0,3

О 5 10 15 20 ь мс

Рис. 5. Высота осцилляций капли для моды 2: 1 — эксперимент [7]; 2 — численно полученные колебания капли

Для изучения особенностей внутренних течений получены мгновенные картины векторов скорости в капле для различных моментов времени (рис. 6).

2.5$ 2.1Ы*

Рис. 6. Изохоры в поперечном сечении капли

ч  1 2

Я    ф и

1 1 и

I     щ

и 1 1

Как видно из рис. 6, течения в капле значительным образом определяются направлением движения подложки. Так в начальные моменты подъема подложки наблюдается перестроение векторов скорости и смена направления движения течений в жидкости. Четко визуализируются парные симметричные вихри. При подъеме подложки наблюдается объединение вихрей в укрупненный симметричный парный вихрь, с четко выраженными двумя линиями стекания. Таким образом, внутренние течения жидкости с опозданием реагируют на смену направления движения подложки и начинают изменение направления движения с нижних слоев жидкости. Отметим также, что в моменты времени, соответствующие началу следующего периода и максимальному положению подложки, в капле наблюдаются безвихревые течения.

ЗАКЛЮЧЕНИЕ

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

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

Работа выполнена при частичной финансовой поддержке РФФИ (грант 16-38-00127 мола, грант 16-37-00060-мол_а) и комплексной программы фундаментальных исследований УрО РАН, проект №15-7-1-11.

СПИСОК ЛИТЕРАТУРЫ

1. Коренченко А. Е., Илимбаева А. Ж., Бескачко В. П. Численное исследование свободный колебаний лежащей капли // Вестник Южно-Уральского государственного университета. Серия: Математика. Механика. Физика. 2011. № 10(227). С. 72-76.

2. James A. J., Vukasinovic B., Smith M. K., Glezer A. Vibration-induced drop atomization and bursting // Journal of Fluid Mechanics, 2003, vol. 476, pp. 1-28. DOI: https://doi.org/10.1017/S0022112002002835

3. Vukasinovic B., Smith M. K., Glezer A. Dynamics of a sessile drop in forced vibration // Journal of Fluid Mechanics, 2007, vol. 587, pp. 395-423. DOI: https://doi.org/10.1017/S0022112007007379

4. Kim H. Y. Drop fall-off from the vibrating ceiling // Physics of fluids, 2004, vol. 16, no. 2, pp. 474-477. http://dx.doi.org/10.1063/L1637352

5. Noblin X., Buguin A., Brochard-Wyart F. Vibrated sessile drops: Transition between pinned and mobile contact line oscillations // European Physical Journal E, 2004, vol. 14, iss. 4, pp. 395-404. DOI: 10.1140/epje/i2004-10021-5

6. Shin Y. S., Lim H. Ch. Shape oscillation and detachment conditions for a droplet on a vibrating flat surface // European Physical Journal E, 2014, vol. 37/8, no. 74, pp. 1-10. DOI: 10.1140/epje/i2014-14074-5

7. Kim H., Lim H.Ch. Mode Pattern of Internal Flow in a Water Droplet on a Vibrating Hydrophobic Surface // The Journal of Physical Chemistry B, 2015, vol. 119, iss. 22, pp. 6740-6746. DOI: 10.1021/acs.jpcb.5b02975

8. Park Ch. S., Kim H., Lim H. Ch. Study of internal flow and evaporation characteristics inside a water droplet on a vertically vibrating hydrophobic surface // Experimental Thermal and Fluid Science, 2016, vol. 78, pp. 112-123. http: //dx .doi.org/10.1016/j. expthermflusci. 2016.05.018

9. Wilkes E. D., Basaran O. A. Drop Ejection from an Oscillating Rod // Journal of Colloid and Interface Science, 2001, vol. 242, iss. 1, pp. 180-201. doi:10.1006/jcis.2001.7729

10. Issa R. I. Solution of implicitly discretised fluid flow equations by operator-splitting // Journal of Computational Physics, 1986, vol. 62, iss. 1, pp. 40-65. doi:10.1016/0021-9991(86)90099-9

SIMULATION OF OSCILATIONS LIQUID DROPLET LYING ON A VIBRATING UNDEFORMABLE BASE

Kuzmin I. M., Sarmakeeva A. S., Chernova A. A.

Institute of Mechanics, Ural Branch of the Russian Academy of Sciences, Izhevsk, Russia

SUMMARY. The behavior of a drop of viscous Newtonian fluid, which is lies on a solid un-deformable rigid base, produces vertical vibration of low frequency (85 Hz) is studied. Numerical simulation of dynamics of a vibrating fluid by Volume of Fluid method (VOF) is produced. Therefore, was discussed how to obtain the equilibrium configuration of drop. It is found that for stationary droplet form, close to the experimental use of a cylindrical column to describe initial situation of liquid is appropriate. Effect of accounting model the contact angle at the triple point of the base liquid gas is investigated. The need to address the dynamic changes in the contact angle is shown. Waveforms Profile drop isotachs in the plane of symmetry, and changes depending on the time of drop height for the second mode of forced oscillations of a drop corresponding to the base of the vibration with a frequency of 85 Hz, also obtained. One ring line of nodes on free surface, corresponding to the second oscillation mode is fixed. Dependence of topology internal flows in a drop from the base movement direction identified and described in the example isotachs. The dependence of the change of height drops resulting from modeling the known experimental data corresponds to. Comparison of numerical results with experimental data is made and allows us to conclude that the method used adequately describes the processes of energy transfer from the base to drop.

KEYWORDS: liquid droplet, Faraday wave, free surface, Volume of Fluid method, contact angle. REFERENCE

1. Korenchenko A. E., Ilimbaeva A. Zh., Beskachko V. P. Chislennoe issledovanie svobodnyy kolebaniy lezhashchey kapli [Numerical study of free vibrations of a sessile drop]. Vestnik YuUrGU. Matematika. Mekhanika. Fizika, 2011, no. 10(227), pp. 72-76.

2. James A. J., Vukasinovic B., Smith M. K., Glezer A. Vibration-induced drop atomization and bursting. Journal of Fluid Mechanics, 2003, vol. 476, pp. 1-28. DOI: https://doi.org/10.1017/S0022112002002835

3. Vukasinovic B., Smith M. K., Glezer A. Dynamics of a sessile drop in forced vibration. Journal of Fluid Mechanics, 2007, vol. 587, pp. 395-423. DOI: https://doi.org/10.1017/S0022112007007379

4. Kim H. Y. Drop fall-off from the vibrating ceiling. Physics of fluids, 2004, vol. 16, no. 2, pp. 474-477. http://dx.doi.org/10.1063/L1637352

5. Noblin X., Buguin A., Brochard-Wyart F. Vibrated sessile drops: Transition between pinned and mobile contact line oscillations. European Physical Journal E, 2004, vol. 14, iss. 4, pp. 395-404. DOI: 10.1140/epje/i2004-10021-5

6. Shin Y. S., Lim H. Ch. Shape oscillation and detachment conditions for a droplet on a vibrating flat surface. European Physical Journal E, 2014, vol. 37/8, no. 74, pp. 1-10. DOI: 10.1140/epje/i2014-14074-5

7. Kim H., Lim H.Ch. Mode Pattern of Internal Flow in a Water Droplet on a Vibrating Hydrophobic Surface. The Journal of Physical Chemistry B, 2015, vol. 119, iss. 22, pp. 6740-6746. DOI: 10.1021/acs.jpcb.5b02975

8. Park Ch. S., Kim H., Lim H. Ch. Study of internal flow and evaporation characteristics inside a water droplet on a vertically vibrating hydrophobic surface. Experimental Thermal and Fluid Science, 2016, vol. 78, pp. 112-123. http: //dx. doi.org/10.1016/j. expthermflusci .2016.05.018

9. Wilkes E. D., Basaran O. A. Drop Ejection from an Oscillating Rod. Journal of Colloid and Interface Science, 2001, vol. 242, iss. 1, pp. 180-201. doi:10.1006/jcis.2001.7729

10. Issa R. I. Solution of implicitly discretised fluid flow equations by operator-splitting. Journal of Computational Physics, 1986, vol. 62, iss. 1, pp. 40-65. doi:10.1016/0021-9991(86)90099-9

Кузьмин Игорь Михайлович, младший научный сотрудник ИМ УрО РАН

Сармакеева Анастасия Семеновна, младший научный сотрудник ИМ УрО РАН

Чернова Алена Алексеевна, кандидат технических наук, старший научный сотрудник ИМ УрО РАН, e-mail: alicaaa@,gmail. com