Численный анализ ветвлений решения задачи Хеле-Шоу о стационарно движущемся пузыре Текст научной статьи по специальности «Математика»

УЧЕНЫЕ ЗАПИСКИ КАЗАНСКОГО УНИВЕРСИТЕТА

Физико-математические пауки

УДК 532.546^519.63

ЧИСЛЕННЫЙ АНАЛИЗ ВЕТВЛЕНИЙ РЕШЕНИЯ ЗАДАЧИ ХЕЛЕ-ШОУ О СТАЦИОНАРНО ДВИЖУЩЕМСЯ ПУЗЫРЕ

М.М. Алимов

Аннотация

Для стационарной задачи о продвижении пузыря в лотке Хеле-Шоу одно из проявлений нерегулярности предела исчезающе малых капиллярных сил состоит в вырождении решения для идеализированной задачи вместо ожидаемого единственного решения получается целое семейство. С помощью численного анализа С. Тапвир (в. Тапуеег) показал, что учет капиллярных сил устраняет вырожденность решения, по по дал этому ясного объяснения. Кроме того, помимо главной ветви решения С. Тапвир получил и другие. С целыо выявления всех ветвей решения задачи о продвижении пузыря в настоящей работе сформулирована модифицированная задача по аналогии с подходом Ж.-М. Вапдеп-Врека (Л.-М. Уаш1еп-Вгоеск) к задаче о продвижении пальца. Проведен ее численный анализ, который показал, что в постановке с заданной площадью пузыря решение задачи единственно и совпадает с главной ветвыо, полученной С. Тапвиром. Никаких других ветвей решения не выявлено. Такое несогласование с результатами С. Тапвира можно объяснить тем обстоятельством, что его методика допускает решения с неоднолистностью физической плоскости. Дапо также объяснение факту устранения вырождения решения: у области течения имеется две характерные точки бесконечности слева и справа, в которых область имеет заданные размеры. Обе эти величины определяются однократным интегрированием основного граничного уравнения. Поскольку решения пеидеализировап-пой задачи по могут быть симметричны относительно вертикальной оси, удовлетворить оба условия с помощью одной копстапты интегрирования невозможно. Возникает условие разрешимости.

Ключевые слова: задача Хеле-Шоу со свободной границей, поверхностное натяжение, стационарное движение пузыря, итерационный метод.

Введение

Течение вязкой жидкости в лотке Хеле-Шоу при наличии второй связной фазы (жидкости другой вязкости) является простой лабораторной моделью многих природных и технологических процессов, отличающихся взаимодействием между собой нескольких фаз с эволюционирующей границей между ними [1]. С использованием процедуры осреднения характеристик процесса по зазору лотка математическая модель такого двухфазного течения сводится к плоской задаче Хеле-Шоу со свободной границей [2]. Различают нестационарную и стационарную задачи. Для стационарных задач Хеле-Шоу вообще и задачи о пузыре, продвигающемся в лотке типа канала, в частности, наибольший интерес вызывают различные проявления нерегулярности предела исчезающе малых капиллярных сил. одно из которых вырождение решения [3]. Если обозначить через Б относительную площадь пузыря (в плане), через и относительную скорость его продвижения, через е относительную величину капиллярных сил, то для каждой заданной пары величин Б и е > 0, как показывает эксперимент, наблюдается пузырь определенной конфигурации и определенной скорости продвижения и(е, Б), в том числе и в случае е ^ 1.

В то же время в идеализированном случае (е = 0) для любой заданной площади пузыря Б вместо единственного решения анализ дает целое семейство точных аналитических решений с параметром семейства и > 1 [4]. Численный анализ, проведенный в работах [5, 6], показал, что учет капиллярности е > 0 устраняет вырожденность решения, поскольку приводит к отбору определенной скорости продвижения пузыря и(е, Б). Однако с математической точки зрения ясного объяснения этому факту до сих пор пет ведь по сути никаких новых условий в постановку задачи не добавляется, просто одно из граничных условий становится нелинейным за счет включения в него кривизны неизвестной границы.

Автор [5] дает такое объяснение: при заданных параметрах е > 0, Б, и итера-

и

то можно добиться их сходимости. Анализ используемого в работе [5] итерационного метода дает более содержательную информацию. Метод включает в себя использование дискретного преобразования Фурье, причем несколько нетрадиционным образом - количество коэффициентов с0, С1, ..., cN_1 разложения функции в ряд Фурье на единицу меньше числа N +1 точек, в которых функция вычисляется. Фактически такое использование преобразования Фурье эквивалентно традиционному использованию дискретного преобразования Фурье [7] при наложении линейного относительно функции условия сN =0. Такого рода условие эффективно используется для подавления роста старшей гармоники в итерационных методах решения задач со свободной границей [8]. В то же время его можно

