Исследование нестационарных плоскопараллельных течений вязкой несжимаемой жидкости (приближение Стокса) методами R-функций и Галеркина Текст научной статьи по специальности «Математика»

УДК517.95:519.63

ИССЛЕДОВАНИЕ НЕСТАЦИОНАРНЫХ ПЛОСКОПАРАЛЛЕЛЬНЫХ ТЕЧЕНИЙ ВЯЗКОЙ НЕСЖИМАЕМОЙ ЖИДКОСТИ (ПРИБЛИЖЕНИЕ СТОКСА) МЕТОДАМИ R-ФУНКЦИЙ И Т А ЛЕРКИН А

АРТЮХА.В., СИДОРОВ М.В._________________

Рассматривается задача расчета нестационарного плоскопараллельного течения вязкой несжимаемой жидкости в конечной односвязной области в приближении Стокса. Для описания течения используется уравнение для функции тока. На основании метода R-функций и проекционного метода (метод Галеркина) строится приближенный метод решения этой задачи. Эффективность разработанного численного метода иллюстрируется вычислительными экспериментами.

Введение

Актуальность задачи. Изучение законов движения жидкости играет важную роль в развитии техники и естествознания. Исследования в этой области стимулируются потребностями авиации, кораблестроения, теплоэнергетики, геофизики, биологиии пр. За последние десятилетия сфера исследования и применения явлений, связанных с движением жидкости, постоянно расширяется и охватывает ведущие направления промышленности (химические технологии, нефте- и газоразработка, металлургия и т. д.) и ряд естественных наук (биология, физика атмосферы и океана и др.). Во многих практически важных случаях жидкость можно с большой достоверностью считать вязкой несжимаемой ньютоновской средой, и проходящие в ней процессы могут быть промоделированы с помощью уравнений Навье-Стокса [1,2]. Различные задачи, возникающие при изучении динамики вязкой жидкости, могут быть исследованы теоретическим путем или с помощью физического эксперимента. Однако с развитием ЭВМ все активнее используется математическое моделирование. Существует множество численных методов, применяемых при расчете вязких течений. Литература по этому направлению обширна [3-5 и др.]. В основном эти численные методы используют метод конечных разностей и метод конечных элементов. Они просты в реализации, но не обладают необходимым свойством универсальности - при переходе к новой области (особенно неклассической геометрии) необходимо генерировать новую сетку, а часто и заменять сложные участки границы простыми, составленными, например, из отрезков прямых. Точно учесть геометрию области можно, воспользовавшись конструктивным аппаратом теории R-функций, разработаной акад. В.Л. Рвачевым и его учениками [6,7 и др.]. Задачи гидродинамики решались в работах С.В. Колосовой, К.В. Максимен-ко-Шейко, И.Г. Суворовой. Т.И. Шейко, М.В. Сидорова и др. [8-11, 17]. Однако в основном рассматри-

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

Цели и задачи исследования. Целью настоящего исследования является разработка новых средств математического моделирования и численного анализа нестационарных плоскопараллельных течений вязкой несжимаемой жидкости в конечных односвязных областях в приближении Стокса на основании методов R-функций и Галеркина. Для достижения поставленной цели необходимо решить следующие задачи:

- получить полную структу р} решения начальнокраевой задачи для функции тока, используя метод R-функций;

- разработать и обосновать алгоритм аппроксимации неопределенной компоненты полученной структу ры на основании метода Г алеркина;

- провести вычислительные эксперименты.

1. Постановка задачи

Нестационарное течение вязкой несжимаемой жидкости описывается хорошо известными уравнениями Навье-Стокса [1]:

dv _ 1 4_

----r(vv)v = —Vp + vAv (])

<5t р

и уравнением неразрывности

divv = 0- (2)

Здесь v - поле скоростей; р - давление; v - кинематическая вязкость; р - плотность. Будем предполагать, что объемные силы отсутствуют.

Решение системы (1), (2) сопряжено со значительными трудностями, связанными, в основном, с присутствием в (1) нелинейного члена (vV)v .

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

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

16

РИ, 2011, №3

