Форсайт Дж., Малькольм М., Моулер К. - Машинные методы математических вычислений (1040536), страница 24
Текст из файла (страница 24)
оо -зооо -гбао -гала -зб оо -1аоо -з ао оао сау Ь Рис. 6.6. Общая ошибка в у(1) при использовании метода Эйлера для решения задачи у'=.у, у(0)=-1. Машина Нопеутчеи 6050, обычная точность 6).=2, 1= =27). для этой задачи, очевидно, равно приблизительно 2-". При решении обыкновенных дифференциальных уравнений обычно нужно значительное машинное время, чтобы ошибки округлений могли ощутимо накопиться.
Например, вычисления для рис. 6.6 потребовали около 16 минут времени процессора машины НопеушеП 6060. Однако, довольно часто приходится интегрировать системы уравнений, для которых необходимы часы машинного времени. ю. зхдхчл коши >за 6А. Методы Хотя для численного решения обыкновенных дифференциальных уравнений существует множество различных пошаговых методов, каждый попадает в одну из четырех общих категорий. Методы рядов Тейлора Если у(г) — гладкое решение, то оно имеет разложение по Тейлору у((+й) =-д(~)+йд'(1)+ — ",, у" (()+... Метод Эйлера можно рассматривать как приближение, использующее первые два члена ряда Тейлора.
Если могут быть вычислены производные функции у более высокого порядка, то можно получить метод порядка р, полагая у у ьг у, .=д (йу' ~ д"+ д" > ( у~о Первый отброшенный член дает оценку локальной ошибки дискретизации и потому может быть использован для выбора величины шага. Производные у'(х), у" (х) и т.
д. можно выразить через частные производные функции ). При атом важно различать функцию )(д, 1) двух независимых переменных и функцию ~(у(1), т) одного независимого переменного„полученную подстановкой решения д в )'. Исходя из основного дифференциального уравнения д' (г) =-) (д ((), т), можно продифференцировать обе части по 1, что дает У" =-)хд'+ )с =- 1„~+)б затем у"' = 1т~'1 + 1р ' ~+ ~р'~>+ 2~о')+)и и т. д.
Этот процесс называется полным дифференцированием. Например, если ~(у, ()=д'+Р, то у' =у'+Р, у" = 2у'+ 2у~Р+ 21, у"' = бу' (-бу'('+ 4у(+ 21'+ 2. Заметьте, что даже для многочлена ) от двух переменных выра- >кения полных производных усложняются по мере роста порядка. Можно составить программу, которая вычисляет полные про- изводные для некоторых классов функций. Для дифференциаль- ных уравнений, определяемых такими функциями, методы рядов, 6.ь метОды 1зт Тейлора могут быть очень эффективны. Однако в общем случае применение этих методов слишком ограничено, чтобы представлять для нас дальнейший интерес.
Методы Рунге — Кутта Метод Рунге — Кутта предназначен для аппроксимации методов, основанных на рядах Тейлора, но без явного вычисления производных, за исключением первой. Приближение получается ценой нескольких вычислений функции. Классический метод Рунге — Кутта четвертого порядка выражается формулой У,+4=-У. + — (64+264+266+66)* где 6,=6)(у„, („), ~,-ч(у.+-,~„~.
—,), 6,=6~ (у„+ —,6„(„+ —,1, 6, =-6((у„+6„4„+6). Заметьте, что на каждом шаге требуется четыре вычисления функции ((у, 4). Классический метод Рунге — Кутта может рассматриваться как обобщение на дифференциальные уравнения квадратуриой формулы Симпсона: Если г ~у, 1) есть функция только от й то обе формулы совпада4отс При 6 0 классический метод Рунге — Кутта четвертого порядка согласуется с рядом Тейлора асимптотически до членов порядка 6'.
Однако здесь нет легко определяемой оценки для локальной ошибки дискретизации, помогающей в выборе величины шага. В З 6.8 будет подробно описана современная модификация метода Рунге — Кутта, которая включает контроль длины шага. Методы Рунге — Кутта имеют несколько достоинств. Их легко программировать; опи численно устойчивы для широкого класса задач. Поскольку для вычисления у„„нужно лишь одно значение решения у„, то метод является гал4остарту4ощил.
Величину шага 6„можно изменить на любом этапе вычислений. б ЗАДАЧА КОШИ 138 Многошаговые методы В методах, рассматривавшихся до снх пор, значение у„„. вычисляется посредством функции, которая зависит лишь от 1„, у, и длины шага и„. Видимо, логично предположить, что можно получить ббльшую точность, используя информацию в предыдущих точках, именно у„„у„„...
и )„ь 1„„..., где ),= =)(уи Г,). Многошаговые методы, основанные на этой идее, весьма эффективны. Если требуется высокая точность, то они обычно более экономичны, чем одношаговые методы, и часто можно тривиально получить оценку ошибки усечения. Запрограммированные соответствующим образом, многошаговые методы могут эффективно выдавать значения численного решения в произвольных точках, не изменяя значение >и Порядок метода может выбираться автоматически и динамически изменяться, тем самым получаются методы, работающие для очень широкого круга задач. В процессе решения могут выявляться некоторые виды ошибок.
Жесткие уравнения (обсуждаемые ниже) могут решаться некоторыми многошаговыми методами, и уравнения можно автоматически классифицировать как жесткие нли нежесткие. Эти достоинства приобретаются ценой усложнения программы н в некоторых случаях за счет возможной численной неустойчивости. Линейные многошаговые методы можно рассматривать как специальные случаи формулы У,«> = Х сб'У«б>->+и Х )>>)«+>-г ~=1 !=б где й — фиксированное целое число и ад либо ))А отлично от нуля.
Эта формула определяет общий линейный и-и>вговый метод. Метод называется линейным, потому что каждое г'; входит в формулу линейно; при этом г может быть, а может не быть линейной функцией своих аргументов. После того как метод «стартовал», каждый шаг требует вы числения у«„ из известных значений: у„, у„ о..., у,-д«>, >, Абь ЕСЛИ и«= — О, МЕтОд НаЗЫВаЕтСя ЯВНЫМ И ВЫЧИСЛЕ- ние проводится очевидным образом. Если же ~,~0, то метод называется неявным, потому что для нахождения у„„требуется Ббльшие трудности при использовании неявных методов окупаются другими их качествами. Обычно на каждом шаге решения совместно используются два многошаговых метода.
Явный метод, называемый предиктором, сопровождается одним илн более применениями неявного метода, называемого корректором — отсюда название «методы предиктор-корректор. » 6. с метОды 439 Примером хорошего метода предиктор-корректор является метод Адамса четвертого порядка; предиктор: у„4, — — - у, +,— (551', — 591„, +37)„4 — 91„4), корректор: д„„=- у„+ — (91,, + 1914 — 51„, +1„,). Обе формулы имеют четвертый порядок. Это частные примеры семейства предикторов и корректоров с различными порядками, Схема предиктор-корректор для вычисления 9„44 такова: 1. Использовать предиктор для вычисления 93, — началь ного приближения к у„е4.
Положить 1=0. 2. Вычислить функцию и положить 1"14=1(у„",'4. 1, 4) 3. Вычислить лучшее приближение у„",")' по формуле корректора, полагая 1,4,=1',",414. 4. Если 1у'„'„" — 9,",,',~) заданного допуска, то увеличить 4 на 1 и перейти к шагу 2; в противном случае положить у„,=-у'„'О'. Здесь три различных процесса: шаг-предсказание 1, который мы обозначим через Р; шаг-вычисление 2, назовем его Е; и шагкоррекция 3, обозначаемый С.
Шаг 4 может быть заменен требованием выполнения ровно т итераций корректора. Получающийся метод обозначают Р (ЕС)'". Имеются веские причины полагать, что разумно заключительное вычисление 7, дающее более точное значение для следующего шага; подобные методы обозначают Р(ЕС) Е. Если порядок предиктора не меныпе, чем порядок корректора, или меньше только на единицу, то более чем двукратная коррекция уже существенно не влияет на результат. Для методов Адамса известно, что если принимать во внимание и устойчивость, и ошибки дискретизации, то метод Р(ЕС)'" уступает методу Р(ЕС) -" Е, а метод Р(ЕС)'Е более устойчив, чем Р(ЕС)'Е.
Таким образом, наиболее разумной схемой предиктор-корректор для методов Адамса является, видимо, РЕСЕ. Методы экстраполяции Формулы предсказаняя, обсуждав4ш4еся в предыдущем разделе, можно рассматривать как методы экстраполяции, где 4ь,„ экстраполируется по известным предыдущим значениям у4 и ~,. Метод, использующий другой тип экстраполяции, предложил в 1927 г. Гонт (Л. А.
С4ацп1); его развил в докторской диссертации 1964 г. Грэг (\Ч. В. С4гаяй). Мы вкратце опишем этот численный метод. По поводу некоторых тонких деталей отсылаем читателя к известной статье Булирш, Стоэр (1966). !4а 6 зкдлчл каши Дифференциальное уравнение интегрируется на заданном интервале, например 0(((Т, посредством одношагового метода для различных длин шага и;; полученные результаты обозначим через У'(1»!).
Это дает дискретное приближение к функции 1'(й), где 1'(О)=-у(Т). По вычисленным узлам строится интерпаляционный полинам пли рациональная функция г'(л). Тогда экстраполированное значение 1'(О) является обычно очень хорошим приближением к 1'(О). Эта еще один пример упомянутого в гл. 5 приема Ричардсона «экстраполяции к пределу».
б.а. Жесткне уравнЕния Постоянная времени для решения дифференциального уравнения — эта время, необходимое для ега уменьшения, выражаемого коэффициентом 11е. Например, уравнение у'== — Ху имеет решение Се ". Если Х положительно, та у уменьшится в !/е раз за время Ц.. Устойчивое дифференциальное уравнение называется жестким, если она имеет частное решение в виде убывающей экспоненты, постоянная времени которого очень мала в сравнении с длиной интервала, иа котором разыскивается решение. Для системы у'=-1(у) скорости убывания связаны, па крайней мере локальна, с собственными значениями матрицы д~/ду=-(дну!). Рассмотрим, например, систему и' =- 998и + 1998о, о' == — 999и — 1999о. Собственные значения матрицы коэффициентов 998 1998~ †9 †19) суть — 1 и †10.