и

решение было эффективно получено в [5], такая стратегия себя вполне оправдала. Вместе с тем объяснить с ее помощью факт устранения вырожденности решения невозможно. Помимо этого решение оказалось неедннственным, а сама методика получения другой ветви решения задачи (будем называть ее ветвыо С. Танвира) изложена весьма сжато и не вполне ясно. В следующей работе [6] без каких-либо подробностей сообщается, что была найдена еще одна ветвь решения, аналогичная ветви Ж.-М. Ваиден-Брека [9] в задаче Хеле-Шоу о стационарной конфигурации пальца.

В соответствии с вышесказанным цель настоящей работы построить отличный от используемого в [5, 6] метод решения задачи о стационарной конфигурации пузыря, продвигающегося в лотке Хеле-Шоу, с учетом капиллярных сил, который позволил бы, во-первых, объяснить качественное отлично неидеализированной математической постановки от идеализированной: во-вторых, прояснить вопрос с неединственностью решения. Основная идея предлагаемого альтернативного метода состоит в переходе к некоторой модифицированной задаче с дополнительным свободным параметром, аналогично подходу Ж.-М. Ванден-Брека к задаче о стационарном пальце [9]. В работе [9] в качестве дополнительного свободного параметра выступала величина угла наклона касательной к поверхности пальца в его вершине. В рассматриваемом случае в качестве такого параметра предлагается взять А - отношение ширины канала на бесконечности справа к ширине канала на бесконечности слева. После построения численного решения модифицированной задачи в широком диапазоне определяющих параметров, все ветви решения интересующей нас задачи о пузыре в лотке типа простого канала получаются наложением ограничения А = 1.

1. Постановка модифицированной задачи и сведение ее к граничному уравнению

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

А

У = 1/2

а)

В

'-..= 1

12

аУ ____„Г

В

С В

-1 -а

И

(

А Б

0 а )

Рис. 1. Вид физической плоскости г и вспомогательной плоскости £

тральной оси канала пузыря воздуха (ввиду симметрии изображена только верхняя половина течения). Дополнительный параметр А > 0 характеризует относительную ширину канала на бесконечности справа.

Определяющие параметры процесса известны [4]: это толщина щелевого зазора Ь и ширина канала та бесконечности слева Ь; вязкость жидкости м и скорость ее движения на бесконечности слева V*; коэффициент поверхностного натяжения на межфазной границе Г; площадь пузыря Б* (в случае А < 1 под Б* понимается удвоенная площадь фигуры ВСС', изображенной та рис. 1; в случае А > 1 необходимо выполнить аналогичное достраивание точки В' - проекции точки В на продолжение оси ДС). Необходимо найти конфигурацию Г пузыря и скорость его продвижения и*.

Реализовать обобщенный случай А = 1 в лабораторных условиях практиче-

А < 1

степку СД надо будет двигать со скоростью пузыря (в случае А > 1 надо будет двигать стенку АВ). Вместе с тем такой процесс вполне представим и ему можно сопоставить корректную математическую постановку. Последнюю выпишем сразу

Ь

и скорость V*. Безразмерные декартовы координаты (ж, у) связаны с пузырем и представляют собой относительную систему отсчета в противоположность абсолютной системе отсчета, связанной с неподвижной стенкой канала, например АД. Как известно, для течений Хеле-Шоу всюду в области течения выполняются уравнения [2, 4]

V = У^, Д^ = 0, (1)

где у(ж, у) - безразмерный потенциал течения. Вследствие этого целесообразно введение комплексной физической плоскости г = ж+гу и комплексного потенциала течения Ш(г) = у(ж, у) + гф(ж, у), где ф(ж, у) - функция тока [10].

Аналогично формулировке стационарной задачи Хеле-Шоу [4, 5] о продвижении пузыря в простом канале, для течения Хеле-Шоу, представленного на рис. 1, к уравнениям (1) необходимо добавить граничные условия

—то

дж

1;

АБ и АВ и СБ : ^ ду

то

0; Г :

дж

¥> =

= V;

—еК, ф = иу,

(2)

где К - кривизна границы Г е

капиллярный параметр следующего вида Ь2

е=

12Ь2 V

Ги

в зависимости от величины е и заданной площади пузыря Б.

