Вычислительная теплопередача Самарский А.А. Вабищевич П.Н. (1013606), страница 104
Текст из файла (страница 104)
1МР1.1С1Т КЕА1 '8 ( А — Н„О-Е ) 01МЕХ$1ОХ ОМО(Х1,Х2), Р31(Х1,Х2) СОММОХ / Т04 / Н131, Н231, Н1412, Н2412, НН81, ТА1Л, РК1, ОК2Н!, Х1, Х2 313М = Н131 + Н1Я + Н231 + Н2Я 707 13.4. Конвенция в нояасгни квадрагннага сечения С С Внутреввве узлы С !ХО 20 Х = 2, Х2 — 1 00 101=2, Х1 — 1 ОМО(!,Х) = Б!ХМ'РБ!(1,Х) — Н1Б!е( РБ!(1 — 1,1) + РБ1(1+1,Х) ) — Н2Б1'( РЯ(1,Х-1) + РБ1(1„!+1) ) 1О СОХТ1ХЮЕ 20 СОХТ1ХУЕ С С Ннжвяя н верхняя гравнцы, краевые условия Тома С !ХО 30 1 = 2, Х1-1 ОМО(1,1) = -2.00'Н2Б1еРБ!(1,2) ОМО(!,Х2) = -2.ОО"Н2Я"РЯ(1,Х2-1) 30 СОХТ!Х!ХЕ С С Левая в правая границы, краевые условия Тома С ХХО 40 Х = 2, Х2-1 ОМО(1,1) = -2.!ХО"Н1Б!'РБ1(2,Х) ОМО(Х1,Х) = -2.00'Н1Б!'РБ!(Х1-1,1) 40 СОХТ1Х!ХЕ С С Углы расчетной областн С ОМО(1,1) = 0.130 ОМО(Х1,1) = 0.00 ОМО(1,Х2) = 0.130 ОМО(Х1,Х2) = 000 13.4.4. Решения сеточных зллнптичеекнх задач Дадим описание используемого итерационного метода решения несамосопряженных эллиптических задач, который используется при расчете температуры н функции тока (см.
(24) и (25)). Для реализации этапа расчета вихря скорости (уравнение (23)) используется подпрограмма Б01УЕ1„подробно описанная в и. 13.1. 708 1Ъава 13. Примеры чиелеииого моделировоиил Рассмотрим систему линейных уравнений с невырожденной матрицей Ау = 7, (27) где А = Ь + В + 17, А ~ А* (17 ~ Ь'), (28) В попеременно-треугольном итерационном методе для приближенного решения задачи (27), (28) выберем В по аналогии с самосопряженным случаем (см. п. 13.1) в виде В = (С+А)С (С+И), (29) где С вЂ” некоторая диагональная матрица с положительными элементами. Итерационный метод запишем в виде: г» = 7 — Ау», в» =В г», во, (30) Р» = ⻠— 73»Р»-и У»+1 = У»+ а»Р».
Из последнего соотношения и определения невязки г» следует формула г»+~ — — 㻠— а»АР», (31) Й= О, й=1,2,..., которая используется для пересчета невязки на новом итерационном шаге. Приведем расчетные формулы для определенна итерационных па- раметров. Пусть С вЂ” самосопряженный положительный оператор, за- дающий энергетическое пространство Н, со скалярным произведением (у, е), = (Су,е) и нормой 8У8~ = (у, у),. В обобшенном итерацион- ном методе минимальных невязок итерационный параметр а» в (30) выбирается из условия минимизации Лг»»Д„что дает (гы Ар») (32) (Ары Ар») Итерационный параметр Д определим по аналогии с самосопряженным случаем из условия ортогональности векторов Ар» и Ар», в Н;.
б»= (Ав», Ар»- ~), (33) (Ар» нАР» 1), В рассматриваемом несамосопряжениом случае (28) подчиним выбор С требованию применимости алгоритма, родственного описанному в п. 13.1, когда итерационный шаг фактически сводится к решению двух систем уравнений с треугольными матрицами и вычислению трех скалярных произведений. С этой целью положим С=((С+ЮГ')'С(С+А) '. (34) 709 13А, 1аонвекиия в новости квадратного сечения чоа — — (б+1/) 'б(б+Ь) 'га, ра = аеа — Дра-ь (б1/2(б+ Ц-! Аща б1/2(б+ Ц-|Ара,) (б1/4(б+ Ц ~Ара ь б1/з(б+ Ц ~Ара 1) ' у,+, -- у, + аарь га 1 = га — ааАр», (б1/2(б+ Ц,г б1/2(б+ г)-1Ара,) (б'/'(б + й)-1Ар„б'/а(б + Т)-'Ара) (35) Дополним соотношения (35) формулой для пересчета вектора Ара.
Ара = Ааеа — ДАра ь (36) Умножим соотношения для аоа, ра и уач.1 слева на б '/ (б + 1/), а соотношения для Ара и га,.1 слева на б'/з(б + Ь) ' и введем обозначения: г = б~/~(б + Ь) ~га = б ~/~(б + 1/)чоа, ра=б '~'(б+1/)рь у,=б '"(б+1/)у„ да = б / (б + Ь) 'Ар».
Вместо чоа вычисляется 1а = б'/г(б+ Ц 'Ааоа. С учетом введенных обозначений итерационный процесс (35), (36) реализуется следуюшим образом: та = б'/ (б+ Ц 'А(б+ У) 'б'/ га, ра = та — Дра ь (ала ) Ь = 1а — /унда-1> /5а =, (Да-и Ча-1) уа+~ = уа+ аара, ('ьча) ге+1 = га — аада, аа = = (Да ча) Так же как и в 13.1 используется масштабирование с диагональными матрицами б'/' и б '/~ и аналогичный прием (см. (30) в п. 13.1) лля вычисления 1».
Ниже приведено описание алгоритма обобщенного метода минимальных невязок — приближенной факторизации. В соответствии с (29)-(33) итерационный шаг реализуется в рамках следуюших вычислений: 71О Бгава 13. Примеры численного моделирования 1. Вычисление зяементов диагональной матрицы П ~/з. В подпрограмме В01УЕ4 реализован вариант поперемеиио-треугольного итерационного метода — метод приближенной факторизации (см. п.4.7). 2. Переход от задачи (17), (18) к задаче (21), (22). ~-~/з27П-~/з 2Н ~-1/2~П-!/2 -ь/з -из 3.
Начальное нриблизкение. Уь:= гг У~ у ии Пуь. 4. Начальная невязка. У:=Уз = О '(П ' 'У вЂ” 27уь — у) — уь, ггО гы (У,,/), тт:= ггб, ааь':= О, игц:= О. 5. Решение системы уравнений с верхне-треугольной матрицей. б. Решение системы уравнений с ниъене-треугольной матрицей. г:=г+м '(/+27г), га:= (ю, д), Ьй:= Гд ь дай 7. Пересчет итерационных векторов а:= à — Ьй ь а, р:= / - Ьй ь о, и:= (й й) гй:= ( И да1:= 1/аа, ай:= га ь аой 13.4. 1(онвекция в волости квадратного сечения 711 р:= у + а!с з р, 1 != у — айза, гг:= (г,г), пВ:= пП+ 1, Ниже приведен текст подпрограммы 801,тЕ4, реализующей описанный алгоритм итерационного попеременно-треугольного метода— сопряженных градиентов в варианте приближенной факторизации лля решения несамосопряженных пятиточечных сеточных задач. Р— С С С С С С С С С С С С С С С С С С С С С С С С 8.
Лересчет итерационного нриблизкения и невязки. лри необтдимости продолжить итерации бо Го 5. 9. Обратное масштабирование (искомое приближенное решение). $1)ВКОЦТ1ХЕ $0(УЕ4 ( Х, 1., А, г', Р, ЕРБ ) 1МР1.1С1Т КЕА1'8 ( А — Н, 0 — Е ) Р1МЕХВ1ОХ А(1), У(Х), Р(Х) ПАРАМЕТРЫ: ИН.чЕ4 — Подпрограмма решения системы линейных уравнений патиточечной ревностной схемы для несамосопряженного двумерного эллиптического уравнения второго порядка с переменкыми коэффициентами в прямоугольнике; реализован попеременно-треугольный метод приближенной факторизации — сопряженных градиентов. порядок системы (число узлов сеточного прямоугольника); полуширина ленты пятидиагональной матрицы системы (число узлов сеточного прямоугольника по первому направлению); массив, в котором задаются элементы несимметричной пятидиагональной матрицы пятью диагоналямя.
Вслед эа коэффициентами в массиве А будут расположены векторы, используемые на итерациях; массив, содержащий ка входе начальное приближение и на выходе — итоговое итерационное приближение к решению; массив, содержащий на входе правую часть системы; 212 Г)ава 13. Примеры численного моделирования ЕР — параметр, задающий отвосителькую точкость итогового итерационного приближения к решению. Описание фрагментов массива А: Элементы главной диагоизли располагаются в компонентах с 1-й по Х-ю; ИЬ = Х ИК = 2'Х + 1 12Ь = ЗеХ 12К = 4"г1 + 1. 1О = 5"Х + Ь Элементы соседней верхней диагонали располагаются в компоиеитах с (11В+1)-й по 111В+г)-1)-ю; Элементы удаленной верхней диагоиали располагаются в компонентах с 112В+1)-й по 112В+Х-Ь)-ю; 1Т = 1О + Х 1ТЬ = 1Т + Ь а = ГГ + )ч) + Ь 1Р = 1О + Х Указатели 10, с7, 1Т1., Щ используются для обращения к Фрагментам массива А, расположенным вслед эа коэффяциеитами системы.
Реализация алгоритма ЯОЬУЕ4. 1. Вычисление элементов диагоизлького множителя матрицы приближеикой факторизации: РО 101= 1, Х О = А(1) — А(ИЬ+1)эА(1Т+1-1)е( А(ИВ+1-1) + А(12К+1 — 1) ) е — А(12Ь+1)еА(!Т+1-Ь)"( А(ИК+1-1.) + А(12К+1-Ь) ) А(1Т+1) = 1,РО / О А()О+1) = РБОКТ(А()Т+1)) 10 СОМГ1Х1)Е 21З 13.4. Конвекиия в новости квадратного сечения С 2. Масштабирование матрицы системы диагональной матрицей: С 330 20 Х = 1, Х А(Х) = А(10+3)еА(3)еА(10+3) — 2.330 А(!11.+3) = А(10+3)'А(!1Ь+3)'А(10+3 — 1) А(ПК+3) = А(КО+3) А(ПК+3) А(10+3+1) А(12Ь+3) = А(10+3)*А(121.+Х)'А(10+Х вЂ” Ь) А(12К+Х) = А(!О+Х)'А(!2К+Х)'А(10+3+1 ) 20 СОХТ1Х1!Е С С 3. Масштабирование начальною приближения С верхне-треугольной матрицей: С 0О 30 3 = Х, 1, — ! А(1Т1.+3) = У(Х) / А(10+3) т(х) = А(1ть+х) — А(ик+х) А(1ть+х+1) — А(12К+3)еА(1ТЬ+Х+Ь) 30 СОМПХЮЕ С С 4.