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

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. №31. С. 141-148. 9. Колосова С.В. Об обтекании невязкой жидкостью цилиндра в трубе // Прикл. мех., 1971. №7. В. 10. С. 100-105. 10. Максименко-Шейко К.В. Исследование течения вязкой несжимаемой жидкости в скрученных каналах сложного профиля методом R-функций // Проблемы машиностроения. 2001. Т. 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 с.

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

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

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

Яловега Ирина Георгиевна, канд. техн. наук, доцент кафедры математики Харьковского национального педагогического университета им. Г. С. Сковороды. Научные интересы: математическое моделирование. Адрес: Украина, 61168, Харьков, ул. Блюхера, 2, тел. 700-35-20.

УДК 517.95 : 519.63

ЧИСЛЕННЫЙ АНАЛИЗ ПРОЦЕССОВ ПЕРЕМЕШИВАНИЯ МЕТОДОМ ^-ФУНКЦИЙ

ГИБКИНА Н.В., РОГОВОЙН.С.,

СИДОРОВ М.В., СТАДНИКОВА А.В.__________

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

Введение

Актуальность исследования. Математическое моделирование и анализ течений вязких жидкостей широко применяется во многих прикладных задачах, в частности в задачах перемешивания. С одной стороны, эта проблема связана с многочисленными применениями в химической, фармацевтической и пищевой промышленности [2, 9, 16]. С другой - перемешивание жидкостей представляет фундаментальную научную проблему, которая тесно связана с современными концепциями хаотической и регулярной динамики [1, 15]. Известно [3, 7], что ламинарные течения при некоторых условиях могут приводить к интенсивному перемешиванию. Такие режимы, получившие название хаотических, являются предме-28

том интенсивного изучения как с теоретической, так и с экспериментальной точек зрения. Большинство методов, используемых при моделировании таких процессов, не обладают свойством универсальности и их сложно применять для «непримитивных» областей. В работах Дж. М. Оттино, Х. Арефа, В.В. Мелешко, Т.А. Дунаевой, Т.С. Краснопольской и др. [3, 5, 7, 9] решалась задача перемешивания для таких простых областей, как круг, полукруг, круговой сектор и т. д., однако для изучения процесса перемешивания в более сложных областях предложенный ими математический аппарат не работает. Точно учесть геометрию области можно, используя конструктивный аппарат теории R-функций, предложенный акад. НАН Украины В.Л. Рвачевым [10]. Поэтому разработка новых методов численного анализа задачи перемешивания, основанных на применении метода R-функций, является актуальной научной проблемой.

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

1) определение поля скоростей течения жидкости (формализм Эйлера):

2) исследование траекторий движения отдельных частиц жидкости (формализм Лагранжа).

Для решения первой части задачи перемешивания необходимо разработать приближенно-аналитический

РИ, 2012, № 3

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

Настоящая работа продолжает исследования, начатые в [4].

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

Рассмотрим плоское квазистационарное течение вязкой несжимаемой жидкости в конечной односвязной области Q с кусочно-гладкой границей дО. Пусть граница dQ области Q состоит из участков SQ, dQ, dQ3, 9Q4, причем 9Q2, dQ4 находятся в состоянии покоя, а , 9Q3 попеременно движутся в противоположные стороны со скоростями vt0p (t) и Ybot(t) соответственно (рис. 1).

9Q4

Решение первой части задачи перемешивания заключается в получении поля скоростей (Vx,vy) в области течения Q . Будем считать, что рассматриваемое течение относится к ползущим, а нелинейными слагаемыми в уравнениях Навье-Стокса можно пренебречь, т. е. можно ограничиться рассмотрением приближения Стокса [6].

Плоское квазистационарное стоксово течение будем описывать с помощью функции тока у (x, y, t), вводимой соотношениями [6]

ду ду

ду , y дх .

vx =

Для функции тока y(x,y,t) можно поставить следующую краевую задачу:

Д2у = 0 в Q , (1)

у»=0, Ш

дО

gtop(t) на 5Qi,

<gbot(t) на 5Q3, (2)

0 на dQ2 U dQ4,

где n - внешняя нормаль к dQ; Д кий оператор,

2

Д2 = _а4_+2 д4 +^

дх4 дx2дy2 ду4

бигармоничес-

Согласно [14], функции gt0p(t) и gb0t(t) имеют вид

gtop(t) “ -(vtop(t)= т1) , gbot(t) = _(vbot(t) т3) ,

где Т1 и Т3 - единичные касательные векторы к дQl

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

Решение задачи (1), (2) можно представить в виде

у (x,y,t) = gtop (t)y 1 (x,y) + gbot (t)y 2 (x,y), (3)

где y1(x,y) - решение задачи

Д2у1 = 0 в Q , (4)