В приближении Стокса уравнения нестационарного плоскопараллельного движения вязкой несжимаемой жидкости в конечной односвязной области Q плоскости хОу имеют вид

3v 1 др (д2\т д2\, N

—~---------I-V —— н-----f-

dt р дх ч дх ду / ’

£v,

а

1 <5р ( д\

—— + v ---f-

р ду ^ дх

fV

ду2 ,

(3)

дх

ду

2. Применение методов R-функций и Г алеркина

Для решения начально-краевой задачи (4) - (6) используем методы R-функций и Г алеркина.

Пусть граница дй1 области Q кусочно-гладкая и может быть описана элементарной функцией со(х. у) согласно методу R-функций [6], причем функция со(х,у) удовлетворяет условиям:

1) со(х,у) = 0 на дО;

2) со(х,у) > 0 в Q ;

СЛ|/ _

3)тлг - - 1 на сії. т.с. о)(\. у) = 0 - нормализованное

Анализ плоскопараллельных течений удобно производить с помощью функции тока \|/(\. у), вводимой соотношениями

уравнение <3Q, п - внешняя к д<} нормаль.

В работе [15] было показано, что краевым условиям (5) удовлетворяет пучок функций

_ А|< _ cty

" ду ' у дх

(уравнение неразрывности при этом обращается в тождество).

Исключая из (3) дифференцированием давление, для функции тока получаем уравнение

5Ду 2

---- = vA viz

dt ' '

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

v. Так, если жидкость примыкает к неподвижной стенке, то в этих точках скорость жидкости обращается в нуль. Это означает, что в нуль обращается нормальная и тангенциальная составляющая скорости (условие прилипания). Если же жидкость примыкает к подвижной твердой стенке, то в таких точках скорость жидкости должна по величине и направлению совпадать со скоростью соответствующей точки стенки. Исходя из этого, на границе дО. области О можно задать значение функции тока ф и её нормальной А|/

производной -Tjr, где п - внешняя нормаль к дО.

