Fletcher-1-rus (1185917), страница 55
Текст из файла (страница 55)
В качестве альтернативного варианта можно воспользоваться схемой первого порядка по времени типа ВВЦП прн условии, что второй слой начальных данных рассчитывается на мелкой сетке, тем самым компенсируя первый порядок точностн. Для этой цели может быть полезной экстраполяция по Ричардсону (и. 4.4.1). 315 $7.4, Метод прямых различными методами решения систем обыкновенныхдифференциальных уравнений 1тхеаг, 1971; 1.ашЬег1, !973; Бежаге( е1. а1., 1984).
Следует, однако, подчеркнуть, что построение полудискретной формы вносит ошибку, связанную с пространственной дискретизацией. В результате этого «наилучшим» вариантом при выборе метода решения получаемой системы оказывается обычно алгоритм более низкого порядка, чем если бы это была система обыкновенных дифференциальных уравнений, не связанная с аппроксимациями.
Прн решении задач с начальными данными, описываемых обыкновенными дифференциальными уравнениями, наиболее эффективными обычно оказываются либо линейные многошаговые методы, либо методы Рунге— Кутты. Последовательное применение соотношения (7А1) ко всем внутренним узлам приводит к системе уравнений, записываемых в виде (7.42) дп/Ж = Ам. Если бы пространственный член уравнения (7.1) был нелинейным, то уравнение (7.42) можно было бы записать в более общей форме е(а/й = Г. (7.43) Скалярные эквиваленты векторных уравнений (7.42) и (7.43) имеют вид Ни/Ж = Аи, (7.44) ни/ж = /. (7.45) Применяя к уравнению (7.45) линейный многошаговый метод общего вида, можно написать (7.46) где ищется решение и"+ .
Если 8 = О, это решение получается в явной форме; если же 6 Ф О, то (7.46) представляет собой неявный алгоритм для определения ия+ . Простейший случай, когда и = 1, приводит к соотношению а,и + аои = б1 (йо/ + й,/ ~ ), (7.47) которое охватывает, в частности, схему Эйлера (а~ — — — ио = 1, ро = 1, ~~ = О), й+' — й = Л1/". (7.48) Последняя формула, если ее применить к линейной системе (7.42) в форме (7.41), совпадает со схемой ВВЦП, соответст- з!е Гл. 7. Одномерное уравнение диффузии вующей уравнению (7.5). Соотношение (7.47) охватывает так- же и трапецоидальную схему ил+! — ил = Л1 (0.57" + 0.51""~).
(7.49) В применении к (7.41) трапецоидальная схема совпадает со схемой Кранка — Николсона в форме (7.22). Двухшаговый вариант схемы (7.46) имеет вид пой+ +а!й+ +пой=И(йод" +И"+ +~Я" ). (7.50) При применении к (7.41) соотношения (7.50) при рз — — 0 схема (7.50) позволяет получить явное решение и~+~ и совпадает со схемой (7.12). Если 5з Ф О, то схема (7.50) в применении к (7.41) охватывает как схему (7.25), так и схему (7.26).
Например, схема ТСЧН (п. 7.2.3) соответствует выбору аз — — 1.5, а! = — 2.0, ао = 0.5, ()о = р! = О, рз = 1.0. Фактически все схемы, разработанные в $ 7.1 и 7.2, могут рассматриваться как принадлежащие классу линейных миогошаговых методов. При решении обыкновенных дифференциальных уравнений используют линейные многошаговые методы более высокого порядка (т. е. включающие большее число шагов), и притом в варианте предиктор — корректор. В этом случае явная схема обеспечивает предиктор, а неявная схема — соответствующий корректор [беаг, 1971).
Однако если стратегия высокого порядка применяется к системе уравнений, полученных в результате полудискретизации дифференциального уравнения в частных производных, например к системе (7.41), то для хранения информации о векторах нл+' и Гл+! необходим большой объем памяти. Кроме того, использование схемы высокого порядка для маршевого продвижения по времени было бы неэффективным, так как если точное решение не является очень гладким, то преобладающее влияние на ошибки решения будет, вероятно, оказывать дискретизация низкого порядка пространственных членов.
Необходимость избежать хранения в памяти большого числа решений, соответствующих предшествующим шагам по времени, побуждает воспользоваться одношаговыми методами решения полудискретных систем, такими, как методы Рунге — Кутты. В применении к уравнению (7.45) схема Рунге — Кутты общего вида, содержащая Я этапов, может быть записана в форме и ил+! ил+ Л! Х~~ с ~е (7.51а) е 1 Згг $7л, Метод прямых где ! =!Ф ! ы~ у ! ЖХ ьг), 1 2. я. !751м г ! а, = ~ Ь„, г = 1, 2, ..., )г. (7.51с) (7.52) Существуют и многочисленные трехэтапные явные схемы Рунге — Кутты третьего порядка 1(.ашЬег1, 1973]. Нижеследующая четырехэтапиая и четырехшаговая схема Рунге— Кутты четвертого порядка получила широкое распространение при решении задач с начальными значениями, определяемых обыкновенными дифференциальными уравнениями ил+1 ия + ((п+ 27 + 27 ~+ ~~~~) (7.53) где и = и" + 0.5 огг'", 7"' = 7 (г'"+и', и'), й'= и" + О.бо!)', г" = г(г"+и", и"), й" = и" + Ьгг", 1*" = 7 (7", и"').
Очевидно, что схемы Рунге — Кутты повышенного порядка требуют большего числа оценок для 1. При решении таких полудискретиых систем, как (?.43), серия повторных оценок всех членов, связанных с пространственной дискретизацией, зачастую играет определяющую роль в подсчете суммарного времени исполнения. Формулы (7.51) охватывают неявные схемы Рунге — Кутты, иеэкономичные в вычислительном плане, так как уравнение (7.51Ь) является нелинейным и должно решаться с применением итераций для каждого из значений гт и на каждом шаге по времени. Поэтому ббльший интерес проявляется к явным схемам Рунге — Кутты, в которых верхние пределы суммирования в формулах (7.51Ь) и (7.51с) заменяются на г — 1. Так, например, схему Эйлера (7.48) можно интерпретировать в качестве одноэтапной (14 = 1) явной схемы Рунге — Кутты при с!=1,Ь! =О,а =О.
Существует бесконечное число двухэтапных явных схем Рунге — Кутты второго порядка. Типичной схемой такого рода является улучшенная схема Эйлера и'=и +огг", и"+' =и" + 0.бог(1" + 1'). зсв Гл. 7. Одномерное уравнение диффузии Если скалярное уравнение (7.44) рассматривается вместе с начальным условием и(О)= ие, то соответствующее точное решение имеет форму и (С) = иве ы = и Чтобы физическое решение соответствовало устойчивому процессу в интервале 0 < 1< оо, необходимо, чтобы Л < О.
Абсолютная устойчивость линейного многошагового метода (7.46) анализируется путем применения указанной схемы к уравнению (7.44) с последующим формированием полинома устой- чивости Х (ас — ЛЫРс) г'=О. с-о (7.55) При заданном ЛА1 абсолютная устойчивость достигается 11атЬег1, 1973] тогда, когда все корни г, полинома (7.55) удовлетворяют условию 1г,1< 1. (7.56) Это условие в принципе эквивалентно ограничению, выполнение которого требовалось для устойчивости решения дискретизированного дифференциального уравнения в частных производных при применении анализа устойчивости по Нейману 6 4.3). В случае схемы Эйлера абсолютная устойчивость достигается, если — 2 < И1 < О. Поскольку Л < О, это соответствует условию Ы < 2/12~. Для трапецоидальной схемы абсолютная устойчивость связана с условием ЛАС < 0; следовательно, трапецоидальная схема устойчива при любом выборе Ы.
Критерий абсолютной устойчивости, эквивалентный (7.56), описан в книге 11аспЬег1, 1973) и приводит к ограничению — 2 <И( <О, относящемуся как к одноэтапной схеме Эйлера, так и к явной двухэтапной схеме Рунге — Кутты (7.52). Однако четырехэтавная явная схема Рунге — Кутты (7.53) имеет слегка расширенный интервал устойчивости — 2.78 ) <ЛМ <О.
При применении метода прямых мы встречаемся с системами обыкновенных дифференциальных уравнений. Следовательно, уместно проверить абсолютную устойчивость уравнений в форме (7.42). Как и прежде, можно получить полином устойчивости (7.55), если вместо Л вводить Ле — собственные значения матрицы А.
В общем случае величины Ла комплексные. Однако чтобы задача имела физически устойчивое решение, необходимо, чтобы вещественная часть Лд, т. е. КеЛю была положительной. 3!9 $7.4. Метод прямых Например, применение схемы Эйлера (7.48) к уравнению (7.41) дает следующий критерий абсолютной устойчивости: — 2((йеАя ЛГ(~0. (7.57) Для системы из М уравнений, соответствующей представлению (7.41), матрица А получается трехдиагональной и ее собственные значения определяются формулами (9.48), т. е.
а ( — 2 + 2 соя (ял)/(М + 1)) Лхх В результате наиболее строгое ограничение на величину А(, вытекаюшее из условия (7.56), имеет форму 4а ат О.за — 2» (— — или Ы ~( — '-3-. 4хя ах (7.59) 1) = ,)' аде е ех + 1)„, «-1 (7.60) где собственные значения ),» определяются формулами (7.58), а величины ея — соответствующие собственные векторы. Коэффициенты ая определяются начальными условиями ()е(х~). Как ясно из формул (7.58), диапазон изменения собственных значений Таким образом, здесь получено то же самое ограничение иа величину И, обусловленное устойчивостью, как и для схемы ВВЦП (7.5).
Этого и следовало ожидать, так как метод Эйлера с шагами по времени в применении к уравнениям (7.41) совпадает со схемой ВВЦП. Применение к (7.41) трапецоидальной схемы (7.49) приводит к ограничению Ке),е.ог (О. Но так как в силу (7.58) )сея,е (О, то ограничение на А( отсутствует. Трапецоидальная схема, примененная к (7.41), совпадает со схемой Кранка — Николсона (7.22). Как можно предположить на основании этих примеров, абсолютная устойчивость алгоритма решения системы обыкновенных дифференциальных уравнений эквивалентна устойчивости алгоритма решения родственной задачи о решении дискретизированного дифференциального уравнения в частных производных ($ 4.3).
Если система, полученная при помощи метода прямых, является нелинейной, например система (7.43), то определение устойчивости алгоритма зависит от собственных значений якобиана дГ/д(). Однако сами собственные значения зависят от решения и могут нуждаться в корректировке в процессе построения решения. Точное решение системы уравнений (7.42), образуемой из (7.41), может быть записано в виде Гл.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
