Эволюция малых деформаций парового пузырька при его сильном сжатии в бензоле

ГИДРОДИНАМИКА, ТЕПЛО- И МАССООБМЕННЫЕ ПРОЦЕССЫ, ЭНЕРГЕТИКА

УДК 532.528.2

А. А. Аганин, Д. Ю. Топорков

ЭВОЛЮЦИЯ МАЛЫХ ДЕФОРМАЦИЙ ПАРОВОГО ПУЗЫРЬКА ПРИ ЕГО СИЛЬНОМ СЖАТИИ В БЕНЗОЛЕ

Ключевые слова: кавитационный пузырек, искажение сферичности пузырька, коллапс пузырька.

Исследовано сильное сжатие кавитационного пузырька в бензоле. Использована математическая модель, в которой движение пара и жидкости разлагается на сферическую составляющую и ее малое несферическое возмущение. Сферическая составляющая описывается уравнениями газовой динамики с учетом нестационарной теплопроводности, неравновесного испарения-конденсации на межфазной поверхности, с использованием широкодиапазонных реалистичных уравнений состояния. При описании несферической составляющей приближенно учитываются эффекты вязкости жидкости и поверхностного натяжения, влияние содержимого пузырька. Показано формирование ударной волны в пузырьке в финале его сжатии. Установлено, что рост несферичности пузырька при его сжатии не превышает 180 раз при начальном радиусе пузырька более 100 мкм.

Keywords: cavitation bubble, sphericity bubble distortion, bubble collapse.

Strong compression of a cavitation bubble in benzol has been investigated. To describe the motion of the liquid and vapor, a mathematical model with decomposing into spherical component and a small non-spherical perturbation is used. The spherical component is described by the gas dynamics equations with the transient heat conduction, nonequilibrium evaporation-condensation at the interface and the realistic wide-range equations of state. To describe the motion of non-spherical component, the effects of viscosity, liquid surface tension and the bubble content are approximately taken into account. A shock wave has been shown to form in the bubble during its final compression stage. It has found that during compression of the bubble the growth of its nonsphericity does not exceed 180 times if the initial radius of the bubble is more than 100 /мп.

Введение

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

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

дейтерированного ацетона [1]. В последнее время опубликован ряд исследований по синтезу

наноразмерных алмазов, где главным действующим

фактором является энергия, выделяемая в результате кавитационных процессов в жидкости. Следуя гипотезе академика Э.М. Галимова [2] о возможности кавитационного синтеза алмазов в природных условиях при схлопывании кавитационных пузырьков в кимберлитовой породе, группа ученых под его руководством произвела синтез наноалмазов кавитационной диссоциацией бензола [3]. С точки зрения приложений успехи в направлении кавитационного синтеза алмазов позволили бы в дальнейшем получать особо чистые, изометрические алмазные порошки заданной узкой фракции, и, как результат, исключить, например, трудоемкие операции дробления и классификации микропорошков. Вполне определенный размер наноалмазов требуется, например, при золочении ультрадисперсными алмазами спиц для компрессионно-дистракционного остеосинтеза [4,5]. В работе [6] проведено математическое моделирование процесса схлопывания

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

В настоящей работе изучается эволюция малых отклонений от сферичности кавитационного пузырька при его сильном сжатии в бензоле в условиях, близких к экспериментальным по кавитационному синтезу наноалмазов [3]. Применяется математическая модель, аналогичная используемой в работе [7], в которой движение пара и жидкости расщепляется на сферическую и несферическую составляющие. Для описания сферической составляющей используется полная гидродинамическая модель как для пара, так и для жидкости. Учитываются нестационарная теплопроводность в паре и жидкости,

неравновесные испарение-конденсация на межфазной поверхности, применяются

реалистичные уравнения состояния, построенные по экспериментальным данным. При описании

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

содержимого пузырька в предположении однородности распределения плотности.

1. Физическая постановка задачи