и*/к

Покажем, что параметр V - величина скорости течения на бесконечности справа не свободен, а находится из балансовых соотношений. Для этого помимо скорости абсолютного движения жидкости V введем и скорость относительного движения жидкости V (движения относительно наблюдателя, связанного с пузырем). Эти скорости связаны очевидными соотношениями

«Х = ух — и, «у = Уу. (3)

В силу несжимаемости жидкости для относительного ее движения суммарные расходы слева и справа на бесконечности должны быть одинаковы: 1 —и = (V—и)А. Тогда

V = А-1 + (1 — А-1)и.

В случае А = 1 (простой канал) получается естественное соотношение V = 1.

Используя традиционный для задач со свободной границей прием параметризации [11], будем искать решение краевой задачи (1), (2) в виде конформного отображения г(£), полагая, что во вспомогательной плоскости £ облает и отвечает область канонического вида. В качестве последней, аналогично [5], выберем верхнюю половину единичного круга (см. рис. 1, б). Конформное отображение г(£) нормируется заданием соответствия двух граничных точек В и С, а также требованием симметричного расположения граничных точек А и П, которым отвечает С = а и £ = — а. Величина а € (0,1) является неопределенным вспомогательным параметром.

Аналогично [5] для функции ) можно выписать соотношение

*(С) = и-1 [ш(С) — Ш(С)] , (4)

вытекающее из формулы (3). Здесь Ш(£) - потенциал абсолютного движения жидкости, Ш(С) - потенциал относительного движения жидкости. Потенциал абсолютного движения жидкости ищется в виде суммы

Ш (С ) = Шс(С) + иш (С), (5)

где Ш0 (С) - потенциал абсолютного движения жидкости идеализированной задачи; ш (С) - регулярный в области добавок.

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

— а1 С + а

(6)

2п (С + а)(С + а-1)' Подстановкой формул (5), (6) в соотношение (4) получим аддитивное представление для искомой функции

*(С) = го(С)+ "(С), (7)

где го (С) - функция, отвечающая решению идеализированной задачи с теми же

Аиа

1 ?и-1 — 1

*о(С) = ^ 1п(С - «) - " ^ 1п(С - а-1) -

- А1п(С + а) + 2^1 2 + А1п(С + а-1). (8)

Эта функция содержит все особенности решения неидеализированной задачи в области , а с учетом ее симметричного продолжения через действительную ось всюду в единичном круге |£| < 1. Регулярный в этом круге добавок ш (£) ищется в виде ряда

ш (0= со + ^ ейСк, 1т ек =0 (9)

к=1

с действительными коэффициентами ск , к = 0,..., то, поскольку ш (£) на СДАВ должна удовлетворять условиям 1т ш (£) =0. Аналогично [5], для определения ш (С) го первого граничного условия на свободной границе Г в постановке (2) можно вывести граничное уравнение

с = е" : Ие ш (С) = —еи-1К(а), (10)

где К (а) — кривизна искомой границы Г. В терминах функции ,г(£) она имеет вид

Z = : K(а) =

dz d(

Отметим, что Im w(Z) восстанавливается по значениям Re w (Z) на единичном круге Z = ега с помощью интегрального оператора Гильберта [12]. Поэтому кривизна K(а) по сути является нелинейным интегро-дифференциальным оператором К от Re w (Z). Соответственно, граничное уравнение (10), (11) может быть записано в операторной форме

Z = : Re w (Z) = - eU[Re w (Z)], (12)

Таким образом, краевая задача со свободной границей (1), (2) свелась к граничному интегро-дифференциальному уравнению (10), (11) относительно функции w (Z) вида (9). Необходимо добавить только ограничения на функцию z (Z) : в силу предположения о конформности отображения Z ^ z производная dz /dZ не должна иметь точек сингулярности в единичном круге |Z| < 1, за исключением точек A и D: Z = ±а в которых dz/dZ так же, как dzo/dZ , должна иметь простые полюса с моментом (2п)-1. Использование регулярной структуры (9) функции w (Z) гарантирует отсутствие новых особенностей у функции dz/dZ в единичном круге IZ| < 1 то сравнению с dz0/dZ , однако не гарантирует отсутствия нулей. Поэтому к задаче (7) (11) необходимо добавить ограничение

dz

ICI<1: ^ф 0. (13)

Отметим, что решения, не удовлетворяющие ограничению (13), будут пефшич-

z

2. Организация итераций

В работе [5] для организации итераций использовалась непосредственно операторная форма (12) граничного уравнения (10), (11). В результате задача сводилась к численному дифференцированию функции ш (£) и применению интегрального оператора Гильберта для восстановления 1т ш(£) по значениям Ие ш (£) на единичном круге С = ега. Однако в случае е ^ 1 при асимптотическом анализе линеаризированного варианта граничного уравнения (10), (11) использовался другой подход, в результате которого задача сводилась к определению из линеаризированного уравнения второй производной функции ш (£) и ее численному интегрированию для получения самой функции ш (£). Именно такой подход предлагается

применить и к исходному интегро-дифференциальному уравнению (10). (11). предварительно преобразовав его подходящим образом.

Аналогично принятой в теории струй [10] практике можно сразу обеспечить безусловное выполнение ограничения (13). если в качестве основной формы представления решения использовать функцию Д(С)

ад (4)-ь(<$)=|><*. <»>

которую будем искать в виде ряда с коэффициентами к = 1,..., N. Соответственно, представление (7)—(9) функции ) будет производным от основной формы в том смысле, что производная ¿^/¿С восстанавливается из (14) путем взятия экспоненты от Д(£), а вид функции ) — последующим интегрированием.