ду = [1 на 5Q1,

У1 Ідо = 0, да = [0 на дQ2 U дQ3 U дQ4, (5)

а у2(х,у) - решение задачи

Д2у2 = 0 в Q , (6)

ду

у2 Ш = 0 , дП

на дQз, дО [0 на дО1U дQ2 U дQ4.

(7)

Решение второй части задачи перемешивания заключается в решении уравнений движения лагранжевой частицы

х (,) = дуХМ, y (t) (8)

ду ’ дх

x(t0) = x0, y(t0) = У0 :

(9)

и в построении и анализе траекторий движения (точка означает дифференцирование по времени).

2. Метод численного анализа

Для решения задач (4), (5) и (6), (7) воспользуемся методами R-функций [10] и Ритца [8].

Известно [14], что структура решения задачи Стокса

у|до=f (s) ду

Д2у = F в Q ,

= g (s^ s єдО;

дО

РИ, 2012, № 3

29

имеет вид

В задачах (4), (5) и (6), (7) сделаем замены

у- f-a(Df + g) + ю2Ф, (10)

Vi = hi + ub (13)

где f - ECf , g = ECg - продолжения функций f , g в Q , оператор Di определяется равенством

Di =

дa д + да д дх дх ду ду ’

ф - неопределенная компонента структуры, а функция а(х,у) обладает свойствами

а) а - 0 на дQ;

где hi - -agi, Ui - а2Фі - новые неизвестные функции, і -1, 2 . После подстановки (13) в (4), (5) и (6), (7) для функций Ui, i -1, 2 , получим задачи с однородными краевыми условиями:

Д2ui - -Д2hi в Q , (14)

Ui ш 0,

дui

дп

=0

дQ

(5)

б) a > 0 в Q ;

да

в) — - -1 на дQ, п - внешняя к дQ нормаль. дп

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

Пусть функции ai(x,y), i -1, 2, 3, 4, таковы, что:

1) ai - 0 на дQi;

2) ai > 0 в Q U (дО \ дQi);

3) "дП"--1 на дQi, п - внешняя к дQi нормаль. Тогда функция

a - a1 А0 a2 Л0 a3 Л0 a4 , (11)

где А0 - знак R-конъюкции:

V2 2

„ U + V ;

обладает свойствами а) - в).

Используя формулу «склейки» [10], для задач (4), (5) и (6), (7) получим

f1 -0, f2 -0,

g a2 Л0 a3 Л0 a4

g1 і ,

a1 +a2 Л0 a3 А0 a4 g a1 Л0 a2 Л0 a4

g2 1 .

a3 +a1 Л0 a2 Л0 a4

Значит, по формуле (10) получим следующие структуры решения краевых задач (4), (5) и (6), (7):

Vi - -agi + a^i, i -1, 2 , (12)

где функция a имеет вид (11), Ф1, Ф2 - неопределенные компоненты структур.

Для аппроксимации неопределенных компонент в (12) воспользуемся методом Ритца [8].

С краевыми задачами (14), (15) свяжем оператор A этих краевых задач, действующий по правилу

A -Д2

на области определения

da -iu

u є c4(q) n c1(q), u| до-Щ

- 0!

дО

Можно доказать [8], что такой оператор будет положительно-определенным. Замкнув множество DA в норме, порожденной скалярным произведением

[u, v] - Ц ДuДvdxdy Q

получим энергетическое пространство Ha . Тогда по теореме о функционале энергии [8] обобщенное решение задач (14), (15) может быть найдено как

ui - arg inf Ц [(Д^2 + 2ДuДhi]dxdy, i -1, 2 .

usHa Q

Согласно методу Ритца [8] приближенные решения этих вариационных задач будем искать в виде

N. N.

ui - a2фi и ui, n - a2фi, N - a2 Z cki)xk - Z ckW ,

к-1

к-1

где i -1, 2 , {тк} - полная в L2(Q) система функций (степенные или тригонометрические полиномы, сплайны и т.п.), фк - a2Tk .

По построению последовательность {фк} является координатной:

1) ф] є DA Vj;

2) VN ф1,ФN линейно-независимы;