Итак, для функции тока г[/(х,у) можно поставить начально-краевую задачу

= уД\|/ , (х.у) є Q, t > 0 , (4)

4'L=f0(s.t),

go(s-t), s є дО, t > 0, (5)

4=0 = Vo(х.У), (х,у)єП. (6)

Методика задания функций f0(s,t) и g0(s,t) рассмотрена в [15].

\|/ = Г — ra(Djf + g) + ю2Ф, (7)

где f = ECf0, g = ECg0 -продолженияфункций f0, g0

в О соответственно,

„ да dv да <5v „ „

D.v =------1-----= (Vco,Vv)

' дх дх ду dy ’

Ф = Ф(х, у, t) - неопределенная компонента, которую будем предполагать достаточно гладкой.

В задаче (4) - (6) сделаем замену

у = ср + и ,

где ср = f -co(D,f +g), и -новаянеизвестная функция. Тогда для функции и получим начально-краевую задачу с однородными краевыми условиями:

д( Au)+vA^ = f^ (х,у)єп, t>0, dt

4 =0.^1

= 0.

(8)

(9)

u

lt-0

= Un

(10)

^ ., ЗДф і

где F = -уД'ф + ——, и0=ч/0-ф| at

Для решения задачи (8) — (10) применим метод Fалеркина для нестационарной задачи [22].

Пусть Т > 0, Н - сепарабельное гильбертово пространство. Символом Ц(0.Т:Н) будем обозначать множество функций u(t), t є [0.11. со значениями в н таких, что

}||u(t£dt<+°0

о

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

т

< u.v >= J(u(t).\(t))„dt

о

РИ, 2011, №3

17

Возьмем H = L2(Q). Пусть u0 є L2(Q), F(t) є є L2(0,T;L2(Q)) .

Введем в рассмотрение операторы А и В, действующие в L,(0,T;L2(Q)) по правилам

Au = A2u , Ви = -Ди на областях определения

D

А

г _

^иІиєСЧоЮСЧО). и|м

D

А

(

■j и | и є С2 (Q) П С1 (Q), и|Л2 І

Аі

Ді

= 0

Можно доказать [19], что операторы А и В будут положительно определены. Ясно, что DA cDb.

Тогда задачу (8) - (10) можно записать в операторной форме

— Bu + vAu = F, (х,у) є Q, t > 0 , (П)

dt

u|t=0=u0. (12)

На Da введем энергетическое произведение [u. v]A по правилу: для любых u, v є DA

[и, v] А = (Au, v)4 (n) = jj Au • Av dxdy ^

Q

а соответствующая энергетическая норма

lull =|[(Au)2 dxdy

Q

Пополняя Da в норме I u|A, получаем энергетическое пространство НА оператора А

На DB введем энергетическое произведение [u,v]B по правилу: для любых u.v є DB

[u,v]B =(Bu.v)4(Q) = j|Au-vdxdy

Q

Применяя формулу Грина [19] и учитывая краевые условия, получаем

[u.v]B =(Bu.v)4(Q) = JJVu-Vvdxdy ^

а соответствующая энергетическая норма

Можно показать, что НА = W2 (Q) с Нв.

Пусть u(t) — классическое решение задачи (11), (12), т.е. для любого t > 0 utt) є Da . utt) непрерывно дифференцируема по t, удовлетворяет уравнению (11) и начальному условию (12).

Пусть v(t) - достаточно гладкая в Пх[0,+х) функция, удовлетворяющая краевым условиям (9) и такая, что при некотором Т>0 \(Т) = (). Умножим (11) скалярно в L2(Q) на произвольную функцию v(t) с указанными свойствами:

(]HBu'v ] +v(Au.v)L (n)=(F.v)L2(n) _

Vul A2(Q)

Интегрируя последнее равенство по t от 0 до т, получаем, что

1 [“Г Bu-v I dt + vJ (Au. v)L (Q) dt = } (F, v)4 (Q)dt

O' - U(Q) 0 0

Если проинтегрировать первый интеграл по частям (по переменной t) и воспользоваться равенством v(T) = 0 , то получим, что

-(Bu0.v(0))U(Q)-J Ви,— I dt + vJ(Au.v), lUldt =

oV ®t/L2(Q) 0

= j(F-v)4(n)dt

о

Учитывая вид энергетических произведений в НА и Нв , последнее равенство перепишем в виде

-(

dv

а"

т

dt + v[|u.v|Adt =

В 0

= [u„.v(0)|B 4-J(F.Y)MQ;dt (О)

0

Последнее равенство возьмем в качестве определения обобщенного (слабого) решения задачи (11), (12) (а значит, задачи (8) - (10)).

Обозначим множество функций

WT ={u|ueL2(0,T;HA),

и' є L2(0,T;L2(Q)), u(T) = 0} .

|u|! =j||Vu|2 dxdy

Q

Пополняя DB в норме |u|B , получаем энергетическое пространство Нв оператора в •

Определение. Функция u(t) называется обобщенным (слабым) решением задачи (11) - (12), если

a) u(t) є L2 (0, Т; W22 (Q));

18

РИ, 2011, № З

б) для любого элемента v(t)eWT имеет место равенство (13).

Справедлива следующая теорема.

Теорема 1. Пусть F(t) є L2 (О, T;L2(Q)) и и0єНв. Тогда существует, и притом единственное, обобщенное решение задачи (11), (12).

Для построения обобщенного решения задачи (11)-(12) воспользуемся методом Галеркина [22]. Приближенное решение задачи (11) — (12) ищем в виде

п

и„0) = ХМ1>Фк , (14)

к -1

где ск (t), к = 1.п - неизвестные пока функции, (срк}

- координатная последовательность, т.е. последовательность (фк} удовлетворяет условиям:

1) для любого к срк є НА ;

2) для любого п ср..фп линейно-независимы;

3) (фк} полна в Нд .

Поскольку из (7) следует, что и = со2Ф, где Ф = Ф(х, у, t) - неопределенная компонента структу ры, то координатную последовательность можно взять

в видефк =со2тк, где (тк} - любая полная bL2(Q) система функций.