Уравнение (10), (11) преобразуем к виду

а;

(15)

С учетом сказанного выше правая часть этого уравнения представляет собой некоторый интегральный оператор А от И,е {('¿Д/йС}- Соответственно, уравнение (15) может быть записано в операторной форме

С— КеК|НА

п

Именно эту форму уравнения (15) будем использовать для организации итераций. Фигурирующая в левой части уравнения производная (■(¿Д/^С)> регулярная всюду в единичном круге |С | < 1, допускает следующее параметрическое представление

¿Д N

= 1т а* =0. (16)

^ и=1

Параметры представления - пока неизвестные коэффициенты аи, к = = 1,..., N,— будут действительны в силу того, что выражение, стоящее в правой части (15), действительно при действительных причем ао = 0 в силу регулярности функции Д(£) в окрестности точки £ = 0.

Непосредственным интегрированием представления (16) получим

N

Д(С)= ¿о + ¿иск, ¿и = аи/к, к = 1,..., N5 (17)

и=1

где с! о — неопределенная вещественная константа интегрирования. Обозначим через с! ± величины суммы ряда с коэффициентами г1 и в лотках £ = ±а:

N

± ¿и (±а)и. (18)

и=1

В свою очередь, представление (17) функции Д(£) порождает новое представление для выражения £ (¿¿/¿С)

„¿г ¿¿о

«С «С

которое отличается мультипликативным выделением особенностей. Полагая, что справедливы оба представления выражения Z(dz/dZ) _ аддитивное, следующее из (7), и мультипликативное (19), найдем представление выражения Z(du/dZ)

Легко видеть, что оно противоречит регулярному представлению

d °° Z k=1

вытекающему из выражения (9). Действительно, за счет множителя Z(dzo/dZ) у представления (20), вообще говоря, есть особенности в единичном круге |Z| < < 1 - простые полюса в точках Z = ±а. Используя имеющуюся свободу в выборе

d0

вой части (20) в одной точке, а именно Z = а. Для этого надо выбрать

d0 = —d+. (21)

В то же время обеспечить исчезновение особенности во второй точке, а именно Z = —а, невозможно, поскольку, вообще говоря, d_ = d + и, соответственно, будем иметь

du A(ed 0+d- — 1) а

Z--а: С—---• 22

S SciZ 2тг Z + а

Следовательно, у функции u (Z), полученной интегрированием соотношения (20), появляется логарифмическая особенность ln (Z + а). Заметим, что согласно представлению (7) функция u (Z) - аддитивная регулярная часть функции z(Z) • Особенность ln (Z + а) уже выделена и входит в другую аддитивную часть функции z(Z), а именно в zo(Z)• Представляется естественным организовать итерационный процесс, в котором возникающая особенность ln (Z + а) функции u (Z) переносится в функцию zo(Z) так, чтобы их сумма — функция z(Z) — не изменялась. Тогда новая функция unew(Z) будет уже регулярна в области а новая функция znew(Z) будет содержать особенность ln (Z + а) с множителем, отличным от (—A/2n), см. формулу (8). В результате параметр A будет "естественным обра-

A

сформулировать модифицированную задачу. Если в результате итерационного процесса A сходится к некоторому А(то), то в итерационном пределе будет выполнено и условие

lim {d+ — d _ } = 0.

Формула перехода от функции u (Z) к функции unew (Z) имеет вид

^(Z)^(Z) + A(ed0+^"1) [ln(C + «)-ln(Z + «-1)]. (23)

Z

11

,dunew du A(ed 0+d- — 1)

Z- = Z--h Z—--

S dC dC 2тг

Z + а Z + а

(24)

Учитывая оценку (22), легко убедиться, что особенности правой части выражения (24) в точке £ = —а скомпенсированы. Второе слагаемое в квадратных

скобках соотношений (23) внутри единичного круга |£| < 1 является регулярным и добавлено только из соображений симметрии по отношению к формуле (8). Это позволяет образовать новую функцию г0е№(С) того же вида (8) с помощью простой замены А ^ Апе™

Апе№ = А ва °+а-,

^(С) = ¿о(С)1л^л- , так, что сумма функций ш (С) и го(С) остается неизменной:

*(С) = го(С) + Ш (С) = ¿^(С) + шпе-(С).

(25)

(26)

3. Формализация итерационного процесса

Введем обозначение Оо (а, А), До (а, А) для вычисляемых в замкнутом виде про-

изводных

= С-^г-

С= е-

С= ei'

а также следующее обозначение

С(а,Х) = (—

С= ei'

Оо(а,А) + ^ Ъив^, Ъ к = кси.

(27)

и=1

Тогда с учетом представления (17) уравнение (15) можно переписать в итерационном виде + 1 — номер текущей итерации)

N

^а^ еов(ка) =

и=1

и

О(а, А(3))

N

„(¿+1)

+ УЗ с^ еов(кст)

и=1

- КеДо(а,А(3)). (28)

Зная А и с и, к = с предыдущей ^'-й итерации, из уравнения (28)

методом быстрого дискретного преобразования Фурье (далее БПФ) найдем а¿+1) для + 1)-й итерации. Копетапта со"^"^, хотя и входит в правую часть уравнения, неизвестна и может быть определена из того условия, что среднее на интервале а £ [0, 2п) от правой части уравнения (28) с необходимостью должно быть равно нулю вследствие тождества ао = 0 (см. пояснения к формуле (16)). Обозначая среднее посредством (...}, найдем

„(3+1) -

N

О(а, А(3))| £ с(3) еов(ка) + е и-1КеДо(а, А(3))

и=1

( |О(а,А(3) )|

(29)

После определения а'(и3+1) сразу же находятся и коэффициенты !j+1), см. формулу (18). Далее, в соответствии с (23) введем обозначения ¿±+1):

¿±+1)

N

Е^ (±«)и

(30)

и=1

и выберем константу интегрирования Я(°+1) согласно правилу (21): Я

-7(^+1)

-Я

(¿+1)

. Учитывая представления (20). (24). можно выписать задачу для опре-

деления коэффициентов Ь(к+1)

N

Еь

к=1

(3 + 1)егка к е

= Со (а, А(5))

ехр ¿0'+1) +

N

Е^

к=1

'0'+1)егкст

- 1

где

в(5) = А(5)

ехр

/?СЯ

2тг

я0^'+1) + Я- (^+1)

+ а

-1

(31)

(32)

Задача (31) решается методом быстрого преобразования Фурье, причем должно

еделс

(5+1) а с учетом формулы (25) - пара-

выполняться условие Ь°"'+1) = 0. Посте определения коэффициентов ьк"'+1) с Уче'

том формулы (27) найдем коэффициенты с к метр А(¿+1):

А(^+1) = А(5)

ехр

7(^ + 1)

+

(5 + 1)^

(33)

Тем самым шаг итераций завершен.

За начальное приближение выбирались пулевые значения для всех искомых коэффициентов: = 1, ск) = 0, к = 1,..., N. Краевая задача (1), (2) решалась по итерационной схеме (28) (33) в полу обратной постановке задавались три параметра: определяющий параметр е и вспомогательные параметры а, и. В ходе итераций определялись все коэффициенты с к = к = 0,..., N, и ширина ка-