3) {фД полна в Ha .

Тогда для определения неизвестных чисел c(i), ..., c(N), i -1, 2 , необходимо решить систему Ритца

N (i)

Z [фъ фj]cki) --(Д^> ^jK j -1, •••, N, i -1,2.

к-1

30

РИ, 2012, № 3

Из теорем сходимости метода Ритца [8] следует, что при N последовательности функций ui N, i = 1, 2,

сходятся к единственным обобщенным решениям краевых задач (14), (15) как в норме L2(Q), так и в норме HA . Это значит, что последовательности функций yA N = hi + ui, N, i = 1, 2 , сходятся в норме L2 (Q) к единственным обобщенным решениям задач (4), (5) и (6), (7). Условия применимости описанного численного метода формулируются в виде условий

Ahi є L2 (Q), i = 1, 2 .

Т аким образом, решение первой части задачи перемешивания получено.

Решение второй части задачи перемешивания представляет собой решение задачи Коши:

x (t)

dy n(x, y,t) 5y

, y (t)

dy N(x,y,t)

dx

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

Рис. 3. Линии уровня функции тока (a = b = 1)

x(t0) = ^ y(t0) = Уо ,

любым численным методом, например, методом Рун-ге-Кутты-Мерсона. Здесь (согласно (3))

УN(x,y,t) = gtop(t)y1, N(x,y) + gbot(t)y2, N(x,y) .

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

Вычислительный эксперимент проводился для прямоугольной полости Q = (0, а) х (0, b), заполненной вязкой несжимаемой жидкостью. Боковые стенки полости находятся в состоянии покоя, а верхняя и нижняя движутся попеременно со скоростями vtop(t) и Vbot(t) (рис. 2). Для модельной задачи можно взять

^x(a -x) л0 -by(b -y)

_ а _ L b J

о>1 = b - y, ®2 = x, ®з = У, ®4 = а - x, при этом свойства а) - в) и 1) - 3) будут выполнены.

v top(t)

b "

v,

bot(t)

Рис. 2. Расчетная область вычислительного эксперимента

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

1.0

0.8

0.6

0.4

0.2

4 < 4 і

4 4 ✓ - - - - - A. 1 1 1

4 * І 4 ¥ 4 4 t f 1

4 4 4 A 1 r 4 t l

4 * \ 4 * r T 4 f J

к 4 4 4 4 4 4 І г