В соответствии с методом Галеркина неизвестные функции с, (t). k = l,...,n , найдем из условия ортогональности невязки, получаемой при подстановке (14) в уравнение (11), первым п координатным функциям

Ф1..фп. Это приводит для определения ck(t),

k = 1,...,п , к системе обыкновенных дифференциальных уравнений

п п

Хск(1)[фк.ф1]в -i-vVck(t)[cpt.cp.L =

к=1 к=1

(к-cpj - j = 1.2..n. (15)

Систему (15) нужно дополнить начальнымиусловиями

ск(°) = ск> к=1,2,...,п. (16)

Начальные условия (16) можно задать различными способами [22], например, решением системы алгебраических уравнений

п

Хск(°ХФк>Ф])мо) =(uo,9j)L,(Q), j = 1.2.n, (17)

к=1

которая получается из условия ортогональности невязки начальных условий (12) первым п координатным функциям cpj срп .

В силу условий, наложенных на координатную последовательность (срк} , система (17) и задача Коши (15), (16) при любом п имеет единственное решение. Справедлива следующая теорема.

Теорема 2. Приближенные решения un(t) задачи (11) - (12), построенные по методу Галеркина, определены однозначно при любом п, причем

о

un(t)->u(t), п->оо,слабов L2(0,T;W22(Q)),

где u(t) - обобщенное решение задачи (11), (12).

3. Результаты вычислительного эксперимента

Рассмотрим задачу (4) - (6) для прямоугольной области £2 = {(х,у)|0 <х<1, 0 < у < і] со следующими краевыми условиями:

дДш . , , . „

—4- = vA-vp, (х,у) є Q, t > 0 , dt

'I'L

О

1-е

О,

у = 1;

х = (). у = 0. X = 1.

= 0.

Решение поставленной задачи найдено с помощью методов R-функций и Галеркина. В структуре решения (7) нормализованное уравнение Q имеет вид

со(х,у) = [ х(1 — х)] ла [у(1 —у)] = 0 ,

где ла - R-конъюнкция [6].

В качестве базисных функций выбирались степенные полиномы, тригонометрические полиномы, полиномы Лежандра. При вычислении интегралов в скалярных произведениях в системах (15)и(17) использовалась формула Гаусса с 16 узлами по каждой переменной.

На рис. 1-4 построены линии уровня функции тока.

Рис. 1. Линии уровня функции тока при t = 0,2

Рис. 2. Линии уровня функции тока при t = 0,5

РИ, 2011, №3

19

Рис. 3. Линии уровня функции тока при t = 1

Рис. 4. Линии уровня функции тока при t = 1.8

На рис. 5-8 приведены линии уровня вихря Q = —Ац/ в разные моменты времени.

Рис. 5. Линии уровня вихря при t = 0.2

Рис. 8. Линии уровня вихря при t = 1,8 На рис. 9 показан график изменения max vi/i х. v.t'i. Каквидно, в задаче существует стационарный режим, ^тах

Рис. 9. Изменение максимума функции тока по времени

сН|/

Нарис. 10 приведены значения скорости ГГ в

UXl у=0,5

разные моменты времени.

в разные

моменты времени t = 0,2; 0,5; 1; 1,8

Полученные результаты хорошо согласуются с результатами физических экспериментов [1, 2] и результатами, полученными другими авторами. Расхождения составили около 3%.

Рис. 6. Линии уровня вихря при t = 0,5 1.Qi

0.8

0.6

J \\ ч:

\ /

у

Я-

0.0 0.2 0.4 0.6 0.8 1.0

Рис. 7. Линии уровня вихря при t = 1

Выводы

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

20

РИ, 2011, №3

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

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