нала А = А(то) (е, а, и) на бесконечности справа. Выход из итерационного процесса

ск

по достижению равенства Я+ = , что эквивалентно установлению величины А. После выхода дополнительно контролировались невязки уравнений (15) и (20). Решение считалось удовлетворительным, если невязки не превышали 10_6. При N от 210 до 213 итерационный метод (28)-(33) дает удовлетворительные результаты решения задачи для а < 0.998 и е > 0.005 (для а, удаленных от предельного значения 1, скажем а « 0.5, метод эффективен и при меньших е).

Отметнм, что для сравнения вычислялись также и повязки непосредственно уравнения (10) с учетом выражения (11), где производные функций вычислялись с помощью конечных разностей на той же равномерной сетке, что использовалась

а

чения, эти повязки на один-два порядка больше повязок уравнений (15) и (20), но

а

ншо единице невязки уравнения (10) растут гораздо быстрее невязок уравнений

(15) и (20), и имеется целый интервал значений а € [0.9, 0.999], где первые уже

10_6

Конфигурацию контура пузыря Г можно найти с помощью параметрического уравнения ж(<г), у(а), которое получается разделением на вещественную и мнимую часть выражения

1 1 — ои _1 А

= Т" 1п(е- - а) +-^-1п(е<<7 - аТ1) - — 1п(е4<7 + а) +