A 4 ( 4 . - % r 4 4 1 i

A 4 4 4 і ** » t 4 4 І i

A / t A і - % 1 4 4 4 ї

A 4 4 4 % V 4 4 4 i

4 і * * A ¥ 4 4 4 4 j

4 4 і V - - - - > r t 4 j

4 > - — — — A i

0.2 0.4 0.6 0.8 — 1.0

Рис. 4. Поле скоростей жидкости (а = b = 1)

Решение второй части задачи перемешивания проводилось в рамках подхода Лагранжа (задача (8), (9)) в несколько этапов: нахождение и исследование типа стационарных точек, построение фазовых траекторий, исследование эволюции линейного и плоского элементов, помещённых в рассматриваемую область.

РИ, 2012, № 3

Рис. 5. Линии уровня функции тока (а = 1, b = 3)

31

V V V V і 4

V V ч ч - -

V V ч -

V V > ч ' -

д V

І 1 4- 4 4 -

к 1 1 т т !

V 1 Ч ч ч *■

1 4 т-

4 4 4 - -

4 4 * -

4 4 4 4 - -

4 4 4 4 ч

0.2 0.4

Jf

ч 1 А А А А т

- - > 4 4 А т

- - - 4 4 А т

- - - 4 А Г 1

4 4 1 1

■ч ч ч ч t ч т

к і 1 к Т т т

- 4 4 4 Т г т

Л 1 к

- - *• Ч ч ч к

- - V Ч ч Л Y

- - ч. Ч ч ч Y

4 Г Ч Л ч ч Y

, ,

0.6 0.8 L0

Рис. 6. Поле скоростей жидкости (а = 1, b = 3)

Уравнения движения частицы интегрировались численно.

Расположение стационарных точек для различных параметров области приведено на рис. 7, 8 ( х -седло; • - центр).

0.0 0.2 0.4 0.6 0.8 1.0

Рис. 7. Стационарные точки (а = b = 1)

Рис. 8. Стационарные точки (а = 1, b = 3)

Были также исследованы эволюции (в течение 10 периодов) линейного элемента - отрезка, образованного четырьмя фиксированными точками, которые соединены между собой. Отрезок помещался в окрестность точки типа «центр» и типа «седло» при а = 1,

32

b = 3. Результаты приведены на рис. 10, 11. Анализируя результаты, видим, что при попадании отрезка в окрестность точки типа «центр» через 10 периодов он переместился незначительно от окрестности этой точки. В то же время отрезок, помещённый в окрестность точки типа «седло», значительно изменил своё положение. Это свидетельствует о наличии хаоса в окрестности гиперболической точки.

Рис. 9. Фазовые траектории маркеров

Рис. 10. Эволюция линейного элемента в окрестности стационарной точки типа «центр»

Рис. 11. Эволюция линейного элемента в окрестности стационарной точки типа «седло»

РИ, 2012, № 3

На рис. 9 приведены фазовые траектории маркеров на сетке с шагом 0,2. Из рис. 9 видно, что в окрестности точек типа «центр» ф азовые траектории представляют собой эллипсы, а вблизи точек типа «седло» наблюдается хаотическое поведение траектории маркера.

Аналогичные выводы можно сделать, проследив эволюцию плоского элемента (прямоугольника, заданного девятью точками). Результаты приведены на рис. 12, 13.

Рис. 12. Эволюция плоского элемента в окрестности стационарной точки типа «центр»

Рис. 13. Эволюция плоского элемента в окрестности стационарной точки типа «седло»

Выводы

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

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

Литература: 1. Андриевский Б.Р., Фрадков А.Л. Управление хаосом: методы и приложения. I. Методы // Автоматика и телемеханика. 2003. № 5. С. 3 - 45. 2. Андриевский Б.Р., Фрадков А.Л. Управление хаосом: методы и приложения. II. Приложения // Автоматика и телемеханика. 2004. № 4. С. 3 - 34. 3. АрефХ. Развитие хаотической адвекции // Нелинейная динамика. 2006. Т. 2, № 1. С. 111 - 133. 4. Артюх А.В., Гибкина Н.В., Сидоров М.В. Об одном методе математического моделирования некоторых процессов перемешивания с помощью метода R-функций // АСУ и приборы автоматики. 2008. Вып. 143. С. 67 - 73. 5. Дунаева Т.А., Гуржий А.А., Мелешко В.В. Перемешивание вязкой жидкости в полукруге при малых числах Рейнольдса / / Прикладна гідромеханіка. 2001. Т. 3 (75), № 2. С. 15 -24. 6. Лойцянский Л.Г. Механика жидкости и газа. М.: Дрофа, 2003. 840 с. 7. Мелешко В. В., Краснопольская Т. С. Смешивание вязких жидкостей. // Нелинейная динамика. 2005. Т. 1, № 1. С. 69 - 109. 8. Михлин С.Г. Вариационные методы в математической физике. М.: Наука, 1970. 512 с. 9. Оттино Дж. М. Перемешивание жидкостей // В мире науки. 1989. № 3. С. 34 - 44. 10. Рвачев В.Л. Теория A-функций и некоторые ее приложения. К.: Наук. думка, 1982. 552 с. 11. Роговий Н.С. Застосування методу R-функцій до чисельного аналізу одного квазістаціонарного процесу перемішування // Наукові праці п’ятнадцятої Всеукраїнської (десятої Міжнародної) студентської Наукової Конференції з Прикладної Математики та Інформатики» (Львів, 5 -6 квітня 2012 р.). Л.: ЛНУ ім. І. Франка, 2012. С. 249 - 251. 12. Роговой Н.С. Математическое моделирование некоторых процессов перемешивания// Материалы XVI международного молодёжного форума «Радиоэлектроника и молодёжь в XXI веке» (Харьков, 17- 19 апреля 2012 г.). Харьков: ХНУРЭ, 2012. С. 548 - 178. 13. Роговой Н.С. Численный анализ одного квазистационарного процесса перемешивания// Научные труды Международной молодёжной научной конференция «XXXVIII Г агаринские чтения» (Москва, 10 - 14 апреля 2012 г.). М.: «МАТИ» - РГТУ им. К.Э. Циолковского, 2012. Т. 1. С. 117 - 119. 14. Сидоров М.В. О построении структур решений задачи Стокса // Радиоэлектроника и информатика. 2002. №3. С. 39 - 42. 15. Табор М. Хаос и интегрируемость в нелинейной динамике. М.: Эдиториал УРСС, 2001. 320 с. 16. Фундаментальные и прикладные проблемы теории вихрей / Под. ред. А.В. Борисова, И.С. Мамаева и М.А. Соколовского. Москва-Ижевск: Ин-т комп. исслед., 2003. 704 с.

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

РИ, 2012, № 3

33

Гибкими Надежда Валентиновна, канд. техн. наук, доц. кафедры прикладной математики ХНУРЭ. Научные интересы: математическое моделирование, оптимальное управление сложными системами, актуарная и финансовая математика. Адрес: Украина, 61166, Харьков, пр. Ленина, 14, тел. (057) 7021436.

Роговой Никита Сергеевич, ст. гр. СА 08-1 ф-та прикладной математики и менеджмента ХНУРЭ. Научные интересы: математическое моделирование, численные методы, теория R-функций и её приложения. Увлечения и хобби: научная фантастика, футбол, бильярд. Адрес: Украина, 61166, Харьков, пр. Ленина, 14, тел. (057) 7021436.

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

Стадникова Анна Викторовна, ассист. кафедры прикладной математики ХНУРЭ. Научные интересы: математическое моделирование, численные методы, теория R-функций и её приложения. Адрес: Украина, 61166, Харьков, пр. Ленина, 14, тел. (057) 7021436.

34

РИ, 2012, № 3