Рассматривается эволюция малых искажений сферичности кавитационного пузырька при его сильном сжатии в бензоле в условиях близких к экспериментальным по кавитационному синтезу наноалмазов [3,6]. Начальная температура пара и жидкости равна 20°С. В начальный момент времени t = 0, за который принимается начало сжатия, кавитационный пузырек имеет радиус R0 (следуя работе [6] R0 изменяется в диапазоне от 10 до 1000 мкм). В начале сжатия пар в пузырьке и окружающая жидкость покоятся, пар находится в состоянии насыщения (с давлением около 0,11 бар). Поверхность пузырька имеет малое осесимметричное отклонение от сферической формы, настолько малое, что, несмотря на свой значительный рост, оно остается малым вплоть до конца сжатия. Следуя работе [6], воздействие на пузырек бегущей ударной волны моделируется мгновенным увеличением давления жидкости при t = 0 до неизменного в последующем уровня рм, которое полагается равным 400 бар (в [6] на предмет организации условий образования алмазов исследуется диапазон 50 < рм < 400 бар). В результате огромного перепада давления между жидкостью и паром, пузырек стремительно и очень сильно сжимается.

2. Математическая постановка задачи

В случае осесимметричных искажений сферичности пузырька уравнение его поверхности в сферической системе координат г, 9, ф представляется в виде ряда по полиномам Лежандра

г = R(t) + Х Zi (t )P (cos 9),

i=2

где t - время, ai - амплитуда отклонения поверхности пузырька от сферической формы радиуса R в виде полинома Лежандра P,(cos9) степени i. Искажения 6, = a/R предполагаются малыми при всех I. | 6 | << L С учетом этого движение жидкости и пара представляется в виде суперпозиции сферического (радиального) движения и его несферического возмущения.

Для описания радиальной составляющей движения пара и жидкости используется система уравнений [8,9]

^ (рг2) + ^ (pwr2 ) = 0, dt дг

|- (pwr 2) + ^ (pw 2r2 + pr2) = 2pr, dt dr

— (per2) + — [wr2(ре + p)] =—\ г2К— |.

dt дг дг V дг)

Здесь p - плотность, w - радиальная компонента вектора скорости жидкости w, р - давление, е -удельная полная энергия, T - температура, к -коэффициент теплопроводности. Для описания свойств жидкости и пара используются уравнения состояния в виде Ми-Грюнайзена, представляющие

собой сумму потенциальной (холодной), тепловой и химической компонент давления и внутренней энергии [9]

р = р

Холодная составляющая давления и энергии описывается потенциалом типа Борна-Майера [9,10]:

p(p) = A

