Fletcher-1-rus (1185917), страница 39
Текст из файла (страница 39)
Распечатка подпрограммы ЗАСОВ. .азото-ох а- - хгхао-ог .59710-09 -.хоота-оз .25550-01 а= .33440-01 .77аао-ог .49930-аа .ззаао-оа а= .173Во-02 .14490-02 .25550-03 .44470-04 а .29250-04 .гааги-04 .44790-05 .4447о-от а- 45звв-оа .55700-оа .ззвво-ов 223 5 6.1. Нелинейные стационарные задачи б.!.3. (с)Е%ТВИ: двумерньсе нестационарные уравнения Бюргерса Нестационарная форма вышеупомянутых уравнений более детально рассмотрена в $ 10.4. Здесь же нам необходимо решить стационарные двумерные уравнения Вюргерса дй дй 1 сдай дкй т й — +б — — — ~ — *, + — *,-)=О, дк ду ке адик ду ) дб дб ! /дкб дтб Х й — +б — — — ~~ — *+ — *) =0 дх ду йе ч дхк дук ) (6.12) где и!+с « — и- к кин « иС с « — 2а «+и! Бккиь « (6.14) Биа Решение ищется в области, изображенной на рис. 6.9.
Граничные значения взяты из точного решения ]Р!е1с)сег, 1983] й= — — — б= — —— 2 дФ 2 дФ (6.15) Кеф дх ' йеФ ду ' где ф=а, + пах+ а,у+а,ху+аз(е"ск-к с+ е «ск к >соз)(ссу), (6.16) тогда как постоянные ас, аз, аз, ас, аз, )с и х, подбираются так, чтобы разнообразить поведение точного решения. Если реализовать представление уравнений (6.13) во всех внутренних узлах, то получится 2(сссХ вЂ” 2)(ст'У вЂ” 2) уравнений, котоРые должны быть РазРешены относительно иь «и бь «длЯ всех внутренних узлов.
Процесс решения с использованием с граничными условиями Дирихле для и и б, применяя с этой целью метод Ньютона после дискретизации с помощью трехточечных центрированных разностных формул. Уравнение (6.12) можно записать в следующей дискретной форме: 1 )ти! « —— и; «Б,и! «+ б! «1.„и; « — — (Баки; «+ 1.уии; «) =О, (6.13) ВУи«=ин«с би«+ бь«(.абн« вЂ” Ж(Ь бйа+ Ьаубй«)=О Гл. 6, Стационарные задачи метода Ньютона может быть символически представлен формулой 1~"~ЛЧ~"~ ~ = (6.17) где составляющими вектора ЛЧ являются Лиь «, Лоь «, а вектор К имеет компоненты Яиь «и Юл м оцениваемые в каждом из и=й(х,у „1, очй(х,у 1 М Ушел и:=й(1,о) о= б(-1,у) ичй(1,8) о=ой(1,у) ц=й(х,о), очй(х,о) х=1 Рнс.
6.9. Вычислительная область для решения урааненнй (6.12). внутренних узлов. Верхний индекс (п) используется для обозначения номера итерации. Якобиан ) содержит такие члены: давид «/доил « — — — 0.5и; «/Лх — 1/(Ке Лх'), дтти /диг « — — (2/Ке) (1/Лха+ 1/Лу«) + (0.5/Лх) (и1+~ « — и1 ~ «), дй1и, «/диг+, „— — 0.5иг «/Лх — 1/(йе Лх'), дйи1 «/до1 «, = — 0.5оГ «/Лу — 1/(йеЛу'), (6.18) дйиь «/до1 „, — — 0.5ог«/Лу — 1/(Ке Лу'), данг «/до1 « —— (0.5/Лу)(иг «„, — и1 «,). Аналогичные члены можно записать для )1о;, «. Вышеприведенная формулировка задачи и ее решение реализуются в программе )ЧЕ%ТВБ. Распечатки этой программы, а также подпрограмм ЕХВ(1К для построения точного решения, КЕБВ(3 для расчета иевязок согласно (6.!3) и )АСВ1) для вычислений по формулам (6.18) приводятся на рисунках 6.10 — 6.13.
В табл. 6.2 расшифровывается смысл параметров, используемых в программе МЕЮТВ(/. 55 56 57 5В 59 60 61 С 62 63 64 65 66 67 ба с 49 ТО '71 72 7З '74 75 С 76 С пс 7В 79 С ВО 01 02 .Вз 54 В5 06 .В7 С ВВ С 09 с ТО с '91 92 5З 94 94 С 97 С ва с 99 100 С 101 С 102 С 103 104 С 105 ЗО6 107 103 С 109 С 110 С 111 112 С ЗЗЗ С 114 С В Ео )О д = З,ИХ во 9 е 1,ит 0(н,д) ОК(К,д) т(к,д) тк(к.д) 9 соитпшк 10 соит1иуе 11 ВО 12 К 1,МТ 12 МВ1ТЕ(6,13) (ОЕ(н,д),д 1,МХ) 13 Роанат(( Ое,7110.4) 00 14 К 1,ИТ 14 мазте(6 ° 15) (те(к,д),д з,мх) 15 Роанат(' ТЕ ',7Р10.4) аи и Хт о 1ТР 0 онн он втз З.вт 16 СОМТ1 МОЕ саьсоьатк акззооаьз саьь акзво(о,т,а) ЗОИ О.
ВО 17 1 з,н Вэ(1) а(и 17 зон зои + а(1)*а(1) 1В Внз Рзяат(33луам) 1Р(1Т ОЕ. 1ТР)МВ1ТК(6,19)1Т,ВНЗ,(В(д),ди1,3) 19 Роанат(' 1Т=',13,' ВНЗ ',011.4,' В ',3011.4) ткзт Роа соиткаокиск оа нипшк шшвка ОР Хтеавт10из ехсеевев Хтознз .ьт. МРЗ)сото Зз 1Р ПТ .Ок. 1ТР)1ТР 1ТР + 1РМ 1Т ЗТ+ 1 1РПТ .ЕО. ХТНХ)0070 23 1Р[ВНЗ .СТ. 1.0Е+03)СОТО ЗО саьсоьате дасОВ1аи саьь дасво(о,т,ад) РвстОВ1зе тне дасов1ам 1мто Ь.н саьь Расто(,ад,дРуть ЗР1дРТТ(н) .Ео; -1)иа1ТВ(6.20) 20 Роанхт(' Хеао Р1тот ветесттв') 11(дРТТ(н) .ЕЯ.
-1)ООТО 30 зоьте РОВ тие соааест1ОИ, вт саь зоьте(я,ад,дРтт,ав) П(СВЕНЕМТ 0 еив Т Рис. 6.10 [продолжеиие). 6 6.1. Нелинейные стационарные задачи РО 22 д 2,КХР ВО 21 К 2,ИТР ИВ 2*(К 2] 4 2*[МУ 2) ° (3-2) РО " — ЯВ(ИВ+1] ВЧ = - ЯР(ИВ+2) 0(К,3) 0[К,д) + РВ*ОК Ч(к,д) Ч[к,д) + РЧ4ОК СОМТ1ИОЕ СОКТ1ИОЕ СОТО 16 115 116 117 118 119 120 121 122 21 123 22 324 125 С 126 С 127 С 128 23 129 24 130 131 1Э2 25 133 26 1Э4 135 27 136 28 137 С 138 139 140 141 14г 143 144 145 ОЕИЕЯАТЕ ООТРОТ КЯ1ТЕ(6,24)1Т,ЯКВ ГОЯКАТ(//,' АГТЕЯ ',13,' 1ТЕЯАТ10КВ ТВЕ ЯИВ ЯЕ$100АЬ 1$ 1И2.
5) ВО 25 К 1,МУ 881ТЕ[6 26) [0(К З) 3=1 ИХ) ГОЯИАТ(' 0 ',7Г10.4) ВО 27 Е 1,КУ ЪЯ1ТЕ(6,28] (Ч(К,З),3 1,МХ) ГОЯИАТ(' Ч ',7Г10.4] 1Г(1ЯР .Ей. 0)СОТО 30 ВО 29 К 1,ИУ УЯ1ТЕ[7,6)(0(К,З],3=1,МХ) ъя1те[7,6)(ч(к,з],а 1,их) 29 СОМТ1ИОЕ 30 СОИТ1МОЕ ВТОР ЕКР Рис. 6.10 (окончание). Выражения (6.13) вводятся в вектор 14 в следующем порядке: тта2,2 4ЪО2, 2 ттаз, 3 4302,3 ° ° ° тта2, ИУ-! 4ЪО2, УУ вЂ” ! 4таз, 2 а порядок введения скоростей в вектор [) таков: 2,2 02.2 а2,3 02,3 а2,МУ вЂ” ! 02,МУ-! М3,2~ 03,2 15" В качестве начального решения относительно й и й было использовано точное решение уравнений (6.12) со значениями параметров а! = а, = 110.13, а, = а4 — — О, а, = 1.О, А = 5, )[е = 10, у „ = и/30 и хо = 1.О.
Однако с помощью программы )ъ[ЕЖТВ[1, примененной на сетке 5)4, 5, т. е, с 18 узловыми неизвестными значениями и и о, не удается получить сходя!пееся решение. Для получения такого решения необходимо Гл. 6. Стационарные задачи 228 Таблица 6.2. Параметры, используемые в программе ХЕ%ТВ() Описание Параметр Число узловых точек в направлениях х и у Максимальное число итераций =О, 1, вводит начальное решение, полученное с помощью подпрограммы ЕХВ())1 =1, 2, записывает окончательное решение в файл ХЕй)ТВ().ЬТА =2, считывает начальное решение из файла ХЕ%ТВ().ЬТА Выдает на печать невязки через каждые 1РХ итераций Максимально допустимая среднеквадратичная невязка уравнений Число Рейнольдса ке, формула (6.13) Показатель ниакней релаксации о(, формула (6.19) Шаг по времени, используемый в псевдонестационарном варианте ($6.4) Число уравнений, подлежащих решению Точное решение й, о аь ..., а, в формуле (6.16) Зависимые переменные в уравнениях (6.13) Навивки уравнений (6.13) =)1 после возврата из подпрограммы кЕЗВ(); содер)кит поправки ЛЧ после возврата из подпрограммы ЗО(АгЕ Якобиан Л в уравнении (6.17) Поправки для и н о; Ьй в уравнении (6.!7) ХХ, ХТ 1ТМХ 1лР 1РХ ЕРЗ кЕ ОМ РТ Х ()Е, т7Е А 1), Н и )(Р АЛ Р(), РТ) подвергнуть нижней релаксации поправку к Ч(л).
Это означает, что применяется следующая формула: (1(л+1) (1(л) + ш Дч(ле1) (6.19) при ш = 0.15. Однако и в этом случае скорость сходимости очень мала (см. данные на рис. 6.14). Такая ситуация не может считаться необычной при применении метода Ньютона к дискретизированным уравнениям, описывающим задачи гидроаэродинамики. А именно поправки, определяемые по методу Ньютона, могут оказаться достаточно большими для того, чтобы текущее решение вышло за пределы радиуса сходимости (см. п. 6.1.1), если только эти поправки вводятся без изменений.
Вид формулы (6.19) подсказывает соответствующую стратегию. Она состоит в том, что величина ЛЧ(л+'), получаемая в результате решения уравнения (6.17), определяет лишь на- правление поиска. Для трех значений о) строятся новые реше.ния Ч'('и=я["]+а йЧ("'и ж вычисляются соответствуюшие им невязки ясй+~~. Для каж.дого значения О) рассчитывается среднеквадратичная невяз- Зовяонтэик кхвоЙ(нк,те,А.Аь] ка ф +(] причем предполагается наличие квадратичной зависимости (р,с „, от оз.
Затем выбирается такое значение (о, (и<-() а следовательно, и Л(]("('>, которое ведет к минимальному эт'„, для заданного направления поиска (з(](аы). Эта стратегия, однако, может оказаться слишком доростоящей в вычислительном смысле, если невязки эсп и зсэг имеют сложные выражения.
1 2 эс бс 5С 6 С 7 1 9 1О 11 12 13 14 15 16 17 1В 19 20 21 22 23 24 25 26 27 2В 29 зо Зз 32 эз 34 35 6 6.1. Нелинейные стационарные задачи СЯЬСОЬАТЕЗ ТИЕ ЕХАСТ ЗОЬОТэои ОР ТКЕ ТМО-В1ИЕИ$1ОНАЬ ВОЙОЕЙЗ' кооятэоиЗ озтио тнк соьк-ноРР ТЙАИЗРОЙиятэои В1НЕИ$10И Я(5),ОЕ(21,21],ТЕ(21,21),РВ(21,21) СОНМОИ ВХ,ВТ,ЙЕ,ИХ,М'Г Р1 "- 3.1415927 Хг 1.О ТИАХ = Р1/6./АЬ и[т ит-1 Вт тмАхгямт яих = их " ВХ = 2.7ЯИХ ВО 2 д 1,ИХ Ад д-1 Х = -1.
+ Яд*ОХ ХВ Х-ХЕ ВЕХ = ЕХР(ЯЬ*ХВ) + ЕХР(-АЬ ХВ) . ВВХ - "КХР[АЬ~ХВ) - ЕХР(-АЬ*ХЮ ВО 1 К = 1,МТ АХ К-1 Т АК*ВТ ЗТ = 51М(АЬ*Т) СТ вЂ” С05 ЬЯЬтТ] РЙ(х,д! А(1] + А(2]*Х + А(3]*Т + Я(4]*Х*Т + Я[5]"ВЕХ*СТ РМХ А(2) + А(4] Т + А(5)тАЬ*ВВХ*СТ Рат А[3] + А[4]"Х " А[5]*АЬ*ВКХ*$7 Ук(к,д] = — 2.*РНХ/РИ(х,д]/ЙЕ ЧК(х,д) = - 2.*РНТ/РН(х,д]ИЕ 1 СОИТ1ИОЕ 2 СОМТ1ИОЕ ЙетоЙМ КМВ Рнс. 6.11. Распечатка подпрограммы ЕХВ[)й.
Гп. 6. Стационарные задачи 230 В качестве альтернативной и более эффективной стратегии по сравнению с применением формулы (6.19) можно указать на сочетание метода Ньютона с псевдонестационарной форму- Ъ 2 вввкооттие ккввщв,удм Э С 4 с еувьщтев кевтввзьв от 20 втезву впвоекз' еоваттоик 5' С 6 7 в 9 С. 10 ИХР ИХ - 1 11 ИТР ИТ" 1 12 СХ 0.5/0Х 13 СУ м 0.5/ЗУ 14 ССХ 1./КЕ/вк/ВХ 15 ССТ 1./КЕ/вт/ВТ 16 ИСТ 1 17 С 10 12 20 21 22 23 24 25 26 27 ЗВ ЗЭ 30 31 32 ЗЭ С 34 КЕТ0КМ 35 еив ведь.в к,в.т,сх,ст,ссх,ест,вви 01ИЕИЗ10И 0(21,21),У(21,21),К(50) соммои вх,вт,кк,их,ит ВО 2 Ю 2,ИХР Л< 3-1 ЗР =Э+ 1 во 1 к = З,итт КМ = К - 1 КР К+1 Ввк СХ*0(К,Ла(0(К,ЗР)-0(К,ЗМ)) + СУ*У(К,З) (0(КР,З)-0(КМ,З)5 к(ист) = ввм — ссх*(0(к.зи)-2.*0(к,з)+0(Х,ЗР)) - сей*(0(км,з) 1 - 2.*0(К,З) + 0(КР,З)) Ввк СХ*0(К, Т) *(т(К,ЗР) -У(К,ЗИ) ) + от*У(К,З) а (У(КР,З)"У(КМ, Т)) ки(ст+ы - вни - ссх.(у(к.зм)-2.