оп оп оп

211-1 -2 + А 2тг

N

1п(ег<т + а"1 )+Е скегка.

к=1

Зная конфигурацию пузыря, можно вычислить и его площадь в плане

Возвращаясь к первоначальной задаче о течении в простом канале (А = 1), покажем, как получить все ее ветви решения. Зафиксируем физический параметр е. В широком диапазоне значений вспомогательных параметров а, и проводим расчеты по методике (28)^(33) с целью вычисления А(е,а, и) и Б(е, а, и). В пространстве переменных и, А для фиксированного е строим поверхность А (5, и). Сечение ее плоскостью А = 1 даст кривую в плоскости переменных и, кото-

е

простому каналу.

Для получения дополнительных точек поверхности А (5, и), помимо описанной выше методики (28) (33), можно использовать еще одну, вспомогательную методику. Обратимся к формуле (8) для ^о(С) • Величины А и 2и-1 входят в качестве слагаемых в множитель при 1п(£ + а-1). Было сделано предположение, что изме-

можно заменить на 1

нение А, согласно формуле (25) равное ДА = А (ей0+ изменение и так, чтобы выполнялось равенство ДА = взять

2и

- и

- 1), Д(2и

). Для этого надо

и

вида (8)

^(С) = ¿о(С

Ди = и п^

2 + и ДА

получить новую функцию г,

new

о

(С) того же

иnew _

2и

2+ иА (е^о+<!- - 1)'

Для получения формулы перехода от функции и (£) к функции и0^ (£) можно использовать то же самое положение (26) о неизменности суммы и (£) + го (С):

а0+^- - 11 С - а-1

- 1п

2п

С + а

-1

_ дш Х(ес1о+с1- - 1)

¿С ~ ¿С +

С - а

-1

С + а

-1

В соответствии с этим в итерационной методике (28) (33) надо заменить задачу (31) для определения коэффициентов Ь<к+1) на следующую

N

^б^е^ = Со (а | А ^)

к=1

N

ехр

-7(^ + 1)

к=1

¿к+1)е4кЧ - 1 +

/3 СЯ еш еш

2тг е}а - а-1 ёа + а.-1

Аи 2и ^

и 0'+1) =

2+ и С^в^ '

(34)

(35)

Численные эксперименты показали эффективность методики (28) (30), (32), (34), (35) итерационного решения краевой задачи (1), (2). Она позволяет проводить расчеты при фиксированных е, а и А с определением соответствующих и = и(то)(е, а, А) и Б(е, а, А). Такая методика удобна для получения точек линий уровня поверхности А (Б, и), а также для решения первоначальной задачи для простого канала. Однако до выяснения характера ветвления решения последней эта методика является все-таки вспомогательной, а методика (28) (33) основной.

Табл. 1

а ет и е А

0.5 0.10 1.9162 0.0013219 1.00000

0.81 0.10 1.8396 0.0036137 1.00025

0.99 0.10 1.7907 0.0055457 1.00001

0.20 0.20 1.9768 0.0004100 1.00118

0.50 0.20 1.8734 0.0027042 0.99997

0.81 0.20 1.7448 0.0076199 1.00005

0.99 0.20 1.6770 0.0118430 0.99999

0.50 0.30 1.8576 0.0040908 0.99998

0.81 0.30 1.6744 0.0119105 1.00002

0.99 0.30 1.5890 0.0187485 0.99980

0.50 0.50 1.8489 0.0068498 0.99999

0.99 0.50 1.4206 0.0349517 0.99998

0.50 1.00 1.8453 0.0137269 0.99999

0.81 1.00 1.6101 0.0412873 1.00002

0.99 1.00 1.3106 0.0757705 1.00000

Табл. 2

а ет и е А

0.5 0.01 1.6155 0.0001568 0.77512

