Бабенко - Основы численного анализа (947491), страница 165
Текст из файла (страница 165)
ЕЕ, ЬН вЂ” следы звуковой поверхности, АВ, С — следы поверхности, свободной от граничных условий) мы едпнствепяости решени Если Л < Л„то существует единственное инвариантное множество полупотока --. стационарное решение гь. Оно абсолютно устойчиво. При возрастапни Л появляются и другие ипвариаптные множества, в частности периодическио решения. Если значения параметров Л, о отвечают точке па нейтральной кривой, то для этих значений происходит бифурьация Хопфа возникает периодическое решение задачи (4) — (6).
Выше точки В на нейтральной крпвой бифуркацня будет докрятнческой. а ниже -- закритической. При фиксированном о множества периодических решений образуют непрерывнью ветви, и в частности. если о ) он, то в докритической области должны иметься как устойчивые, так и неустойчивые периодические решения, однако обнаружить численно устойчивые периодические решения пока не удается. Эти периоднческне решения прн увеличении Л должны терять устойчивость, н в принципе должны возникнуть двояко-периодические решения, которые. к сожалению, пока не наблюдались.
При численном решении задачи (4)- (6) удается наблюдать апериоднческие по времени решения, н возникает подозрение, что полу- 801 Э 5. Могоод установления поток У~ имеет инвариантные множества сложной природы, возможно, странные аттракторы. 1,0 0,9 0,8 0,7 б 8 !О Я 10 Рис. 11. Нейтральная кривая для течения Пуазейля ф(ж у, !) = ~~~ аь(у, 1) ехр(йол). ~ь ви Г!одставляя это выражение в уравнение (4), мы пе получим точного равенства. Потребуем, чтобы коэффициенты Фурье левой и правой частей совпадали при ехр(йак) для Й( < и.
Положим в!ам, о! уь ехр(дсох); тогда получим уравнение !ми о — ! 2 д! — !ьаь -'- )ь — —. В !ьаы '!с < и, (7) где !ь = д' /г(у — (гас) . Это уравнение можно записать в виде !ь( — ХГ !ьаь) = — 1ь, д! Опишем численный алгоритм решения задачи (4) — (6). Очень часто алгоритм строится так, что обыгрывается тот либо иной метод решения задачи Дирихле лля уравнения Пуассона. Применять разностный метод или метод конечных элементов в данном случае нецелесообразно, поскольку можно ожидать, что решение будет иметь большую гладкость, несмотря на наличие малого параметра Л ~ при старших производных. В принципе решение задачи аналитично по переменным х, у, но возможны зоны больших градиентов -- погранслоев, где этой аналитичностью воспользоваться затруднительно.
Однако решение спектральной задачи Орра-.Зоммерфельда показывает, что собственные функции не имеют больших градиентов, и поэтому весьма вероятно сделанное выше предположение о гладкости. Апостериори расчетами это пре,лположение подтверждается. Итак, целесообразно построить алгоритм без насыщения. Приближенное решение будем искать в виде 802 Глава 10. Некоторые вонрост численного решснил нраеоъсс задач откуда дат д1 — В 1ьаь = 1 )~,. — сц,(1)сМоу+ (дь(1)айну, (8) где под / мы понимаем оператор, обратный к дифференциальному опе— 1 ратору (ь с нулевыми граничными условиями; функции аь(1), рь(1) будут определяться нз условий (5).
Пусть Т„,(у) — многочлен Чебышева первого рода н у1, ..., у его нули; положим 1,(д) —. Т (у) ((д — гд)Т (уг)), и, (у) =- С (у)(1 — дз)(1— — уг) "л. Аппраксимируем аь(у, ъ) интерполяциопным многочлепом ы ааг(1)ид(р), й т': О, 1=1 и ае(р, 1) ннтерполяционным многочленом р — — + ~ ~ае,(1)иг(у).
у 1=1 Тем самым первое граничное условие (5) будет удовлетворено. В уравнения (8) подставим вместо аь интерполяционный многочлен н потребуем, чтобы невязка обращалась в нуль на множестве узлов интерполяции. Пусть матрица Яь иъиеет внд $ь = ( 1„— (о1с)~и,(у,)) Если 1 оператор ограничения на множество узлов, то, как доказано ВЬ|ШЕ, у(Ь уъ = ЯЬ ГЬ ГДЕ ГЬ = (ую, ..., Тань), уаг = ТЬ(ун 1). ПОЛО1КНЫ 5ь = (сЬ а)гу1....., сЬ ойу ), уь = (сЬ айу1, ..., сЬ акд ) .
Тогда из (8) имеем ЛЛ„ Ж вЂ” дг Буль — — 81 Рк+ оь(1)чь+ пь(1)ую !11 < и, (9) где А1,- .=- (аы, ..., аг, )'. До сих пор величины оь(1), гь(С) остаются неопределенными. Они находятся из требования выполнения второго граничного условия (5). Чтобы уменьшить объем вычислений, система обыкновенных дифференциальных уравнений преобразуется.
Пользуясь симметрией задачи по переменной у и с.сигая т четным, разложим векторы Ль на четную и нечетную компоненты, введя величины аь~ = (аьд+ан~,тг-д)/2 (1 < д < тг12) ъъ 1В ~вам,1 )ъъв . Несложно доказать, что (8~,1) = (зь ) . В самом де- ЛЕ, Зги СООТНОШЕНИЯ ВЫтЕКаЮт ИЗ РаВЕНСтВ В, = В .1 го,т1 . (1, 1 = ь ь 1, 2, ..., т) и аналогичных равенств, справедливых для злементов 803 З 5.
Мегвод установлении матрицы Я ~. Поэтому из системы (9) получаем две системы поряд- ка ш/2: АА' и — Л, 'Я~~А~ — — — Я') 'Е,~ + аьд,, А1 дА ' — В 'да Аь =--Фь) 'Рь +В 1 пг (10) Й~ < и, — В 'Л=В=, = С" + 31 и ь ь+ где~,~=(Р ) ггь.д =(Р ) 'д Уберем верхние индексы х и запишем систему (11) покомпонентно: 1Вьу (о,.6, — рьуВиз = дь~ —: (~%9ьэ~ (12) где дь компоненты вектора Сы Система (12) решается методом предиктор-корректор, причем и па этапе прогноза, и на этапе коррекции производится удовлетвореяие граничным условиям (5).
На этапе прогноза используются формулы Вь — - . ((1 -'; 8ыи.)Вь — шь Вь -01 1 1, (е) ( — ц 1 — бшьу 12 — -- (23дг. — 16дь + бди ) — . 'тодуьу), (13) где ыь — В (Л вЂ” (ой)в)т/12. В последней формуле величины В,1, оь прогнозируемые величины на временном шаге с номером 1 определяются по величинам на времонных шагах с номерами — 2, — 1, О. Уравнение, аналогичное (13), имеет место и для величин (В„), Д, Очевид— % но, что второе граничное условие (5) имеет вид ~в т ав и (1) = О, ~~~ аьуи ( — 1) = О. э=1 1=1 где Р определяется по вектору Рь аналогичным образом, ~„", и„ векторы половинной размерности, образованные из первых т/2 компонент векторов ди, щ.
Несложно доказать, что матрицы яь приводятся к диатональному виду в обп~ем базисе. Пусть Я~ =. Р+Л (Р+) ~, где Льт .= г11а8(Л~ — (ко)з), Л, -- собственные значения матрицы Ве . Положим В = (Рт) 'А„", С„= =- — (Л~ ) ~(Рв — ) Рь . Тогда из (10) получаем диагональную в линейных членах систему 804 Глава 10. Ненашарые аанраеы шмленнога речаенил нраееыв,чадач Легко проверить, что ипа (1) = — и',( — 1).
Поэтому, положив 1 (1 — дг)т' (д,)(1 — д,) ' последние соотношения можно записать в эквивалентном виде: т/г (6 ш6а, )а„=. О. (14) Пусть г, = 6 .Е 6,„те Ц = 1, 2, ..., го~2), г= = (г~, ..., гт, ). Положим ов = геР='-, от = (ее=, ..., ет ): тогда соотношения (14) принимают вид (от, В~) — — О, где (, ) -- скалярное произведение в Итчг. Учитывая (13), имеем трг о = ( -',(Вт)') = ~ ' ~(1--~ лз)В„",'- (15) чп/г — аът В,, — (23д — 1бд + 5д )~ + гон ~ М т ~е>, ~0 „е гч 1 5ьааз 12 ги ьг "а ) 1-- бага.' — г из которого находим тол.
Аналогичное уравнение записывается и для В из которого находим тД. Коррекция выполняется по формулам 1 — бшвз .~. — (5д. ч-8д„- д. ) +тоьетг) 60 ( — П 12 где д„вычисляется по величинам В . Аналогичное уравнение записы- , -% ...,, . 00 1 ввется и для (В, ) . Граничные усчовия удовлетворял>тся аналогичным образом, и они позволяют найти величины тоги тдь. 61ы детально слпюали формулы, за исключением формул, которые служат для определения величин Д ® < и). Эти величины находятся с помощью метода быстрого преобразования Фурье, и детали может восстановить сам читатель, Приведем результаты двух расчетов, вьпюлненных по описанному алгоритму.
В первом расчете отыскивалось периодическое решение при значении параметров В = 9000, о = 0.85 (точка С на рис. 11). В расчете принималось и .= 7, т = 30; он начинался с довольно хаотичных начальных данных. В процесса решения задачи по времени поле скоростей сглаживалось, функция ае(у. 1) — у+уз,г3 стремилась к пулю, аз (д, г) принимала вид ехр( — расее)иг(д), где р1(у) первая собственная функция 805 8 5. Метод установления спектральной задачи Орра — Зоммерфельда, а функции а (д, 1) ('у~ < 2) резко убывали. Затем они возрастали, происходила перестройка функций аь(д, 1) и решение приобретало вид бегущей волны: Полученный в расчете период оказался достаточно большим: Т =- 27,6.
Формула (16) выдерживалась с замечательной точностью на протяжении белое 150 периодов. Подробные таблицы коэффицнеггтов приведены в (22), а здесь приведем только таблицу коэффициента ае(д), пуазейлевского профиля йо(у) н их производных. Из табл. 1 видно, что !и ао/ау ( > ~д 1У,7ад ~ ), т.е. сила трения, испытываемая стенкой, больше в случае периодического решения. О достаточности взятого числа гармоник (и = 7) говорит тот факт, что шах аг(д)~ = 4 10 у а гпах ат(д) ~ = 4,1 10 ". Рис. 13. Траектория жидкой частицы 0,6 0,4 0,2 0 — 0,2 Рис.
12. Линии уровня функции тока й(т, у, 0) !5(л, у, 1) = ~ ~ехр(Ыо(я — с!))ои(у). !ь!<о ь~ььь~йА141111! 18 30 42 54 66 78 90 102 1!4 !26 !38 150 !62 174 1 ь ! тг1Г 806 Глава 10. Ненатерта вопроси численного решения нраееът задач аблица 1. Коэффициент аа(у) и его производная цля периодического решения Т а 4ц* за е у аа ао 0,0523 0,1551 0,2530 0,3431 0,4228 0,4908 0,5462 0,5893 0,6207 0,6420 0 6552 0,6624 0,6655 0,6665 0,6667 0.3960 0,3750 0,.2966 0,2841 0,2061 0,2051 0,1284 0,1358 0,.0670 0,0751 0,0244 0,0291 0,0028 0,0051 0,0523 0,1564 0,2588 0,3584 ОА540 0.5446 0,6293 0,7071 0,7771 0,8387 0,8910 0,9336 0,9659 0,9877 0,9986 0,6215 0,6418 0,6546 0,6619 0,6653 0,6665 0.6667 0,0531 0,1583 0,2579 0,3490 0,4290 0,4961 0,5504 0,5917 0,9973 0,9755 0,9330 0,8716 0,7939 0,7034 0,6040 О,а000 1.0182 0,9941 0,9469 0,8788 0,7924 0,6917 0,5833 0,4761 5571,92 < 1 < э715,92 0,7 0,5 0,3 0,1 6 18 30 42 54 66 х Рис.
14. Траектория жидкой частицы Второй пример относится к значениям параметров Д = 104, о = = 1. При о = 1 бифуркация Хоцфа докритическая, но есть надежда, что в фазовом пространстве можно обнаружить устойчивое периодическое течение, лежащее далеко от течения Пуазейля (см. рнс. 11). Расчет велся при п = 15, ш = 30. Отправляясь от хаотических начальных данных (на рис. 12 приведены линии уровня функции ф(я, у, 0): масштабы по осям координат выбраны различными, и значения величин уровня не указаны, чтобы не загромождать рисунок), расчет велся довольно большое время (О, 1,(, 1, = 6 10з, Естественная едняица времени для данной задачи равна 2.~г/0,2375 = 26,455.
Заметим, что 0,2375 = 1шЛг, где Л1 собственное значение линеаризированной задачи, имеющее хгвксиыальную вещественную часть. Вначале течение быстро приближалось к пуазейлевскому течению, а затем очень медленно от него уходило, оставаясь очень гладким по пространственным координатам, причем не наблюдалась какая-либо периодичность либо регулярность. <1тобы проиллюстрировать сказанное, на рис. 13, 14 представлена траектория жидкой частицы на довольно больших временных промежутках. Создается впечатление, что при выбранных параметрах мы натолкнулись на аттрактор, движение по которому апериодично. Заключение Численный анализ столь бурно развивается и стал столь обширной теорией, а точнее совокупностью теорий, что невозможно в одной книге охватить все его аспекты.