Fletcher-1-rus (1185917), страница 50
Текст из файла (страница 50)
Поэтому логическим усовершенствованием алгоритма (7.6) является схема Те+ Те — и (Ти 2Тп+ Т ) 2ДГ Дх (7.8) предложенная Ричардсоном (см. рис. 7.2). Однако, несмотря на то что схема имеет порядок О(Ь(з, Лхз), анализ ее устой- Как можно видеть, имеет место согласованность (см. $4.2) уравнений (7.5) и (7.1). Исходя из структуры главного члена в выражении для Еу схема ВВЦП будет считаться имеющей первый порядок точности по времени и второй порядок точности по пространству.
Однако следует помнить, что это утверждение, строго говоря, справедливо только в пределе Д(-з- О, Лх-~0. Практически решения строятся на сетке конечных размеров и величина таких членов, как дз7(д(з, заранее неизвестна. Как показывает анализ устойчивости по Нейману (2 4.3), коэффициент усиления 6 выражается формулой 287 $73. Явные методы чивости по Нейману 1Хоуе, 1983] показывает, что при з,') 0 схема безусловно неустойчива.
Таким образом, она не представляет практической ценности. Можно отметить, что замечание о неустойчивом поведении относится к уравнению в целом. Если аппроксимация с центральной разностью для производной по времени применяется к уравнению конвекцни (9.2), то может быть получен устойчивый алгоритм (9.15). Схему Ричардсона (7.8) можно модифицировать так, чтобы получить устойчивый алгоритм. Это достигается за счет замены Т~ в уравнении (7.8) на 0.5(Т~ '+Т~+ ).
В результате получается уравнение гл+! Тл-1 [гл (Тл-! 1 го+1) + гл Уравнение (7.9), известное как схема Дюфорта — Франкела (см. рис. 7.2), можно преобразовать так, чтобы получился явный алгоритм Т,"+'= ( — ) (Т";, +Т)+1)+ ( — ) Т~ (7 10) Схема Дюфорта — Франкела является трехслойной по времени во всех случаях, кроме з =0.5, а в последнем случае она совпадает со схемой ВВЦП. Прн использовании трехслойной схемы нужно хранить в памяти данные на двух временных слоях и возникает необходимость в альтернативной двухслойной схеме для реализации первого шага по времени.
Применение к уравнению (7.10) анализа устойчивости по Нейману (см. $4.3) дает коэффициент усиления 6, приведенный в табл. 7.1. Так как при любом значении 8 н при з,) 0 имеем 16( ( 1, то схема Дюфорта — Франкела устойчива прн любом значении Лй Однако имеется некая цена, которую нужно уплатить за этот весьма благоприятный результат, касающийся устойчивости.
Если в (7.9) подставить точное решение, разложенное в ряд Тейлора относительно (1, п)-го узла, то в результате получится %-аз„', +а( —;Л) Я+О(Л',Лх')=О. (711) Следовательно, для обеспечения согласованности необходимо, чтобы Л(/Лх-+.0, когда Л(-в.О н Лх — О, т. е. для наличия согласованности требуется выполнение условия Лг(( Лх. Однако а(Л1/Лх)т =зЛй а при решении задач о диффузии ожидается, что з имеет порядок О (1). Поэтому схема Дюфорта — Франкела согласуется с уравнением (7.1), но будет неточной, если Ги. 7. Одномерное уравнение диффуани величина зМ велика.
Альтернативная форма выражения ошибки аппроксимации (см. табл. 7.1) показывает, что если з = =(1/12)ые, то схема Дюфорта — Франкела дает ошибку аппроксимации 0(Лх'). Соответствующая этому точность решения указывается в табл. 7.3. Все же с практической точки зрения использование схемы Дюфорта — Франкела влечет за собой серьезное ограничение на величину М, даже несмотря на то что это ограничение связано с требованиями точности, а не устойчивости, как это было со схемой ВВЦП. 7.1.3. Трехслойная схема Явная форма трехслойной дискретизации уравнения (7.1) в общем случае записывается в виде аТ1+'+ ЬТ~ + сТ~ ' — (сИ.„„Т~ + ет.а,Т~ ) = О, (7.12) где 1.„Т1 = (Т~ ~ — 2Т1 + Т~+ 1)~бхе.
Параметры а, Ь, с, а и е могут быть определены путем разложения каждого из членов уравнения (7.12) в ряд Тейлора относительно узла (1,п) и требования о том, чтобы уравнение (7.12) было согласовано с уравнением (7.1). Примеры подобных процедур даются в п. 3.2.2 и 3.2.3. Эта процедура позволяет переписать уравнение (7.12) в форме, содержащей лишь две варьируемые константы у и р вместо пяти.
В результате уравнение (7.12) заменяется на В+у)(т",+' — т",) у(т,"— т",-') де Ы вЂ” а [(1 — 5) 1.аиТ," + $1.иаТ~ 1 = О. (7 13) Разложение членов уравнения (7.13) относительно узла (1, и) указывает на согласованность данного уравнения с уравнением (7.1), причем ошибка аппроксимации выражается формулой Е~ =азбх ех, ~0 5+ у+ р — 1з )+О (Лх"), (7 14) где з = аб(/Лха. В выражении (7.14) все производные по времени заменяются производными по пространству посредством использования исходного уравнения, как это сделано в табл. 7.!.
Ясно, что уравнение (7.13) имеет ошибку аппроксимации четвертого порядка, если постоянная (1 выражется по формуле (1 = — 0.5 — у+ 1/(12г). (7.15) $7.1. Явные методы Из уравнения (7.13) получается алгоритм +(, + ) 1; ((! — ЯТ1+[)Т1" ) (7.16) где 7.„'„Т =Т>, — 2Т +Т+и Оказйвается, что алгоритм (7.16) обладает лишь условной устойчивостью, причем максимальное значение я, обеспечивающее устойчивость, является функцией у. Это можно установить, применяя к уравнению (7.13) анализ устойчивости по Нейману (см. $4.3). При этом возникает необходимость решения следующего квадратного уравнения относительно 0: 0е(1+ у) — 0 [1+ 27+ 2я (1 — 6) (соя 8 — 1)[+ + [у — 2[)я (соя 9 — 1)) = О. (7.17) Для обеспечения устойчивости необходимо, чтобы было [О[ ( ( 1 при любых значениях 6. В результате можно построить график, показанный на рис. 7.3 и демонстрирующий границу 05 03 О 2.0 4.0 6.0 у Рис.
7.3. Граница области устойчивости для схемы (7.13) и формулы (7.15). области устойчивости в координатах у и я. При у = 0 уравнение (7.17) приводит к обусловленному требованием устойчивости ограничению на я( я ( 0.34), более суровому, чем это делает схема ВВЦП. Точность трехслойной схемы (7.13) исследуется в п. 7.1.4. 7.1.4. ИгЕХ численные резулетаты применения явных схем В данном пункте сравниваются результаты применения схемы ВВЦП (п. 7.1.1), схемы Дюфорта — Франкела (п. 7.1.2), а также трехслойной схемы (п. 7.1.3). Все три метода отражены 19 К. Флетчер, т.
1 60 61 62 63 64 65 66 67 С 6$ С 69 С 70 71 72 73 С 74 С 75 С 76 77 7$ 79 ВО $1 $2 $3 $4 С $5 С Вб С В7 ВВ $9 90 С 91 С 92 С 93 94 95 96 97 9$ 99 С 100 101 102 103 104 105 С 106 С 107 С 10$ С 109 110 111 С 112 С 113 С 114 тоь<1) = 1. тоь<днлх> Тм(1) = 1. тн(днлх> ТР(1) 100.~ТН(1) ТР(дНЛХ) 100.~ТН<днЛХ) Н 0 ЕЛСН Т1НЕ ЗТЕР ЗТЛЯТЗ ЛТ ЗТЛТЕНЕНТ 10 10Н Н+1 1Г(нк .КО.
1)ООТО 15 1Г(НЕ .ЕО. 2)СОТО 13 не 3, З-ьет,, 4тн ОЛОВЕ ВО 12 д 2,3НЛР РРН(д) ЛЗеТНИ) + ВЗ*ТОЬ(д) ВО11Е 13 Кд д- а+К РОН<Л - РОН(д) + СЗ*КЬ(К)*(Об*ТЕ(хд) + ЕЕ*ТОЛ(кд)В 11 СОМТ1МОЕ 12 СОНТ1НОЕ СОТО 17 Нк 2, РОГОЯТ-РЯЛНЯЕЬ 13 ВО 14 д 2,4НЛР 14 ОРн<л - лз <тн(д-1) + тн<д+1)) е взятоь(д) ООТО 17 НЕ 1, ГТСЗ 15 ОО 16 д 2,дНЛР РОВ<4) (1. 2.*$)*ГНИ) + 3 (Тн(д"1) + ТНИ+1)) 16 сомт1Я>е 17 ОО 1$ д 2,3НЛР 1Г(НЕ .ОТ. 1)ТОЬ(д) ТН(д) 1В ТН(д) РОН(д) РО 19 д 2,ЛЧЯР 19 ТР<Л = 100.*ТН(д) Т Т + ОЕЬТ ЗГИРВ .ко.
1>яяхтк<6,20)т, <тн(л,д 1.днлх) 20 ГОЯНЛТ(' Т ',Г5.2,' ТО ',11Р6.2) 1Г НЛХ1НРН Т1НЕ ОЯ НЯХ1НОН НОНВЕЯ ОГ Т1НЕ-ЗТЕРВ ЕХСЕЕРЕО ЕХ1Т ГЯОН ЬООР 1Г<н .ск. Ннлх>сото 21 1Р(Т .ЬТ. ТНЛХ)Сото 10 ОВТЛ1Н БХЛСТ ЗОЬРТ10М ЯНР СОНРЛЯЕ 21 Зпн " О. Рис. 7.4 (продолжение). Га. 7. Одномерное уравнение диффузии ВМЗ 1$ тне йиз ейной ауз Знн/(1 ° + йдн) йиз ВЗЯВТ[ЗТЗ) НВ1ТЕ(6,26)ВИЗ РОВМВТ(' ВМЗ 91$ ',$11А.//) ЗТОР ЕМИ Рис. 7 4 (окончание).
в тексте программы 131РЕХ, распечатка которой приведена на рис. 7.4 и которая является обобщением программы Р1РР (рис. 3.13). Решения ищутся в вычислительной области 0 (к( 1.0 и 2.00 (1 < 9.00 при начальных условиях, заданных при = 2.00, с помощью точного решения (3.42), деленного на 100. Граничные условия имеют вид Т (О, [) = Т (1, /) = 1.0. (7.18) Представление о точности различных схем можно получить путем оценки среднеквадратичной разности между численным н точным решениями при Т = 9.00.
Точное решение рассчитывается с помощью подпрограммы ЕХТКА (рис. 7.6). Схема Дюфорта — Франкела и трехслойная схема нуждаются в задании начальных данных на двух слоях при 1= 2.00 и /= 2.00 — И. Основные параметры, используемые в программе Р1РЕХ (см. рис. 7.4), приводятся в табл. 7.2. Характерная форма представления решения, построенного с помощью трехслойной схемы на сравнительно грубой сетке, показана на рис. 7.6.
Точности различных вариантов решения, полученных при использовании программы Р1РЕХ, сравниваются в табл. 7.3. Чтобы получить возможность сравнения трехслойной схемы четвертого порядка в виде (7.13) со схемами ВВЦП и Дюфорта — Франкела, использовалось два значения 3: 6 =0.3 и 0.41. Кроме того, в таблицу включены два частных случая — ВВЦП с 3 = 116 и схема Дюфорта — Франкела с $ =(1/12)'/'. Оба 115 С 116 111 С 11$ 113 12О 121 22 122 123 126 23 125 С 126 С 123 С 12$ 123 130 131 26 132 13И СВЬЬ ЕХТВВ [дийХ, И1ХЕХ. ВЕ)Л, Р1, ВЬРИ,Т,ТЕ) ВО 22 д - ),дийХ $НР ТЕ(д) — ТН(д) ЗОН ~ ЗНН + ОНР*ОНР СОИТ[НОЕ 1Р(1Рй .ИЕ. 1)НВ1ТЕ(6,20)Т, (То(ш,д 1,дийх) ИВ1ТЕ(6,23)Т,(ТЕ(д),д 1,дНВХ) РОВИХТ(/,' Т ',Р5.2,' ТЕ ',11$6.2,//) 293 З 7.!.
Явные методы этих случая соответствуют резкому уменьшению ошибки аппроксимации. 1 2 зввкоитхне кхтла(знах,илхкх,вкьх, Рт,льва,т,тк) 3С 4 С ЕХЛСТ ВОЬОТ10Н ОР ТИЕ ТААИ$1ЕИТ ВЕАТ СОНВОСТ1ОИ РАОВЬЕН 5 С С тс В ВО 2 д 1,ЗИЛЕ 9 Лз 2-1 10 Х(Ш ВЕЬХ*ЛЗ 11 тк ю) 100. 12 Во 1 И 1,ИАХКХ 13 ЛИ Н 14 ВЛН (2.*ЛН " 1.) 15 ВХН ВвиьР1ьх Ю) 16 ВТН -АЬРИ*ВЛН*ВЛН*Р1*Р1ьт 17 С 1$ С ЬХН1Т ТИК АКОИНЕИТ Е12Е ОР ЕХР(ВТН) 19 С 20 21 22 23 24 25 26 В1ИКИ$10И Х( ° 1),ТЕ(41) 1Р(втн .ьт, - 25.штн - -25.о ВТИ КХР (ВТН) 1Р(ВТН .ЬТ. 1.0К"10)ООТО 2 1 Тк(0) Тк(д) " 400./ВЛН/91*51И!ВХИ)*ВТН 2 СОИТ1ИНЕ АНТОНИ ЕИВ Рис.
7.5. Распечатка программы ЕХТКА. 3 Ькткь, 4тв-Олвка $снене ЛНАХ 11 )ВХЕХ 10 ННАХ 500 ТНАХ 9.00 Т$Т 1.00 ОИГ ОАО Л .)00 ЛЬРН .100Е-01 ВЕЬТ .ЗООЕ+00 ВЕЬХ .100Е+00 Т 9.20 Тв 100.00 $4. 11 69.$0 5$.45 51.16 4$.65 51.16 5$.45 69.$0 $4, 12100.00 Т 9.20 ТЕ 100.00 $4.12 69.ВО 5В.45 51.17 4$.66 51.17 5$.45 69.$0 $4.12100.00 Виа 01Г .4!690-01 Рис. 7.6. Типичная выдача результатов по программе Р)РЕХ (начало).
При к=О.З схема ВВЦП обнаруживает скорость сходи- мости второго порядка. Приближенное значение скорости сходимости было получено по данным об отношении величин среднеквадратичной ошибки на двух сетках: Лх = О.! и 'Ах = 0.05, т. е. г = ) И (КМ8аз-о.) /)д (КМЗал-слв). Схема Дюфорта — Франкела значительно точнее схемы ВВЦП. Таблица 7.2. Параметры, используемые в программе 01РЕХ Параметр Опнсаане МЕ Это связано с близостью данного значения и к специальному з = (1/12) и' = 0.289, для которого достигается скорость сходи- мости четвертого порядка. Трехслойная схема четвертого порядка (при Т = О) оказывается менее точной, чем схема Дюфорта — Франкела на грубой сетке (Ах = 0.20). Обычная тенденция для схем повышенного порядка состоит в том, что для выявления своих преимуществ в точности они нуждаются в более мелких сетках.
По мере возрастания Т точность трехслойной схемы четвертого порядка на заданной сетке уменьшается, хотя четвертый порядок скорости сходимости приближенно сохраняется. При з = 0.41 и Т = 1.0 трехслойная схема четвертого порядка обнаруживает точность, сравнимую с точностью схемы Дюфорта — Франкела на грубой сетке, но является существенно более точной на мелкой сетке. Результаты, показанные в табл.
7.3, свидетельствуют о том, что построение схем повышенного порядка путем уменьшения 1РК аМАХ МАХЕХ )ЧМАХ А1.РН Ь ТМАХ ТЬТ СаАМ ТВ т)ч Т01. ТЕ ПЕ(.Х РЕЕТ Е). Т Х йМБ Гл. 7. Одномерное уравнение диффузии =1: схема ВВ((П =2: схема Дюфорта — франкела =3: трехслойная схема четвертого порядка Управление формой выдачи; 1РК=1 для решения ТП Число точек в направлении х Число членов в выражении точного решення Максимальное чвсло шагов по времени Коэффициент тепловой диффузии сс = стах/ЬГа Максимальное время Начальное время Параметр у в соотношении (7.!6) Размерный температурный массив Безразмерный температурный массив Т" в (7.16) Безразмерный температурный массив Т"-' в (7.16) Размерный температурный массив точного решения Лх 61 l Коэффициенты в выражении разностного оператора Ь „ Время Пространственная координата, О ( х ~ 1.О Среднеквадратичная ошибка решения 6 7.2.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
