Вычислительная теплопередача Самарский А.А. Вабищевич П.Н. (1013606), страница 86
Текст из файла (страница 86)
С учетом граничных условий (19), (22) имеем уз — 4уз+ 5у~ йх 1=1, 9.5. Переменные ефункция тока, вихрь скорости, температура» 473 Это соответствует записи разностного уравнения (29) и в граничных узлах сетки с учетом (19), (20). 9.5.3. Аппроксимация по пространству в задаче конвенции После рассмотрения одномерного примера нам можно перейти к формулировке дифференциально-разностной задачи, соответствующей аппроксимации по пространству задачи конвекции (3)-(5), (8) — (12). Пусть теперь П вЂ” прямоугольник, а й — равномерная сетка в нем.
На множестве сеточных функций П, обращающихся в нуль на границе сетки, определим сеточный оператор Лапласа обычным образом: 2 Лу = ~~! Л~у, Л~у = -у;;,„а = 1,2. (34) а=! Для конвективных слагаемых используем аппроксимации второго поряд- ка на основе дивергентной записи (6), (7). В безындексных обозначениях получим !г — Лф + У(т)Лф + и(Лгф + р(к)ф) — 13и; = О, кЕы, 0<1~(Т. Здесь сеточная функция р(к) определена в соответствии с (31), (32) и имеет вид (37) Р(к) — ~ Ра(ки), а=! (38) О, Ра(ка) = 2 Ь4! Ла ( аа С 1а Ла~ о=1,2.
1е Ла~ ка — Ла~ г У(т) = Ч~; У.(е.), а ! У!(е!)го (Фадт)аа У2(е2)и! (фа,и!)я ~ причем на основании условий прилипания и непротекания компоненты скоростей, вычисляемые по центральным разностям функции тока ф, равны нулю на дш. При выборе (35) имеет место разностиое тождество (К(т)и!, ф) = О, (36) так как функция тока обращается в нуль на ди!. Уравнению (13) с учетом однородных граничных условий (8), (9) для функции тока ставится в соответствие дифференциально-разностное уравнение 474 Егава 9. Конвекн!нвный нгенлообмен Здесь учтено, что при аппроксимации конвективных слагаемых в соответствии с (35) значение вихря на границе несущественно. На расширенной сетке для разностного бигармонического оператора Лз имеем з Ле г7гхгхгвгхг + 2Рхгхгвгхг + Рхгхгвгхгг У(т) — ~, У (е ) а=! 1/ да д У,(е )а = — ~е,— + — (е,а) 2 ), д*.
дан (39) а = 1, 2, и поэтому (40) Аппроксимация (39), (40) дает У(т) — ~~г У,(е ), агп (4 1) где 1 У!(е!)а = -(г)га а;, + (г)гз,а)а,), (42) Уз(ез)а = — -(Фз,аа, + (Фз,а)й,). При такой аппроксимации конвективных слагаемых оператор Р(т) является кососимметричным и поэтому (ср. с (36)) (х'(т)а, а) = О. (43) С учетом (41) — (43) и граничных условий (11) дифференциальноразностное уравнение для температуры имеет вид оа — + Р(т)а+ нЛа = уг(н, 1), ю Е нг, 0 < 1 ( Т. (44) Правая часть уравнения (44) поправляется в приграничных узлах с учетом неолноролности граничных условий (11) (см.
п.4.7). При аппроксимации уравнения теплопроводности (5) используется обычный (вместо (6), (7)) выбор конвективных слагаемых в виде полусуммы конвективиого переноса в недивергентной и дивергентной формах. Положим 9.5. Переменные «функцин тока, вихрь скорости, теанература» 475 (47) 9.5.4. Рааностные схемы Приведем некоторые разностные схемы для сформулированной дифференциально-разностной задачи (37), (44) — (46). Естественно использовать чисто неявную разностную схему, которая имеет вид Лань! Лфн 2 +Р(тны)Лф„ь(+и(Л ф„+~+р(х)ф„+~) —,0(а„+~)4,=0, (50) Уравнения (37), (44) дополняются начальными условиями ф(х, О) = ф (х), х б ы, (45) и(х, О) = и (х), х б ы, (46) которые соответствуют (10) „(12). Приведем простейшие оценки решений дифференциально-разностной задачи (37), (44) — (46), которые согласуются с соответствующими оценками (см.
п. 9.4) решений исходной дифференциальной задачи. Определим сеточное гильбертово пространство Нн, в котором скалярное произведение определяется выражением (и, а)а = (Ли, а), а Ц Ц вЂ” норма в Н. Домножим скалярно в Н уравнение (37) на ф, а уравнение (44)— на а. С учетом (36) следует ((г(т)Лф, ф) = О. Принимая во внимание однородность граничных условий, имеем 1 (ие, ф) = — (и, ф; ) = — -(а, ф„+ фе,) < ЦиЦ ЦфЦн. Поэтому из уравнения (37) получим 4 г(1 — ЦфЦ, < ДЦаЦ. Из уравнения (44) с учетом (43) будем иметь — ЦиЦ ( Цф~. 4 (48) Неравенства (47), (48) полностью согласуются с неравенствами (15) из п.9.4 для решения дифференциальной задачи.
Из (47), (48) непосредственно вытекают следующие простейшие априорные оценки: Ца(х,з)Ц < Ци (х)Ц+ Т шах Ц(о(х,з)Ц, Фе(О,Т( (49) Цгр(х, 1)Цн ( Цф'(х)Цн+,0ТЦи~(х)Ц+ ВТ' шах Ц(о(х,з)Ц. ме(ь,т( При исследовании разностных схем мы будем ориентироваться на получения подобных априорных оценок. 476 Пгава 9. Конеектиеный тенлообмен (53) (56) у(х,г) = — ~~1 Ьлр (х )Р(х,е). (57) и»+1 - и» + ч (ч„+1)а„+1 + кйа,ц.1 = 1»», (51) хб1о, п=0,1.... В силу нелинейности при реализации разностной схемы (50), (51) используется тот или иной итерационный процесс.
Можно применять и простейшую линеаризованную схему йлр»1-1 - й171» 3 + г'(ч„)йгрлч1 + о(Л 3(1»ч1 + р(х)лрл+1) -)у(а»11)Е, =0> (52) а»+1 ал + У(Ч»)а»1.1+ Кйи»Ч1 = 'р т х 61о, а=0,1 когда конвективный перенос осуществляется по полю скоростей с пре- дыдущего временного слоя.
Но даже в этом случае для нахождения Ф„„.1 необходимо обращать сеточный эллиптический оператор четвертого порядка. Для разностной схемы (50), (51) нетрудно получить соответствующие оценки устойчивости. Для этого скалярио домножим уравнение (50) на лр„+1, а уравнение (51) соответственно на а»+,. По аналогии с (47), (48) получим 11111р»+11111л ( 11111Е4л+ тйи 1.1~1, !!и„+1~~ < ~~а !1+ тЬЛ. Дальнейшее упрощения удобно пояснить, возвращаясь от уравне- ния (52) к системе уравнений для сеточной функции тока и сеточному вихрю скорости. Пусть йлрл+1 =го»+1, х Еы а=О 1, (54) тогда уравнение (52) примет вид ГО»Е1 ГО» т + 1г(чл) 1о»11 + о(йи1»~1+ р(х)гр„+1) - р(а +1) Š— — О, (55) хЕы, п=0,1,.... В уравнении для вихря (55) член р(х)лр»+1 обусловлен тем, что граничное условие для вихря не задано.
Как и в одномерном примере (см. (27), (33)), операторное уравне- ние (55) можно записать в виде обычного уравнения и1»+1 гол + ч (чл)и1»+1+ ойи1»1-1 13(а»ч1)е, = 0> хны, п=0,1,... с некоторыми неоднородными граничными условиями. Определим функ- цию 9.5. Переменные «функция тока, вихрь скорости, температура«477 Тогда граничные условия для уравнения (5б) можно записать в в1ше У»+1(х1+ "1, х2), У»+1(х1 «21~ х2)~ д «1 (хи х, + л2), Ул«1(х1 х2 222) х1=0, х2 =О, х2 12ю (58) 2О»+1(Х) «л так что, например, 2 гол+1(О, х,) = — — 2Ф(Ь1, х,), О < х2 < 12. Ь1 2 дл,(х) = — ~~1 ь~ Р„(х )фл(х).
»=1 Это соответствует определению функции тока из уравнения (54), а вихря из и2»+1 го» + ~(тл)~~~ + (™~~ + р(~)Ф~) Ф(~~~ )Е, = 01 хбь2, п=0,1,.... Схема (54), (59) очень удобна для реализации (решаются два эллиптических уравнения второго порядка), но она не является абсолютно устойчивой. Ее условная устойчивость обусловлена именно граничными условиями.
Второй пуп построения безусловно устойчивых разностных схем связан с итерационным уточнением граничного условия при реализации схемы (54), (55) ((54), (5б) — (58)). Пусть р — итерационный параметр, тогда при итерационном уточнении граничного условия на кюкдой итерации решается система уравнений »+1 л+1 л + р(т') О 1+ т + и(Лш„"', + р(х)(рф„"«1+ (1 — р)«р„+11)) — )5(и„+1);, = О, »+1 »«1 йф„1 = юл«„х б ь2, и = О, 1,.... Условия (57), (58) есть не что иное, как хорошо известные граничные условия Тома для вихря скорости. Еще раз подчеркнем, что уравнение (55) эквивалентно уравнению (56) с граничными условиями (57), (58).
Это есть просю операторная запись (5б) — (58), переход от задачи с неоднородными граничными с условиями к задаче с однородными граничными условиями. На основе (54), (55) можно строить разностные схемы более удобные для вычислительной реализации. Простейший подход связан с тем, что граничные условия для вихря берутся с предыдущего временного слоя, т.е.
в (58) 478 Егава 9, Конеекгнненый теплообмен Понятное дело, что такой подход ведет к существенному увеличению вычислительной работы. Кроме того, встает проблема выбора итерационного параметра,и. 9.б.б. Беаытерациоииая реалиаация граничных условий для вихря скорости Построим безусловно устойчивые разностные схемы в переменных «функция тока, вихрь скорости», которые не связаны с итерационным уточнением граничного условия.
Для зтого будем исходнп из дифференциально-разностной задачи (37), (44) — (46). Уравнение (37) запишем в виде Н вЂ” ЛлР+ А ~Л3а+ Азй1а =,Оати х Е ы, 0 < Е < 2". (60) где А~ = У(ч) + иЛ, Аз = ир(х)й. (61) Тем самым, выделение граничных условий для вихря осушествляется оператором Ат. Для уравнения (60), (61) можно использовать те или иные аддитивные разностные схемы. Например, схема суммарной аппроксимации будет иметь вид йрл+цз-ЛРл + р(ч„+Ч~з)Л4„~0ц+ий Р„лез-— Яаны)аи (62) +ир(х)$„+1=0, хны, п=0,1,.... (63) т Полагая ЛФп+лД глп+о1п о 1~ 2г первый полушаг разностной схемы (62), (63) можно рассматривать как уравнение для вихря скорости: лил+1/3 чел + Р(ч„+цз)ле„~11т+ ийы„лцг = 13(а„ы)ти (64) а второй полушаг (63) — как уравнение для функции тока: йгр„+1 + тир(х)лр„+, — — юл+,р, х Е ы, и = 0,1,....
(65) Заметим, что уравнение (64) соответствует определению вихря скорости при однородном граничном условии первого рода. Учет неоднородности проводится в уравнении (65), которое отличается от обычного уравнения для функции тока (54) в приграничных узлах. Для того, чтобы показать абсолютную устойчивость схемы (62), (63), достаточно домножнть скалярно уравнение (62) на Ф„~цм а уравнение (63) — на Р„„.~.
Непосредственно следуют оценки 113рлмуг1!л < Ь 11л+ тР11ан ы!1 !!аде м11л < Ьл+цз11» т. е. !1чр 111л < 11лр 11л+ то!!а М. 9.5. Переменные »функция я!ока, енлрь скороен/и, н/анкерна/ура» 479 йф„,ц, — лф„ + 7/(тп+!/2)йфп+2/2+ + Рй 'фп!.!/2 + н~о(к)фп = /3(и~»!)» (66) + л (тп»!/2)йфл»1/2+ + ий~фп»!/2+ ир(к)фп»! =/3(ял»!)», (67) 0,5т йфл+1 йфп+!/2 0,5т Устойчивость разностной схемы (66), (67) исследуется так же, как устойчивость схемы переменных направлений лля уравнения теплопроводности (см.