0.78 0.01 1.2046 0.0005117 0.37686

0.80 0.01 1.1530 0.0005623 0.29990

0.82 0.01 1.1270 0.0006045 0.25878

0.84 0.01 1.1109 0.0006436 0.23533

0.50 0.004 1.7953 0.0000564 0.89209

0.50 0.100 1.7142 0.0014776 0.88001

0.50 0.150 1.5684 0.0024224 0.77383

0.50 0.180 1.3878 0.0032852 0.65342

0.50 0.200 1.3072 0.0038754 0.61249

0.50 0.220 1.2615 0.0044174 0.59599

0.50 0.240 1.2322 0.0049334 0.58968

4. Анализ результатов расчетов

В первую очередь было проведено сравнение с приведенными в работе [5] ро-

А

схеме (28) (33). для главной ветви решения (табл. 1) и для ветви Танвира (табл. 2). При этом необходимо пересчитать е по формуле е = ета2 и-1п_2 , поскольку в работе [5] для капиллярного параметра ет использовалась другая нормировка. Как видно, для главной ветви решения условие А « 1 воспроизводится с хорошей точностью. В то же время ни в одном из расчетных вариантов для ветви Танвира условие А « 1 не воспроизводится.

Чтобы получить полную картину ветвей первоначальной задачи о стационарА = 1 е проводона серия расчетов краевой задачи (1), (2) по методике (28) (33) в широком

аи

(34), (35). На рис. 2, а, б в двух ракурсах представлен вид поверхности А(и, 5) для е = 0.02. Поверхность раскрашена в тона серого цвета: чем больше вели-А

А = 0.98 (штриховая линия), А = 1 (сплошная линия), А = 1.02 (пунктирная линия). Поверхность характеризуется наличием небольшого гребня, приводящего к немонотонной зависимости и(5) для всех линий уровня (см. рис. 2, о).

Рис. 2. Вид поверхности Л(и, Б) для £ = 0.02. На поверхности отмечены три линии уровня: Л = 0.98 (штриховая линия), Л = 1 (сплошная линия), Л = 1.02 (пунктирная линия)

Для других значений е качественные отличия вида поверхности А(и, Б) можно продемонстрировать на проекциях линии уровня А = 1 на плоско с ть и, Б (см. рис. 3). С ростом е гребень отодвигается в сторону больших значений Б и уменьшается по величине, в частности, для е > 0.054 от вовсе исчезает. С уменьшением е к гребню добавляется впадина, однако их размер весьма мал. Соответственно, ко-

и(Б)

однако они становятся практически неразличимыми (см. кривую, отвечающую е = = 0.005, на рис. 3).

2 1.9 1.8 1.7

и

1.6

1.5

1.4

Н 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9

Рис. 3. Проекции линии уровня А = 1 на плоскость и, 5 для значений £ = 0.005, 0.01, 0.02, 0.035, 0.054

Теперь можно ответить на вопрос о количестве ветвей решения первоначальной задачи о стационарном пузыре в простом канале, когда А = 1. Если решать эту задачу при фиксированных и и г, то в области малых г, а именно £ < 0.054, у решения помимо главной ветви (участок монотонного убывания и, примыкающий к области малых 5) будут и другие ветви - отдельные участки монотонного

и

[9] решения задачи о стационарном пальце в канале, которая также решается при фиксированных и и г. В отличие от задачи о стационарном пальце, задачу о стационарном пузыре в простом канале можно также решать при фиксированных 5 г

решение. Можно высказать предположение, что найденные в работах [5, 6] ветви решения, отличные от главной, представляют собой те самые нефизичные решения, которые не удовлетворяют условию однолистности физической плоскости (13).

Наконец, вернемся к вопросу о том, как с математической точки зрения объяснить факт устранения вырожденности решения задачи при учете капиллярных сил. Следуя логике разд. 2, легко убедиться, что при фиксированных г > 0, а, и решение задачи можно получить только тогда, когда выполняется два условия