Литература: 1. Ландау Л. Ф., Лифшиц Е.М. Теоретическая физика. В 10 т. Т. VI. Гидродинамика. М.: Физматлит, 2003. 736 с. 2. Лойцянский Л.Г. Механика жидкости газа. М.: Дрофа, 2003. 840 с. 3. Роуч 77. Вычислительная гидродинамика. М.: Мир, 1980. 616 с. 4 .DoneaJ., Huerta A. Finite Element Methods for flow problems. London: Wiley, 2003. 350 p. 5. Zienkiewicz O.C., Taylor R. L. The finite Element Method. Vol. 3: Fluid Dinamics. Oxford: BH, 2000. 334 p. 6. Рвачев В.Л. Теория R-функций и некоторые ее приложения. К.: Наук, думка, 1982. 552 с. 7. Колодяжный В.М., Рвачев В.А. Структурное построение полных последовательностей координатных функций вариационного метода решения краевых задач: Препр. АН УССР. Ин-т пробл. машиностр. Харьков, 1975. 75 с. 8. Суворова И. Г. Компьютерное моделирование осесимметричных течений жидкости в каналах сложной формы // Вести. НТУ ХПИ. Харьков, 2004. №3]. С. 141-148. 9. Колосова С.В. Об обтекании невязкой жидкостью цилиндра в трубе //Прикл. мех., 1971.№7.Вып. 10.С. 100-105.10.Максименко-Шейко К.В. Исследование течения вязкой несжимаемой жидко-

сти в скрученных каналах сложного профиля методом R-функций // Проблемы машиностроения, 200 Г Т. 4, № 3 -

4. С. 108 - 116. 11. Рвачев В. Л., Корсунский А.Л., Шейко Т.И. Метод R-функций в задаче о течении Гартмана // Магнитная гидродинамика. 1982. № 2. С. 64 - 69. 12. Ладыженская О.А. Математические вопросы динамики вязкой несжимаемой жидкости. М.: Мир, 1972. 588 с. 13. Лионе Ж.-Л. Некоторые методы решения нелинейных задач.М.:Мир, 1972. 588с. 14. Темам Р. Уравнения Навье-Стокса. Теория и численный анализ. М.: Мир, 1981. 408с. 15. Сидоров М. В. О построении структур решений задачи Стокса //Радиоэлектроника и информатика. 2002. № 3. С. 52 - 54. 16. Сидоров ХІ.В. Применение метода R-функций к расчету течения Стокса в квадратной каверне при малом числе Рейнольдса // Радиоэлектроника и информатика. 2002. № 4. С. 77 - 78. 17. Колосова С.В., Сидоров М.В. Применение метода R-функций // Вісн. ХНУ. Сер. Прикл. матем. імех. 2003. №602. С. 61 -67.18. Слободецкий Л.Н. Обобщение пространства С.Л. Соболева и их приложение к краевым задачам для дифференциальных уравнений в частных производных. //Уч. зап. Ленингр. гос. пед. ин-та им. А.И. Герцена. 1958. Т. 197. С. 54-112. 19.Мишин С.Г. Вариационные методы в математической физике. М.: Наука, 1970. 512 с. 20. ФедотоваЕ.А. Атомарная и сплайн-аппроксимация решений краевых задач математической физики: Дис. ... канд. физ.-мат. наук: 01.01.07. Харьков, 1985. 170 с. 21. Федотова Е.А. Практические указания по использованию сплайн-аппроксимации в программирующих системах серии «Поле»: Препр. АН УССР. Ин-т пробл. машиностр. ;202. Харьков, 1984. 60с. 22 .МихшнС.Г. Численная реализация вариационных методов. М.: Наука, 1966. 432 с

Поступила в редколлегию 24.08.2011

Рецензент: д-р физ.-мат. наук, проф. Колосов А.И.

Артюх Антон Владимирович, аспирант, ассистент кафедры Прикладной математики ХНУРЭ. Научные интересы: математическое моделирование, математическая физика, численные методы. Адрес: Украина, 61166, Харьков, пр. Ленина, 14, тел. 702-14-36.

Сидоров Максим Викторович, канд. физ.-мат. наук, доцент кафедры Прикладной математики ХНУРЭ. Научные интересы: математическое моделирование, математическая физика, теория R-функций и ее приложения. Увлечения и хобби: история культуры. Адрес: Украина, 61166, Харьков, пр. Ленина, 14, тел. 702-14-36.

РИ, 2011, №3

21