(P ^ -p+1 ( (P ^ -p> ( P ^

_P exp b 1 - P - K P

vP° J Ч 4P° J J 4P° J

\'+і

e( p) = •

A

PPob

exp

1 -

( Yp^

_P_

Po

K

'Po

(4

vP0 У

+ eo

Тепловая составляющая давления и энергии принимается в виде:

р(Т) = рГ(р)с„, е(Т) = ^Т,

где Г(р) - функция Грюнайзена, су - теплоемкость при неизменном объеме, которая полагается постоянной.

Молекулярная газовая и жидкая фазы описываются единым по давлению уравнением состояния и различаются теплоемкостью и функцией Грюнайзена [10]:

Rg

Г(р) = -g-cV

+а(1) exp

а(0) +(l - а(0) jexp

( \ P р°) ЧР У

,

+ а(2) exp

( \ P

v p(°) ,

ЧК

( \ P

P(2),

,

PS

с

где Rg - удельная газовая постоянная.

Зависимость давления насыщения аппроксимируется в соответствии экспериментальными данными [10].

( т ^

Ps (T) = Р* exp - .

V ' ' 2 У

Для аппроксимации теплоты

парообразования l будет использоваться несколько модифицированное по сравнению с [11] уравнение, удовлетворительно согласующееся с

экспериментальными данными.

l (T) = 4.45 • 105

-TL -1

Граничные условия в центре пузырька (г = 0), вдали от него (г = гт гт >> Я) и на межфазной границе (г = Я(ї)) имеют вид [8,9] дТ

г = 0 : ш = 0, — = 0; г = гш : р = рш, Т = То; дг

г = Я : Я = +-^~ = ^

Pi = p

'dT

Pi

4|a,w, 2ст

~R ’

Pg

R

dT

IF )rKg IsT1 = jl (pg >’T = Tg ’

где ц - коэффициент вязкости жидкости, ст -коэффициент поверхностного натяжения, 7 -

интенсивность фазовых превращений (испарения или конденсации), отнесенная к единице времени и единице поверхности. Нижние индексы I и д относятся соответственно к параметрам жидкости и газа (пара). Интенсивность фазовых превращений 7 определяется равенствами [8,9]

J =

X =

Ps(Т) xpg

л/f ЛІТ

(

Q =

exp (q2 )-Q/tc 1 ^-7=j exP (~x2 )dx

V2p,

g

Здесь а' - коэффициент аккомодации.

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

Процессы диссоциации и ионизации в паре не учитываются в силу чрезвычайной кратковременности их протекания.

Несферическая составляющая движения пара и жидкости определяется без учета влияния сжимаемости жидкости и неоднородности распределения плотности пара. Для описания отклонения а, применяется модель [12], в которой учет влияния плотности пара производится согласно [13]

(i+q¡ )а, +

3 — + 2(i +1)(/ + 2)-^-R R

2 4v; (i2 - 1)R R

2 + ~---------(/ -1)(1 - Qi)—

i

R

R3 .............R

VQi

i +1

R

• Qi

Ra

ai + R Pi

= 0,

— f Q^dr, pi = f

i+1R ґ R

2,

я и + 1)рд _2 ст

Я = —.------, , =----------—3—.

,Ро Ро^3

Здесь рд - величина, выбирающаяся из стремления

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

что плотность жидкости равна своему начальному значению р0. Функция О(г, Ї) вводится следующим образом [12]

V X «7 = V X Р (СОв0 ^ег,

где ег - единичный вектор вдоль радиальной координаты г. Функция О,(г, Ґ) определяется из уравнения

дО . Я я 2 АГ О

dt

Г i (i + 1)Qi

d2Q ^ dr 2

: 0

и граничных условий

О, Л) = т+- ((/ + 2)а - (/ -1) Л а, + (2/ + 1)Л -1а,

I + 1 I Гл

О (<», t) = о.

В качестве рд в настоящей работе

выбирается граничное значение плотности пара на границе пузырька, полученное из решения системы (1)-(4), в отличие от, например, [13], где выбиралась среднее по пузырьку значение плотности пара. Проведенное в [14] сравнение результатов применения аналогичной математической

постановки с результатами прямого численного моделирования для случая сильного сжатия кавитационного пузырька в другой углеводородной жидкости (ацетоне) [7] показало, что

математическая постановка настоящей работы дает удовлетворительную оценку максимального роста искажения сферичности пузырька для наиболее важного диапазона длинноволновых возмущений сферичности (2 < У < 10).

3. Входные данные

Константы для уравнения состояния, Г(р), р3 и I [10,11]: Лд = 106,6 Дж/(кг К), су = 960 Дж/(кг К) для пара и су = 1194 Дж/(кг К) для жидкости, невозмущенная плотность жидкости р0 = 879 кг/м3, е(Щ = 1,15105 Дж/кг, ео = 0,392-106 Дж/кг, А = 7107 Па, Ь = 37,35, К = 0,38 109 Па, 4 = 1,2, в = 0,3333,

a(0) = 14,57919,

a(1) = 6,04958,

a(2) = 24,62614,

р(0)/р0 = 1,29314, р(1)/р0 = 1,06761 , р(2)/р0 = 1,3439, а(0) = 1,75, а(1) = -2,4, а(2) = -3,1, р* = 2,68-109 Па, Тсг = 562,6 К, Т = 3450 К, Т2 = 15 К. р0 = 879 кг/м3.

Зависимость параметров ц, ст, к от температуры не учитывается, т. е. их величины принимаются неизменными и равными начальным значениям при Т = 293,15 К: ц = 6,52-10-4 Па-с,

ст = 0,02888

Н/м,

к, = 0,009315 кг-м/(с3-К),

кд = 0,1465 кг-м/(с-К). Коэффициент аккомодации а' выбирается равным 0,8 [15]. Скорость звука в жидкости С0 = 1324 м/с.

При t = 0 полагается:

R0 = 10, 100, 1000 мкм, R0 = 0, е/ = е/,0, éi,0 = 0 , Qi = 0.

Давление в жидкости рм = 400 бар.

Решение отыскивается численно методом С.К. Годунова [16] с применением подвижной неравномерной сетки - 1000 ячеек в жидкости и 300 в паре (размер ячеек увеличивается по

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

4. Результаты

Радиальное движение. Из-за большой разности давлений в окружающей жидкости и в покоящемся в начальный момент времени паре пузырек начинает очень сильно сжиматься. Ряд важных особенностей радиального движения пара в пузырьке, окружающей жидкости и поверхности пузырька при его сильном сжатии иллюстрируют рис. 1, 2. Для трех случаев R0 = 10, 100, 1000 мкм на рис. 1 показано изменение радиуса пузырька R, отнесенного к R0, и скорости его изменения R в ходе всего сжатия в зависимости от параметра t c0/R0. Радиус пузырька в этих случаях уменьшается соответственно до Rmin = 0,096, 3,61, 51,1 мкм. Изменение R(t)/R0 для всех трех случаев практически не отличается. Можно отметить, что время сжатия удовлетворительно согласуется с рэлеевским временем полного схлопывания пустой полости [17] tR = 0.915R0^/р;/рш , что было ранее

проиллюстрировано в [6] для более простой математической постановки.

На довольно продолжительном начальном

участке величина скорости сжатия | R | относительно невелика (| R | < 300 м/с при т < 3,83, где R/R0 > 0,6).

Однако в конце сжатия | R | сначала быстро растет (для R0 = 10, 100, 1000 мкм до 13,2, 3,8, 2,8 км/с

соответственно), после чего еще быстрее уменьшается до момента, когда скорость R скачком меняет отрицательное значение на положительное. На рис. 1 эволюция радиуса показана до момента

экстремального сжатия содержимого пузырька tc (R = Rc = 0,1, 3,76, 53,3 мкм при R0 = 10, 100, 1000 мкм соответственно), который наступает несколько раньше, чем момент достижения Rmin. При t = tc

скорость сжатия для R0 = 10, 100, 1000 мкм

оказывается равной R = -11, -3,1, -2,3 км/с

соответственно.

Рис. 1 - Изменение радиуса пузырька Я при Я0 = 1000 (1), 100 (2), 10 (3) мкм и скорости Л при тех же Я0 (соответственно 4-6) при его коллапсе. Точками отмечены значения в момент экстремального сжатия пара в пузырьке Ъ

Финальная стадия сжатия пузырька характеризуется большими радиальными скоростями. Из-за этого в пузырьке образуется сходящаяся сферическая ударная волна. Ее интенсивность по мере схождения быстро нарастает. На рис. 2 показано ее формирование,

распространение и фокусировка в центре пузырька для R0 = 100 мкм. Можно отметить, что первой появляется простая волна сжатия небольшой амплитуды, которая не успевает

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

10* і о8 ю7 106 і о5 ю4 103 ю2 ю1 10°

1<Ґ

0 30 60г’мкм90

Рис. 2 - Радиальные распределения давления в шесть последовательных моментов времени финальной стадии сжатия

Следует отметить, что при измельчении ячеек сетки в паре в 10 раз представленная на рис.2 как эволюция межфазной границы, так и числовые характеристики формирования и распространения ударной волны изменяются незначительно, за исключением центра пузырька (r = 0), где решение в рамках принятой постановки задачи имеет особенность.

Сжатие кавитационного пузырька в бензоле продолжается и после момента фокусировки ударной волны в центре пузырька (R = Rc), пока расходящаяся ударная волна, возникающая в центре пузырька в момент фокусировки там сходящейся волны, не дойдет до поверхности пузырька

(R = Rmin) и не развернет его движение. Разница между радиусами Rc и Rmin незначительна (не более 5% для рассматриваемых R0). Поскольку наибольшее внимание вызывает экстремальное сжатие пара при R = Rc, то на отрезке Rmin < R < Rc не представляют большого интереса как изменение термодинамических параметров, так и эволюция несферичности пузырька.

Рост несферичности пузырька. Сферическая форма пузырька при его сжатии, как известно [18,19], неустойчива и малые отклонения формы пузырька от сферической в ходе сжатия нарастают (более точно, испытывают колебания с возрастающей амплитудой). За меру роста несферичности удобно принять величину |е,- / б,-,0|. Типичное изменение величины |е,- / б,-,0| в ходе сжатия для R0 = 10, 100, 1000 мкм иллюстрирует рис. 3. Видно, что чем более сильное сжатие испытывает пузырек, т.е. чем больше отношение R0/Rc, тем сильнее вырастает искажение к концу коллапса. В начале сжатия искажение во всех трех случаях изменяется подобным образом. При этом можно увидеть, что в случае при R0 = 10 мкм амплитуда колебаний искажения в процессе колебаний (при R/R0 > 0,3) оказывается чуть меньше, чем при других R0. Это обстоятельство объясняется уже заметным влиянием вязкости при R0 = 10 мкм, тогда как для R0 = 100, 1000 мкм влияние вязкости на эволюцию искажения при i = 7 еще незначительно. Как известно, чем меньше пузырек, тем большее влияние на его динамику оказывает вязкость.

Рис. 3 - Эволюция относительной амплитуды

возмущения |в,- / S/,o| при сжатии пузырька для / = 7 при R0 = 10 (3), 100 (2), 1000 (1) мкм. Точкой отмечен момент tc

Рост амплитуды несферичности пузырька при сжатии для фиксированного i будет определяться величиной |є/>с / є/о|, где |є/>с| - максимальное значение величины |s/|, которое (за исключение сверхвысокочастотных искажений сферичности) достигается либо в момент tc, либо незадолго до него. На рис. 4 проиллюстрирована зависимость величины |s/,c / s/,o| от номера гармоники i при R0 = 100 мкм, а также огибающие таких зависимостей при R0 = 10, 100, 1000 мкм. Т.к. при сжатии |s/| испытывает колебания, зависимость |є/,с / є,0| от номера / из-за

влияния фазы колебаний носит колебательный характер (рис. 4). Огибающие зависимостей |eic / si,0| от номера i, следуя [7], используются в качестве приближений зависимостей от номера i максимальных значений степени роста амплитуды несферичности пузырька |eic / si,0|, реализующихся при вариации начальной скорости é¡,0, которая в

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

Видно, что с увеличением номера i от i = 2 в пределах определенного интервала длинноволновых возмущений сферичности, где влияние вязкости незначительно, величина |s¡,c / s¡,0| возрастает. Начиная с некоторого номера i, частота и, как следствие, скорость колебаний возмущения сферичности возрастают настолько, что из-за достаточно возросшего влияния вязкости величина |s¡,c / s¡,0| начинает быстро уменьшаться. Как отмечалось выше, чем меньше пузырек, тем сильнее влияние вязкости на его форму. Поэтому номер i, с которого начинает проявляться влияние вязкости на эволюцию формы пузырька, тем меньше, чем меньше R0. Как следствие, максимальное значение |s¡c / s¡,0| для R0 = 10, 100, 1000 мкм достигается при i = 7, 13, 130 и оказывается равным 950, 180, 160 соответственно. Для R0 = 10 и 100 мкм огибающие показывают монотонный рост (значений локальных максимумов зависимости величины |s¡,c / s¡,0| от i) при i < 7 и 13 соответственно, после чего монотонное уменьшение. Для R0 = 1000 мкм огибающая при i < 130 показывает немонотонный рост, после чего монотонное уменьшение. Для R0 = 100 мкм также заметно некоторое колебание огибающей на интервале изменения i, где нарушается монотонность роста огибающей для R0 = 1000 мкм.

Рис. 4 - Огибающие зависимостей величины |б/,с / е#,о|, характеризующей рост несферичности пузырька при его сжатии, от номера / для Я0 = 10

Если ограничиться наиболее опасными (с точки зрения устойчивости сферичности пузырька при его сжатии) длинноволновыми возмущениями

сферической формы пузырька, а именно диапазоном изменения номера i от 2 до 10, то для R0 = 10, 100, 1000 мкм, судя по огибающим, величина |е,,с / s¡,0| меняется в интервале от 700, 115 и 75 (при i = 2) до 950 (при i = 7), 180 и 115 (при i = 10) соответственно. Как уже отмечалось выше, косвенное подтверждение удовлетворительности полученных с использованием системы (1)-(5) результатов для огибающей зависимости величины |е,,с / s¡,0| от номера i на интервале 2 < i < 10 можно найти в [14], где рассмотрен случай сжатия слабонесферического пузырька в ацетоне.

Заключение

Выполнено исследование роста возмущений сферичности отдельного кавитационного пузырька в бензоле при его сжатии в условиях экспериментов [3] (температура жидкости Т0 = 20°С, давление жидкости рм = 400 бар; в начале сжатия пар в пузырьке находится в состоянии насыщения при температуре 20°С). Начальный радиус пузырька выбирается равным 10-1000 мкм. Тогда радиус Rmin, до которого сжимается пузырек, меняется в диапазоне 0,1-51 мкм соответственно. В финальной высокоскоростной стадии сжатия в пузырьке формируется сферическая ударная волна, сходящаяся к его центру. Ее интенсивность по мере схождения быстро нарастает.

Установлено, что амплитуда осесимметричных возмущений сферичности пузырька в виде полиномов Лежандра P,(cos9) при R0 = 10, 100, 1000 мкм в силу разной степени сжатия пузырька может возрасти к концу коллапса не более чем в 950, 180, 160 раз соответственно. Максимум роста достигается для указанных R0 при i = 7, 13, 130. Меньший рост амплитуды коротковолновых возмущений для меньших R0 объясняется тем, что демпфирующее влияние вязкости тем сильнее, чем меньше пузырек.

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

радиусом R0 > 100 мкм к концу сжатия, независимо от /, она возрастет не более чем 0.02 R, т.е. пузырек и в конце сжатия останется довольно близким к сферическому.

Работа выполнена в рамках программ РАН и при поддержке РФФИ (проекты 12-01-00341-a, 12-01-З1410-мол_а).

Литература

1. R.P. Taleyarkhan, C.D. West, J.S. Cho, R.T. Lahey

2. E.M. Galimov, Nature, 243, 389-391

3. Э.М. Галимов, А.М. Кудин, В.Н. Скоробогатский и др., ДАН, 395, 2, 187-191

4. Р.А. Кайдриков, Я В. Ившин, А.И. Зильберг, Вестник Казан. технол. ун-та, 22, 137-140

5. Я В. Ившин, Р.А. Кайдриков, А.И. Зильберг, Вестник Казан. технол. ун-та, 15, 4б-48

6. А.Ю. Днестровский, С.А. Воропаев, Е.А. Пономарева, ДАН, 436, 5, б11-б14

7. А.А. Аганин, М.А. Ильгамов, Р.И. Нигматулин, Д.Ю. Топорков, Известия РАН, МЖГ, 1, 57-б9 (2010).

8. R.I. Nigmatulin, I.Sh. Akhatov, A.S. Topolnikov, R.Kh. Bolotnova, N.K. Vakhitova, R.T. Lahey (Jr.), R. Taleyarkhan, Phys. Fluids, 17, 10, 10710б (2005).

9. Р.И. Нигматулин, Динамика многофазных сред, т. 1 и 2. Наука, Москва, 1987.

10. Р.Х. Болотнова, В сб. Труды института механики Уфимского научного центра РАН. Уфа, 2007. вып. 5. С. 113-120.

11. Р.Х. Болотновой. Дисс. докт. физ.-мат. наук, Ин-т механики УНЦ РАН, Уфа, 2008, 312 с.

12. A. Prosperetti, Quarterly of Appl. Math., 34, 339-352 (1977).

13. H. Lin, B.D. Storey, A.J. Szeri, J. FluidMech., 452, 145-1б2 (2002).

14. А.А. Аганин, Д.Ю. Топорков, Т.Ф. Халитова, Н.А. Хисматуллина, Математическое Моделирование, 23, 10, 82-9б (2011).

15. М. Баранаев, ЖФХ, 13, 11, 1бЗ5-1б41 (1939).

16. С.К. Годунов, А.В. Забродин, М.Я. Иванов, А.Н. Крайко, Г.П. Прокопов, Численное решение

многомерных задач газовой динамики. Наука, Москва, 197б, 400c.

17. Г. Флинн, В кн. Физическая акустика, Мир, Москва, 19б7, С. 7-138.

18. M.S. Plesset, T.P. Mitchell, Quart. Appl. Math., 13, 4, 419-430 (195б).

19. G. Birkhoff, Quart. Appl. Math, 13, 451-453 (195б).

© А. А. Аганин - д.ф.-м.н., проф., зав. лаб. Вычислительной динамки сплошной среды ИММ КазНЦ РАН, aganin@kfti.knc.ru; Д. Ю. Топорков - к.ф.-м.н., с.н.с. той же лаборатории, top.dmtr@gmail.com/