Численные методы. Соснин (2005) (1160462), страница 17
Текст из файла (страница 17)
То есть надо выбирать очень мелкий шаг интегрированияв соответствии с поведением функции f.Пример 5.9. Одношаговый метод Гира первого порядка аппроксимации.yn+1 − yn= f (tn+1 , yn+1 ).τСогласно модельной задаче, f (tn+1 , yn+1 ) = λyn+1 , откудаyn+1 (1 − τ λ) = yn .Для устойчивости метода мы требуем, чтобы |yn+1 | 6 |yn |. Это неравенство выполнено всегда,так как τ > 0 и λ < 0, то есть сеточная функция не возрастает при любых шагах интегрирования.Мы получили, что данный метод Гира является абсолютно устойчивым. Неоспоримым достоинствомэтого метода по сравнению с методами Адамса является то, что шаг интегрирования мы выбираемс оглядкой лишь на требуемую точность, а к недостаткам можно отнести то, что на каждом шагеприходится решать нелинейную систему уравнений.При рассмотрении модельной задачи у нас возникают так называемые сеточные уравнения. В общемвиде их можно записать так:yn+1 = an yn + an−1 yn−1 + . .
. + a0 y0 .В правой части стоит линейная комбинация значений сеточной функции. По сути дела, это некоторое рекуррентное соотношение. В общем виде оно решается с помощью характеристического многочлена, но мы ограничимся случаем, когда решение можно задать в виде yi = q i .На основе такого представления решений аппарат исследования на устойчивость таков.
Если удается показать, что |q| 6 1 для любых параметров метода (от них будет зависеть рекуррентное соотношение), то общее решение будет монотонно убывать. Если же существует хотя бы один набор параметровтакой, что |q| > 1, то условие устойчивости будет нарушено для данного метода.5.5. УСТОЙЧИВОСТЬ ЧИСЛЕННЫХ МЕТОДОВ РЕШЕНИЯ ЗАДАЧИ КОШИ87Рассмотрим пример. Проверим на устойчивость явные одношаговые методы Адамса.yn+1 − yn= f (tn , yn ).τРешение сеточного уравнения будем искать в виде yn = q n . Так как схема сводится к видуyn+1 = (1 + τ λ)yn , то q будет равным 1 + τ λ.Для устойчивости метода необходимо, чтобы |1 + τ λ| 6 1, но это условие мы уже получали выше.Таким образом, пока наш метод проверки на устойчивость ничего нового не дал, но мы увидим вполученных формулах закономерности.Если применить этот аппарат к методам Гира, опять получаем знакомую формулу:|q| =1.|1 − τ λ|В общем случае для сеточных уравнений, в которых задано более двух узлов, основание решенияq может быть, вообще говоря, комплексным числом.
Таким образом, проведем некоторое обобщениеаппарата исследования устойчивости — будем считать, что в модельной задаче λ ∈ C.Еще раз применим метод Адамса к решению модельной задачи, но уже с комплексным λ. Обозначивв получившихся выкладках τ λ = µ ∈ C, получим:|1 + µ| 6 1.=µ−1µ<µДля всех точек, лежащих внутри области, соответствующие методы Адамса устойчивы.Аналогично для методов Гира:|1 − µ| > 1.Как нетрудно заметить, метод устойчив во внешней области.Определение.
Назовем областью устойчивости численного метода решения задачи Коши дляОДУ ту область значений µ, в каждой точке которой |q| 6 1.Определение. Численный метод называется А-устойчивым3 , если его область устойчивости содержит отрицательную полуплоскость µ.В соответствии с этим определением, метод Гира — А-устойчив, а метод Адамса не является Аустойчивым.Переход с оси на комплексную полуплоскость привел к тому, что можно показать, что явных Аустойчивых методов не существует, а среди неявных нет А-устойчивых методов выше второго порядкааппроксимации.3 отабсолютно устойчивый.88Глава 5. ЧИСЛЕННЫЕ МЕТОДЫ РЕШЕНИЯ КРАЕВЫХ ЗАДАЧ=µµ<µ1Все это вытекает из «плохого» определения устойчивости. Попробуем немного исправить эту ситуацию.Определение.
Линейный многошаговый метод называется А( α )-устойчивым, если область егоустойчивости содержит угол |arg(−µ)| < α.=µµαα<µ−µ = rei(π−ϕ)=⇒Заметим, чтоµ = reiϕ=⇒arg(−µ) = π − ϕ.Можно доказать, что для введенного определения устойчивости выполняется следующее утверждение.Утверждение 5.1. Среди явных линейных m-шаговых методов нет А(α )-устойчивых.Это следует учитывать при решении химических, биологических задач, где часто решение имеетэкспоненциальный вид.Среди неявных методов существуют А(α )-устойчивые методы. Рассмотрим пару примеров такихметодов.Пример 5.10. Двухшаговый метод Гира второго порядка аппроксимации.13yn − 2yn−1 + yn−2 = τ f (tn , yn ).22Пример 5.11. Четырехшаговый метод Гира четвертого порядка аппроксимации.1(25yn − 48yn−1 + 36yn−2 − 16yn−3 + 3yn−4 ) = τ f (tn , yn ).12Найдем область устойчивости для рассмотренного выше двухшагового метода Гира второго порядка аппроксимации.
Как уже делали ранее, зададим f (tn , yn ) = λyn и, обозначим, τ λ = µ, тогда5.5. УСТОЙЧИВОСТЬ ЧИСЛЕННЫХ МЕТОДОВ РЕШЕНИЯ ЗАДАЧИ КОШИ89сеточное уравнение будет выглядеть так:31yn − 2yn−1 + yn−2 = µyn .22(5.16)Как мы уже договаривались, нас интересуют только решения вида yn = q n . Подставив это значениев (5.16), получим:13 2q − 2q + = µq 2 .22или:13 2µ = − + 2.2 q2qНайдем те µ, для которых |q| < 1. При этом комплексная плоскость разбивается на две области —ту, где метод устойчив, и — где неустойчив. При этом, так как на границе |q| = 1, то q = eiϕ иуравнение границы примет вид:31311− 2e−iϕ + e−i2ϕ = − 2 cos ϕ + 2i sin ϕ + cos 2ϕ − i sin 2ϕ =2222231= − 2 cos ϕ + (2 cos2 ϕ − 1) + i(2 sin ϕ − sin ϕ cos ϕ).22µ(ϕ) =Сделаем замену переменной x = cos ϕ (x ∈ [−1; 1]), получим:pµ = 1 − 2x + x2 ± i(2 − x) 1 − x2 .=µ|q| < 1|q| > 14 <µ|q| = 1Вообще то неплохо было бы отметить здесь еще область определения, а то, например, при µ = 23итерационный процесс неопределен.
Проверить область определения мы√ предоставим читателю. Под1±i1/2ставим точку µ = − 32 , она лежит во внешней части плоскости, а q =, |q| < 1 и, таким образом,3мы получили, что внешняя область является областью устойчивости данного метода Гира.Хотя при реализации неявных методов требуется больше времени на просчет одного шага, но ихплюс в том, что мы можем выбирать произвольный шаг — тот, который нам нужен.Численные методы решения задачи Коши для систем ОДУИсследуем методы решения задачи Коши для систем ОДУ.
Будем искать приближенное решениедля такой системы:(ut = f (t, u), t > 0;u|t=0 = u0 ,где u = (u1 , u2 , . . . , un ) и f = (f1 , f2 , . . . , fn ).90Глава 5. ЧИСЛЕННЫЕ МЕТОДЫ РЕШЕНИЯ КРАЕВЫХ ЗАДАЧМетоды Рунге-Кутта и линейные одношаговыеЧтобы стало яснее, приведем пример задачи:ut =vt =u(0)=v(0) =методы легко переносятся на системы уравнений.f (t, u, v);g(t, u, v);u0 ;v0 .Покажем, как к этой задаче применяется метод Рунге-Кутта четвертого порядка аппроксимации(и точности соответственно). Обозначим u(tn ) = yn , а v(tn ) = zn , тогда формулы для подсчета следующих сеточных значений будут таковы:1yn+1 − yn= (K1 + 2K2 + 2K3 + K4 );τ6 zn+1 − zn = 1 (M + 2M + 2M + M ).1234τ6Параметры Ki , Mi , i = 1, 2, 3, 4 вычисляются по следующим формулам:K1K2K3K4= f (tn , yn , zn ),M1 = g(tn , yn , zn );= f (tn + τ2 , yn + τ2 K1 , zn + τ2 M1 ), M2 = g(tn + τ2 , yn + τ2 K1 , zn + τ2 M1 );= f (tn + τ2 , yn + τ2 K2 , zn + τ2 M2 ), M3 = g(tn + τ2 , yn + τ2 K2 , zn + τ2 M2 );= f (tn + τ, yn + τ K3 , zn + τ M3 ), M4 = g(tn + τ, yn + τ K3 , zn + τ M3 ).Мы долго и упорно говорим, что условно устойчивые методы хуже устойчивых, но не привели неодного примера, объясняющего разницу между ними.
Исправим ситуацию. Представим, что мы моделируем два процесса, причем один из них протекает существенно быстрее другого. Система уравненийдля такой задачи будет иметь вид: 0u + a1 u1 = 0, a1 > 0; 10u2 + a2 u2 = 0, a2 > 0;u1 (0) = u0 ;u2 (0) = u0 ,а точные решения соответственно:u1u2= u0 e−a1 t ;= u0 e−a2 t .Как мы уже говорили, a1 a2 — константы суть характерное время протекания первого и второгопроцессов. Нас интересует описание системы при достаточно больших t > t∗ .
Пусть мы используемдля решения этой задачи условно устойчивый метод (например, какой-нибудь из методов Адамса).При этом на шаг τ мы должны сделать ограничение τ λ > −2, где роль параметра λ играют a1 и a2 ,то есть λ = −a1 или λ = −a2 . Неравенства должны выполняться одновременно, то есть τ 6 a21 (всилу того, что a1 a2 ).Теперь представим, что нам важно исследовать первый процесс в то время, когда второй уже невлияет на систему. Тем не менее, от поведения второго процесса зависит шаг аппроксимации (он будетуменьшаться), поэтому весьма вероятно, что мы будем выполнять «лишнюю» работу. Отсюда можносделать вывод, что в случае решения систем ОДУ, описывающих разномасштабные процессы, надоиспользовать абсолютно устойчивые методы (например, методы Гира).Опишем еще одно свойство систем уравнений. В качестве примера системы возьмем такую: u∂u∂tt=0= f (t, u),= u0 .0 < t 6 T;(5.17)5.6. ИНТЕГРО-ИНТЕРПОЛЯЦИОННЫЙ МЕТОД ПОСТРОЕНИЯ РАЗНОСТНЫХ СХЕМ91uu0|t∗u2 (t)u1 (t)tПоставим ей в соответствие матрицу из производных (якобиан):∂fj (t, u)A(t, u) =∂uiи обозначим ее собственные значения за λk (t), t ∈ [0; T ].Определение.
Система (5.17) называется жесткой, если выполнены два условия:Re λk (t) < 0 ∀ k, t ∈ [0; T ];max |Re λk (t)|supt∈[0; T ]kmin |Re λk (t)| 1 — много больше единицы.kТакие системы чаще всего решают неявным абсолютно устойчивым методом.5.6Интегро-интерполяционный метод построения разностныхсхемПерейдем к рассмотрению более сложных краевых задач. Для начала исследуем применение численныхметодов для решения такой модельной задачи: (k(x)u0 (x))0 − q(x)u(x) + f (x) = 0, 0 < x < l;−k(0)u0 (0) + βu(0) = µ1 ;(5.18)u(l) = µ2 .Это краевая задача для обыкновенного дифференциального уравнения второго порядка с заданными функциями k, q, f и неизвестной функцией u.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