с = ± а : Я(() = 0. (36)

Относительно функции Я(() основное граничное уравнение (15) является дифференциальным уравнением 1-го порядка. Произвол при его интегрировании определяется только одной константой интегрирования. Соответственно, из двух условий (36), вообще говоря, можно удовлетворить только одно. Единственно возможный случай одновременного удовлетворения сразу двух условий (36) - это наличие у решения симметрии относительно вертикальной оси, проходящей через середину пузыря, и, соответственно, представимости функции Я(() в виде разложения в ряд (17) только по четным степеням Но, в отличие от идеализированной задачи, у решения неидеализированной задачи такой симметрии быть не может [5].

Следовательно, одно из условий (36) становится условием разрешимости, которое, собственно, и отбирает определенное значение U(е, S) в неидеализированном случае.

Работа выполнена при финансовой поддержке РФФИ (проект Х- 10-01-00629).

Summary

М.М. Alimov. Numerical Analysis of Solution Branchings of the Hele-Sliaw Problem for Steadily Moving Bubble.

It is well known that a free boundary problem for the bubble steady motion in a Hele-Sliaw cell is lionregular in the limit of zero surface tension. Through this lionregularit.y a degeneracy of the solution appears: for a given area of the bubble P.G. Saffman and G.I. Taylor found a family of exact solutions. S. Tanveor showed that the solution degeneracy is removed by taking into account the effect, of surface tension, but he gave no clear mathematical explanation for such removal. In addition. S. Tanveor found several branches of the bubble solution. To find all solution branches, we have defined a modified Helo-Sliaw problem by analogy with J.-M. Vanden-Broock's approach to the problem of steady fingers. Numerical solution of this modified problem has boon found. A unique solution has boon obtained for a given area of the bubble. This solution coincides with the main branch of S. Tanveor's solution. No other solution branches have boon found. An explanation for this disagreement, is that. S. Tanveor's solutions may include the lionunivalent. physical piano, while we have found only univalent, solutions. In this paper we give a clear explanation for the reasons of the degeneracy removal when surface tension is introduced. In the physical plane the flow domain has two characteristic points: at. infinity on the left, and at. infinity on the right., at. which the domain width is assigned. Both these values have to be defined by single integration of the main boundary equation. With only one integration constant, the two conditions cannot, be satisfied because the bubble contour shape has no fore and aft. symmetry. Thus a solvability condition appears.

Key words: free boundary Helo-Sliaw problem, surface tension, bubble steady motion, iterative method.

Литература

1. Океидои Дон:.P., Ховисои С.Д. П.Я. Кочипа и Хеле-Шоу в современной математике, естественных пауках и технике // Прнкл. матем. и механика. 2002. Т. 66, Вып. 3. С. 515 524.

2. Saffman P.G., Taylor G.I. The penetration of a fluid into a porous medium or Helo-Sliaw coll containing a more viscous liquid // Proc. Roy. Soc. London. Ser. A. 1958. V. 245, No 1242. P. 312 329.

3. Saffman P.G. Selection mechanisms and stability of fingers and bubbles in a Helo-Sliaw coll // J. Appl. Math. 1991. V. 46. P. 137 145.

4. Taylor G.I., Saffman P.G. A note on the motion of bubbles in a Helo-Sliaw coll and porous medium // Quart.. J. Mecli. Appl. Math. 1959. V. 12. P. 265 279.

5. Tanveer S. Tlio effect, of surface tension on the shape of a Helo-Sliaw coll bubble // Pliys. Fluids. 1986. V. 29, No 11. P. 3537 3548.

6. Tanveer S. New solutions for steady bubbles in a Helo-Sliaw coll // Pliys. Fluids. 1987. V. 30, No 3. P. 651 658.

7. Самарский ^4.^4., Николаев E.G. Методы решения сеточных уравнений. М.: Наука, 1978. 592 с.

8. Алимов М.М. Определение конфигурации фронта насыщения слоистого пористого материала // Матем. моделирование. 2010. Т. 22, Л' 6. С. 69 83.

9. Vanden-Broeck J.-M. Fingers in a Hele-Sliaw cell with surface tension // Pliys. Fluids. 1983. V. 26, No 8. P. 2033 2034.

10. Lamb H. Hydrodynamics. Cambridge: Cambridge Univ. Press, 1924. = Ламб Г. Гид-родипамика. М.; Л.: Гостехиздат, 1947. 928 с.

11. Гуревич М.И. Теория струй идеальпой жидкости. М.: Наука, 1979. 536 с.

12. Гахов Ф.Д. Краевые задачи. М.: ГИФ МЛ, 1963. 640 с.

13. Лаврентьев М.А., Шабат Б.В. Методы теории функции комплексного переменного. М.: Наука, 1973. 736 с.

Поступила в редакцию 21.12.11

Алимов Марс Мясумович кандидат физико-математических паук, ведущий научный сотрудник НИЦ «НИИММ им. Н.Г. Чеботарева» Казанского (Приволжского) федерального университета.

E-mail: Mars.AlimovQksu.ru