Ортега Дж., Пул У. Введение в численные методы решения дифференциальных уравнений. Под ред. А.А.Абрамова (1986) (1095855), страница 16
Текст из файла (страница 16)
Рассмотрим уравнение у'= — 1ооу+1оо, у(о)=у . что метод Эйлера должен дать достаточную точность при сравнительно большом шаге л. Однако нз (2.5.30) видно, что если Й > 0,02, то ! 1 — 1006 ! > 1 и значения у„с увеличением номера шага начинают быстро расти, свидетельствуя о неустойчивости. Если сравнить точные решения (2.5.29) и (2.5.30), то видно, что частные решения уравнений (2.5.25) и (2.5.27) тождественно совпадают (и равны 1). Величина (1 — 1006)" служит аппроксимацией экспоненциального члена е 'ОО и действительно является хорошим приближением при малых 6, но это приближение быстро становится неудовлетворительным, когда л достигает значения 0,02, И хотя этот экспоненциальный член после х = 0,1 фактически не вноснт никакого вклада в решение, для сохранения устойчивости метот Эйлера по-прежнему требует, чтобы этот член аппроксимировался достаточно точно.
Эта ситуация характерна для жестких уравнений: решение содержит член, вклад которого очень мал, однако обычные методы для сохранения устойчивости требуют, чтобы этот член аппроксимировался достаточно точно. Эта проблема очень часто возникает при решении систем уравнений.
Рассмотрим, например, уравнение второго порядка у" + 101у'+ 100у = 0. (2.5.31) Мы можем, как это показано в приложении 2, преобразовать (2.5.31) в эквивалентную систему двух уравнений первого порядка, но для наших целей достаточно рассмотреть его в исходной форме. Общее решение (2.5.31) имеет вид У(х) = с1 Е 1ООх + С1Е Так что решением, удовлетворяющим начальным условиям у(0) = 1,01, у'(О) = — 2, является функция 1 у(Х) Е-1ООх + Š— х 100 (2.5.32) После того как х достигнет значения порядка 0,1, вклад первого члена в решение будет очень мал. Если, тем не менее, применим к соответствующей уравнению (2.5.31) системе первого порядка метод Эйлера, то мь1 столкнемся с той же самой проблемой, что и в предыдущем примере: придется выбрать шаг достаточно малым, чтобы точно аппроксимировать член е ' О", несмотря на то что его вклад в решение очень мал. Приведенный пример демонстрирует сущность проблемы жесткости, возникающей прн решении систем уравнений.
Обычно в таких задачах независимой переменной является время, а в моделируемой физической проблеме возникают быстро затухающие переходные процессы, с которыми численной схеме приходится иметь дело н после того, как они уже практически ничего не вносят в решение. Общий подход в решению проблемы жесткости заключается в использовании неявных методов. Подробное рассмотрение этих вопросов выходит за рамки настоящей книги, и мы здесь ограничимся демонстрацией значения неявнь1х методов для решения жестких задач, применив один нз простейших таких методов к задаче (2.5.25) . Метод решения общего уравнения у = ~(х, у), состоящий в использовании формулы у„„=у„+Анхо ),у„„), (2.5.33) назьвается обратным методом Эйлера.
Этот метод по форме совпадает с обычным методом Эйлера, за исключением того, что значение т вычисляется в точке (х„+,, у„„), а не в (х„, у„), так что метод является неявным. Применяя (2.5.33) к задаче (2.5.25), получаем разностное уравнение у„+, =у„+ Ц вЂ” 100У„, + 100), (2.5.34) которое можно переписать как у„„=(1+ 1006) '(у„+ 1001). (2.5.35) Точное решение уравнения (2.5.35), как легко проверить (см. упражнение 2.5.10), задается выражением Ул = (Уо — 1) (1 + 100п) + 1 (2.5.3б) которое при начальном условии у, = 2 принимает вид (1+ 100() — и+1 (2.5.37) Отсюда видно, что решение ведет себя устойчиво при любой величине шага. Обратный метод Эйлера, как и сам метод Эйлера, имеет только первый порядок точности, так что к лучшим результатам привело бы использование метода Адамса — Моултона второго порядка Ь у„+, у„[~„+ Дх„+, „у„„, )1, (25.38) который называют также правилом трапеции.
В упражнении 2.5.14 предлагается применить этот метод к задаче (2.5.25) . Использование неявного метода для решения уравнения (2.5.25) выглядит обманчиво простым потому, что это уравнение линейно и соответствующее разностное уравнение (25.34) легко разрешается относительно ул~,. Если бы дифференциальное уравнение было нелинейным, то для определения у„+, потребовалось бы на каждом шаге решать нелинейное уравнение.
В общем случае системы дифференциальных уравнений на каждом шаге. приходится решать систему уравнений (линейную или нелинейную, в зависимости от вида дифференциальных уравнений). Конечно, это связано с большими затратами машинного времени, но для эффективного решения жестких уравнений численный метод должен включать те или иные неявные элементы.
Дополнительные замечания и ссылки 2.5 По теории устойчивости решений дифференциальных уравнений имеется обширная литература. В качестве хорошо написанного вводного курса можно рекомендовать книгу [39]. Мы привели основные результаты теории линейных разностиых уравнений только для случая, когда корни характеристического уравнения различны. При наличии кратных корней в решение входят полиномиальные члены, совершенно аналогичные соответствующим членам для дифференциального уравнения.
Более детальное и зложение теории линейных разностных уравнений приведено, например, в книге (921. Метод 12.5.6) возникает естественным образом в результате дифференцирования интерполяционного полинома; вывод имеется в книге (92,с. 219[. 5.
Дж. Ортега Основные результаты по теории устойчивости многошаговых методов получены Далквистом в 1950-х годах; подробное изложение этой теории имеется в работе [92 [ С тех пор был предложен ряд новых определений устойчивости. В частности, были введены понятия жесткой устойчивости и А-устойчивости, характеризующие те типы устойчивости, которыми должны обладать численные методы для успешного решения жестких уравнений. Дальнейшие сведения о различных определениях устойчивости можно найти, например, в книгах [11, 37 [.
Хотя к настоящеагу времени предложен и действительно используется целый ряд методов решения систем жестких уравнений, вопрос о наиболее эффективном способе решения таких систем все еще остается открытым. Обсуждение многих методов и дальнейшие сведения о жестких уравнениях имеются в книгах [11, 37, 38] и в обзорной статье [102[. УПРАЖНЕНИ>/ 2.5 2.5,1. Положив г = у ', покажите, что задача (25.1) — (2.5.2) эквивалентна системе первого порядка у' = г, г' = 10г + 11у с начальными условиями у(О) = 1, г(0) = — 1.
Попробуйте решить эту систему численно, воспользовавшись любым методом из этой главы, и проанализируйте полученные результаты. 25.2. Попробуйте решить задачу (2.5 5) численно, воспользовавшись каким-либо методом из разделов 2.2 и 2.4, и проанализируйте полученные результаты. 2.5.3. Покажите, что метод (2.5.6) имеетвторой порядок точности. 2.5.4. Проведите численные расчеты по алгоритму (2.5.9) — (2.5.10) прн разных значениях й.
Проанализируйте получающиеся результаты. 2.5.5. Выпишите решение разносгного уравнения у„+1 = (5/2)у„+ У„1, У, =У, = 1 через корни его характеристического уравнения. Исследуйте поведение последовательности(Ул) пРн л 2.5,6. Выберите значение у, таким образом, чтобы соответствующее решение разностного уравнения (2.5.9) с у, = 1 стремилось к нулю при л — '. Составьте программу, реализующую алгоритм (2.5.9) как с найденным значением у,, так и со значением у,, задаваемым формулой (25.10) .
Проанализируйте получаюшнеся результаты. 2.5.7. РассмотРите метод У„+1 = Ул 3 + (4»/3) (2/и — /'„1 + 2/я а) известный как метод Милна. Исследуйте этот метод иа устой ивость и строгую устойчивость, 2.5.8. Составьте программу, реализующую метод Эйлера (2.5.27), и проведите расчеты при различных значениях Ь, как меньших, так и больших, чем 0,02. Проанализируйте полученные результаты. 25.9.
Система у' = г, г' = — 1ООУ вЂ” 101г является системой первого порядка, эквивалентной уравнению второго порядка (2.5.31). Примените метод Эйлера к этой системе с начальными условиями у(0) = 2, г(0) = — 2 и определите экспериментально, насколько малым должен быть шаг», чтобы счет был устойчивым. Попытайтесь подтвердить ваши выводы о необходимой величине шага й аналити юскн. 2.5 ЛО. Найдите решение разностного уравнения у„+1 = сУл + Ы как функцию начального значения у,.
Используйте полученный результат для доказательства того, что формула (2.5.36) дает решение уравнения (2.5.35). 2.5.11. функция у (х) = е " является решением задачи у" = у, у(0) = 1, у '(0) = — 1. Устойчиво лн это решение? Докажите ваше утверждение. 25Л2. Установите, какие нз следующих методов устойчивы и какие строго устойчивы: а) у»+ У» + (6/» — 3/» 1+ 3/»-з)' И б) у»+ =у, + — (/»+ + 2/' +/» ); » в) У»+1 =ЗУ» -2У» — 1+ (Г»+1+2/»+/»-1) 2 2.5.13.
Рассмотрите многошаговый метод у»+1 =. (1/2) (у» + у» 1) + (3л/4) (3/»вЂ” — 1) . Является ли этот метод устойчивым? Строго устойчивым? Каков его порядок? 2.5.!4. Примените правило трапеции (2.5.38) к уравнению (2.5.25). 66 Глава 3 ЗАКРЕПЛЕНИЕ НА ОБОИХ КОНЦАХ: ДВУХТОЧЕЧНЫЕ КРАЕВЫЕ ЗАДАЧИ 3.1. Задача диффузии Многие задачи научного программирования представляют собой нелинейные дифференциальные уравнения в частных производных вида (рих)х +ь(» и) (3.1.1) с соответствующими начальными и граничными условиями, такими как и(О,х) =а(х), 0-'х <1 (начальное условие), (3.1.2) и(г, О) = а, и(г, 1) = Д, г ) 0 (граничные условия) .
с граничными условиями и(0) = а, о(1) = 11. (3.1.5) Так как и — функция только переменной х, уравнение (3.1.4) является обыкновенным дифференциальным уравнением, которое мы перепишем в виде (ро ) + я(х, и) = О. 67 Здесь р и д — заданные функции х, 8 — заданная функция двух переменных, а и р — заданные константы, и — решение, которое требуется найти; через и, обозначена частная производная ди/дг. Мы вернемся к изучению зависящих от времени дифференциальных уравнений в частных производных вида (3.1.1) в гл. 7. В этой же и двух последующих главах будем заниматься соответствующими стационарными задачами. Во многих ситуациях, описываемых уравнениями типа (3.1.1), в которые в качестве независимых переменных входят пространство и время, главный интерес представляет поведение решений при ~ — '.
Если решение задачи имеет предельное значение при г —, т.е. и(г,х)- и(х) при г-, хЕ10,1], (3.1.3) то говорят, что функция о является стационарным решением. Один из подходов к нахождению такого стационарного решения состоит в вычислении зависящего от времени решения и(г, х), соответствующие методы рассматриваются в гл. 7. Такой подход, конечно, является неизбежным, если нас интересует и поведение решения как функции от г до достижения стационарного состояния. Другой подход заключается в следующем: мы замечаем, что стационарное решение должно удовлетворять уравнению (3.1.1), в котором частная производная по времени положена равной нулю, т.е. уравнению (ри„)„+ я(х, о) = 0 (3.1.4) Поскольку граничные условия (3.1.5) задаются на обоих концах интервала, уравнение (3.1.6) вместе с граничными условиями (3.1.5) называют двухточечиой краевой задачей.
Численное решение таких задач, формулировка которых более подробно рассматривается в приложении 2, и будет предметом этой и двух последующих глав. Отметим, что отрезок [0,11 выбран исключительно для удобства, а не по каким-либо иным специальным соображениям. Если уравнение (3.1.6) определено на произвольном конечном интервале, то задачу всегда можно свести к отрезку [0,11 (см. упражнение 3.1.1) . Рис. 3.1. Сферическая клетка В качестве примера уравнения типа (3.1.6) рассмотрим следующую задачу из биофизики.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
