Ортега Дж., Пул У. Введение в численные методы решения дифференциальных уравнений. Под ред. А.А.Абрамова (1986) (1095855), страница 32
Текст из файла (страница 32)
например, д,/; 1 (х) = — ~Ях,,...,х,, х, +й,х...,...,хн) — Д(х)]. (4314) а., ) При этом подходе требуются выражения только для самих функций,Гт, 9. дж. Ортега 129 Таблица 42 Сходнмость метода Ньютона в случае двумерной снстемы х', + х,* — 1 = О, х — х =0 а 1 2 0,5 0,87499999 0,79067460 0,78616432 0,78615138 0,0 1,2 2,4 4„8 8,8+ 0,5 0,62499999 0,61805555 0,61803399 0,61803399 которые нужны в любом случае. Но фактическое вычисление матрицы Якоби либо с помощью выражений для частных производных, либо через аппроксимации типа (4.3.14) может оказаться дорогостоящим с точки зрения машинного времени. Поэтому на практике часто используют модификацию метода Ньютона, при которой матрица Якоби пересчитывается лишь время от времени, а не на каждой итерации.
Такой алгоритм может выглядеть следующим образом. 1. Вычислить Е (хе). 2. Провести итерации хг +' =х' — (Г'(хе)] ' Г(х'), ~' = О, 1, й. 3. Вычислить Г (х +'). (4.3.15) 4. Провести итерации х'+' =х' — [Е'( х ')] 'Р(х'), ~' = й + 1,..., 2/с. 130 Модифицированные итерации Ньютона могут быть также полезны, чтобы ослабить второй недостаток метода Ньютона, связанный с необходимостью решения на каждой итерации системы линейных уравнений.
Преимущество алгоритма (4.3.15) в этом отношении заключается в том, что при фактической его реализации приходится решать ряд линейных систем (как на шаге 2) Р'(хе)у' =Г(х'), 1'= 0,1,... „Ус, где матрица коэффициентов Р'(хс) не изменяется. Как указывалось в разделе 3.3, это позволяет запомнить в ходе гауссова исключения множители Е Уразложения матрицы,Г (хс) и использовать их для всех /с+ 1 правых частей.
Третьей и наиболее серьезной трудностью, возникающей при использовании метода Ныатона, является то, что при заданном начальном приближении хе итерации могут расходиться; локальная теорема сходимости гарантирует сходнмость только в том случае, если точка хе (или какая- либо другая итерация) окажется "достаточно близкой" к х'. Одно из средств преодоления этой трудности заключается в построении наилучшего возможного приближения исходя из физических (или каких-либо инь1х) соображений.
Однако это не всегда бывает достаточным. Более менее систематическим подходом, который часто, но отнюдь не во всех случаях приводит к успеху, является метод продолжения по параметру, к изложению которого мы теперь переходим. Во многих научных задачах требуется решить определенные уравнения при различных значениях одного или нескольких параметров. Предположим, что в задаче имеется один параметр а, и запишем систему уравнений в виде г (х; а) = О. (4.3.16) Допустим, нам надо найти решения хг' при значениях ао < а, «... ад, причем прн а = ао задача решается тривиально или достаточно просто; например, при а = а, система уравнений может оказаться линейной. Если хс может быть вычислено и величина ] а, — ао ] мала, есть основания надеяться, что хс' достаточно близко к х*,, чтобы быть подходящим начальным приближением для решения уравнения Г(х; а,) = О.
Продолжая этот процесс, используем каждое найденчое решение предыдущей задачи в качестве начального приближения для решения следующей задачи. Если уравнения, которые нужно решить, не содержат такого параметра, мы всегда можем ввести его искусственно. Пусть, например, Г(х) =0— система уравнений, подлежащая решению, и х — наша лучшая аппроксимация решения (но недостаточная для сходимости метода Ньютона) . Введем новую систему уравнений, зависящую от параметра а: Г(х; а) = Г(х) + (а — 1) Г(х') = О, 0 < а < 1. (4.3.17) л Тогда при а = 0 имеем систему Г(х; О) = Г (х) — Р (хо) с известньгм решением хо, а лрн а = 1 получаем исходную систему Г (х; 1) = Г (х) = О.
Следовательно, мы можем действовать, как в предыдущем абзаце для последовательности О = ас < а, «... агу = 1, Дополнительные замечания и ссылки 4.3 Более подробное обсуждение богатого многообразия методов численного решения систем нелинейных уравнений имеется в книгах [54, 55]. В частности, в этих книгах описываются различные дискретные формы метода Ньютона, в которых матрица Якоби каким-либо образом аппроксимируется. Некоторые из таких аппроксимаций приводят к естественным обобщениям метода секущнх на системы уравненгй, другие ведут к так называемым квазиньютоновским методам, которые в настоящее время выглядят наиболее многообещающими.
Обзор состояния этих квазиныотоновских методов можно найти в статье [18] . Многие системы уравнений возникают при минимизации или максимизации функций и переменных. Из анализа известно, что если функция я непрерывно дифференцируема, то необходимое условие локального минимума имеет вид т,е.
градиент функции я должен обратиться в нуль. Следовательно, решив эту систему уравнений„получим точку возможного локального минимума функции я, причем часто бывает известно, что в такой точке минимум действительно доспп ается. Обратно, если нам дана система уравнений ); = О (г = 1...
гг), то, вводя функцию л е(х) = Х [Ях)] ', г=! можно свести решение системы к задаче минимизации функции я. Ясно, что минимум функции ь, равный нулю, достигается тогда, когда все ~г обрацгюотся в нуль. Для 131 получения численного решения такое преобразование обычно не рекомендуется, так как при этом ухудшается обусловленность задачи. С методом продолжения по параметру тесно связан метод Лавиденко. Для его иллюстрации мы рассмотрим систему (4.3.17) и предположим, что при а Е [ О, 1[ она определяет непрерывно днфференцируемое по а решение х(и).
Дифференцируя тогда тождество Е(х(э)) +1а — 1) Е(х') = О по а, по правилу дифференцирования сложной функции получаем Е'(х(е))х'(а) + Е(хз) = О. Отсюда, предполагая невырожденность матрицы Якоби Г'(х(а)), приходим к дифференциальному уравнению х (а) = — [Г (хба))[ ' Е(хз) с начальным условием х(О) =х'. Решение х(а) этой задачи Коши при а =1 даст, как мы надеемся, желаемое решение исходной системы уравнений Е(х) = О. На практике нам придется решать дифференциальные уравнения численно, и для этого можно, в принципе, использовать любой из методов гл. 2.
Хотя метод Давиденко и метод продолжения по параметру выглядят весьма привлекательно, в практических задачах они нередко оказываются не столь надежными, как хотелось бы. Так, в частности, возможно, что при некотором зз ( 1 матрица Якоби Г'(х(зз)) станет вырожденной, или, например, кривая решения уйдет в бесконечность, не достигнув значения а = 1. За последние несколько лет были достигнуты значительные успехи в преодолении некоторых нз этих трудностей; обзор недавних результатов смотрите в работе [51].
УЛРА ЖНЕНИЯ 4.3 4.3.1. Покажите графически, что система уравнений хз' + х,' = 1, к,* — х, = О имеет ровно два решения. 4.3.2. Покажите, что если Е имеет вид (4.3.3) и В = А ', то (4.3.5) сводится к (4.3.4). 4.3.3. Вычислите матрицу Якоби С' (х) для вектор-функции хз +х!хзхз +хз С(х)- х, +х,х, ~! з 4.3.4. Покажите, что если С (х) =х — ВГ(х), то С' (х) =1 — ВГ'(хз, и выведите отсюда, что условия (4.3.9) и (4.3.10) следуютиз (4.3.8). 4.3,5. Для функций из упражнения 4.3.1 найдите уравнения касательных плоскостей в точкех, = 2,х, = 2. 4.3.6. Выпишите формулы итераций Ньютона для уравнений из упражнения 4.3.1. В каких точках х матрица Якоби не вырождена? 4.3.1.
Составьте программу для решения системы и уравнений с л неизвестными по методу Ньютона. Для решения линейных уравнений используйте гауссово исключение с частичным упорядочиванием. 4.4. Конечно-разностные методы для нелинейных краевых задач Теперь займемся распространением на нелинейные краевые задачи метода конечных разностей, рассмотренного в гл. 3. Сначала обратимся к до вольно простому уравнению и =я(х, и), ()(хз..1, 14.4.1) я(х,, и,) 0 Н(~с) = 6 (4.4.4) 0 я(х„, оя) !! Тогда система (4.4.3) примет вид Е(ч) — = Ат+ Н(~) = О. (4.4.5) В качестве примера рассмотрим задачу и (х) = 3 и (х) + х а + 10 ~ о (х) ~ ~, 0<х <1, и(0) = и(1) = О. (4,4.6) Здесь ,е(х, и) =Зи+х' + 10 о', (4.4.7) и при х,=~И, 1 — О, !»...
)7+1> (4.4.8) разностные уравнения (4.4.3) имеют вид 2о + „йг(3„+ Р~,~ + 10оз) ! ! = 1, ..., л, (4.4.9) где из ганичных условий о~ = и„~ ~ =О. Для этой задачи матрица А и 133 с граничными условиями и(О) =а, и(1) =Д. (4.4.2) Здесь х — заданная функция двух переменных, а а и !1 — заданные константы.