Федоренко Р.П. Введение в вычислительную физику (Федоренко Р.П. Введение в вычислительную физику.djvu), страница 17

При вычислении линейно-независимых решений х'(1) используем для них напрашивающиеся данные Коши: х'(0) = е', 1= 1, 2, ..., р, где е' = [О, ..., О, 1о О, ..., 0), т.е. 1-й орт р-мерного пространства. Итак, хо = е', и далее х~+~ х«~ + тВАьхь~, 8 = О, 1, „., К вЂ” 1. Отметим, что такой способ решения краевой задачи «стоит» (р + 1)-кратного решения задачи Коши. Однако часто объем работы можно сократить. Это относится к очень распространенному типу краевых задач: г < р компонент х задано при г = 0 и р — г компо- 81 гешеиие кгяевых зАдАч для систем одг 881 мент — при ! = Т, т.е.

краевые условия имеют вид х,(0) =у!, х (0) =/м ..., х„(0) =У, (здесь нижний индекс — номер компоненты). Правые краевые ус- ловия произвольные; например, Рх(Т) = Ь или Вх(0)+Рх(Т) = а!, где В, Р— прямоугольные матрицы р- (р — г) (р столбцов, р — г строк), Ь вЂ” (р — г)-вектор.

В этом случае для решения хв(1) берем данные Коши: х(в!=0 т=г+1 г+2 р хв(0) =Ум !'=1, 2, ..., г; а решение краевой задачи ищем в виде х(!) = х~(1) + ", а! х!(!) ! г+! (х1(!) находятся так же, как и раньше). легко видеть, что такая конструкция при любых а! удовлетворяет уравнению х = Ах + а и левым краевым условиям, а свободных параметров а! как раз столько, чтобы за их счет выполнить р — г условий на правом конце интервала времени. Нелинейные краевые задачи.

Метод «стрельбы». Перейдем теперь к нелинейным краевым задачам. Как всегда, в нелинейной ситуации лучше говорить о возможном подходе, чем о методе. Итак, пусть требуется найти решение уравнения — ",=Т(х,г), Омтмт, Когда мы говорим «определим траектории!», это означает, что при каждом заданном значении вектора параметров а мы можем с какой-то точностью численно проинтегрировать задачу Коши. Введем функцию г"(а) к«Ф(а, х(Т, а)). Решение краевой задачи свелось к решению системы нелинейных уравнений р(а) = 0 (р уравнений с р неизвестными). Еще раз под- при общих, например, краевых условиях Ф(х(0), х(Т)) = О.

Используем умение достаточно надежно решать задачу Коши. Введем данные Коши х(0) в качестве искомых неизвестных. ОбозначаЯ их чеРез а= (а,, ам ..., а ), опРеделим тРаектоРию х(0 а) задачи Коши: х Лх~ !)~ х(0) = а. основы вычислитвльной млтвмлтики зз 1ч.| черкнем, что функция Р задана нам достаточно сложным алгоритмом, позволяющим для любого а вычислить вектор Р; такое вычисление «стоит» одного численного интегрирования задачи Коши с начальными данными а. Решение системы можно осуществить методом Ньютона и его модификациями, Конечно, в этом случае вычисление матрицы Р,(а) проще всего выполнить численным дифференцированием, хотя есть и более аккуратные методы, используемые в вариационном исчислении (они предполагают нсиользование так называемой системы уравнений в вариациях; см.

з 27„ 28). Спектральная задача Штурма-Лиувилля. Специальный, но очень важный класс краевых задач связан с определением точек спектра для уравнения Штурма — Лиувилля. Рассмотрим простейший случай. Задано линейное однородное дифференциальное (само- сопряженное) уравнение ф р(1) ф +о(Г) х(Г) = ьг(Г) х(г), Оч Гч Т, содержащее параметр Х (функции р(г) > О, д(г), г(г) заданы). Уравнение дополнено простыми краевыми условиями (тоже линейнымн однородными), например х(0) = О, х(Т) =О. При почти всех Х краевая задача имеет тривиальное решение х(1) ив О, но при некоторых специальных значениях )ь, называемых точками спектра, появляются и нетривиальные решения. Они-то (и соответствующие им значения Х) представляют основной интерес в приложениях, Соединим технику решения задачи Коши и решение нелинейных уравнений.

Поставим для уравнения условия Коши х(0) = О, х(0) = 1. (Нетрудно видеть, что вместо х(0) = 1 можно взять в краевом условии любое число.) После этого, если Х задано, определяется траектория х(0 )ь) — решение задачи Коши. Разумеется, при произвольном )ь зта траектория не удовлетворяет второму краевому условию, и теперь надо подобрать Х так, чтобы на траектории х(0 Х) было выполнено второе условие.

Другими словами, определим функцию Р(Х) яв х(Т, Л) (опять-таки, напомним, что вычисление Р при заданном Х требует численного интегрирования задачи Коши) и станем решать уравнение Р(Х) = О. Корни этого уравнения — суть точки спектра задачи Штурма — Лиувилля (разумеется, приближенные, коль скоро функцию Р мы вычисляем лишь приближенно). Самое грубое решение задачи можно представить себе так. Заменив функции р, й, г(г) на постоянные, равные, например, средним значениям по интервалу 1О, Т1 вычислим спектр такой модельной задачи.

Тем самым будет получена ориентировочная информация о 8З 881 гешзние кгяззых зядяч вля систзм оат расположении точек спектра исходной задачи — о расстояниях между ними. Исходя из зтого, выберем некоторый шаг Ь, заметно меньший расстояния между собственными числами, но не слишком малый, н вычислим значения и"(хЬ) в точках сетки хй, х = О, ж1, й2, ...

(Вспомним, что каждое значение г(хй) — зто интегрирование задачи Коши.) ПоР(ХЭ строим «график» г"(Л) по точкам г(хп). Он будет выглядеть примерно так, как показано на рис. 8. Проведя через полученные точки (хотя бы с помощью лекала) гладкую кривую, найдем приблизительные значения корней г(Х) О. Потом, если нужно, нх можно уточнить. Заметим, что все это сравнительно Рис. 8 просто в самосопряженной задаче, когда нам точно известно, что корни Г(Х) = О находятся на вещественной оси.

В общем случае они комплексны н ситуация заметно осложняется (см. 8 15, 16). Решение нелинейных краевых задач. Метод Ньютона. Закончим зтот параграф описанием еще одной популярной алгоритмической конструкции, предназначенной для решения нелинейных краевых задач. Общую идею поясним на следующем примере. Пусть требуется решить краевую задачу для системы уравнений х=У(х, г), Оп1ат, с краевыми условиями хотя бы общего вида Ф(х(О), х(т)) = О (где Ф вЂ” р-вектор). имеется некоторая функция хз(1), не удовлетворяющая ни краевым условиям, нн уравнению. Используем ее в качестве начального приближения и построим алгоритм типа метода Ньютона (в функциональном пространстве). Это обобщение и соответствующая теория разрабатывались Л.

В. Канторовичем в начале сороковых годов. В соответствии с обшей схемой метода Ньютона следующее (первое) приближение ищем в виде х'(г) = хс(С) + Ьх(г), где Ьх(г) — «малая» поправка, В результате получаем уравнение для Ьх: хо+ Ьх — У(хо + Ьх, 1). Линеаризуя его (отбрасывая малые второго порядка), имеем х'+Ьх=У(хс(1),1)+У„( (г), г)Ь .

1Ч.1 ОСНОВЫ ВЫЧИСЛИТЕЛЬНОЙ МАТЕМАТИКИ Чтобы функция х'(Г) удовлетворяла краевым условиям, потребуем выполнения условий Ф[х»(0) + Ьх(0), хО(Т) + Ьх(Т) ] = О. Они тоже линеаризуются: 0 = Ф [хо(0), х~(Т) ] + Ф„<О>[х~(0), хо(Т) ] Ьх(0) + Ф„<г1 Ьх(Т). Итак, Ьх(г) находится решением краевой задачи для системы линейных неоднородных уравнений ~" — У,[хо(Г), Г] Ьх = /[хО(1), 1] — хО(1) с известной правой частью.

Далее процесс повторяется до получения нужной точности, если он сходится (что требует выбора не слишком случайного начального приближения). Ограничимся общим описанием и укажем, что в последнее время этот метод стали называть методом «квазилинеаризацни» (по инициативе Р. Беллмана). Пример решения краевой задачи. Покажем, как фактически реализуется алгоритм. Задача заимствована из американской литературы, где она характеризуется как »неустойчивая». Сложность решения этой задачи существенно зависит от длины интервала [О, Т].

Поэтому ее решение, которое мы будем считать «точным», было получено методом продолжения по параметру, Сначала нашли решение при Т = 10, затем, используя его как начальное приближение, нашли решение при Т = 11.6, далее при Т = 13.2, 2' = 14.8, Т = 16.4, Т = 18.0, Т = 19.6 и,наконец, прн Т = 20, Будем решать задачу модифицированным методом Ньютона в функциональном пространстве сразу на интервале [О, 20]. Сформулируем краевую задачу. На интервале 0 <1И 20 ищется решение х(1) = (х', хз, ..., хз), удовлетворяющее системе уравнений х'=хз, хз= хз, хз = -1.55 х'хз + 0.1 (хз)з — (х4) з + 0.2 хз + 1 х4 х', хз = 1.55 х'хт + 1,1 хзх4 + 0.2 (х4 — 1). Краевые условия; х~=хз=х4 0 4 О.

хз=О, х4=1, г=т= 20. РЕШЕНИЕ КРАЕВЫХ ЗАДАЧ ДЛЯ СИСТЕМ ОДУ 5Е1 В качестве начального приближения берется функция, удовлетворяющая краевым условиям к'(Г) = хз(г) = О, х'(Г) = Г/Т, [0.5г/Т, ' (г) = ]0,5(т — 1)/т, г и 0.5т, г н0.5т г и 0.1Т, 0.1Т и1и 0.5Т, 0.5Ти / и 0.9Т, 1> 0.9Т.' 1/Т, О.

15 — 0.51/Т, 0.51/Т вЂ” 0.35, 1 — //Т, Это приближение очень грубое, оно выбрано без использования известного решения задачи. Теперь поясним детали технической реализации алгоритма. Сетка, сеточная функция. На интервале [О, Т] вводится равномерная сетка (/л)л, Гл = пт, т = Т/М.

В узлах сетки определяется сеточная функция (хл]~ В, х„Е йз. Невязка. Это очень важный объект, определяющий (наряду с числом узлов М) точность приближенного решения. Обозначая систему х = /(х), определяем невязку г(1) как кусочна-постоянную на сетке функцию: Г(Г) / З ' 1 ~ (1л' 1л+1)' Итерационный процесс имеет целью свести норму невязки к нулю, т.е. приводит (при успешных вычислениях) к решению системы разностных уравнений, аппроксимирующих задачу со вторым порядком.

нормой невязки Я считаем (] гз ю) '/т, добавляя еще нормы невязок в краевых условиях. Я, Ь) х(/) + Ь Ьх(г), где Ь вЂ” подлежащий определению шаг. Для определения Ьх(1) решается краевая задача, получающаяся линеаризацией исходной нелинейной задачи на имеющемся уже приближении х(г). Эта краевая задача имеет вид Ьх /,[х(г)] Ьх(~) — г(8). Уравнение в вариациях. Напомним, что процесс итераций по схеме модифицированного метода Ньютона состоит в вычислении поправки Ьх(1) н образовании однопараметрнческого семейства функций основы вычислительной млтемлтнкн зб Линеаризация краевых условий строится очевидным образом.

Пояснения заслуживает способ вычисления матрицы /„(г( вв /,(х(1) ). В расчетной схеме в качестве этой матрицы использовалась кусочно-постоянная матрица у»~ 1 = у» 2 Решение линейной краевой задачи осуществляется методом ортогоналъной прогонки, описанным в 5 18. Не следует думать, что шаг численного интегрирования совпадает с шагом сетки т, он в несколько раз меныпе. В процессе интегрирования запоминается ограничение функции Ьх(г) на сетку (1„) — сеточная функция (Ьх„). Определение итога. Определяется функция й(Ь) — норма невязки сеточной функции (х„+ Ь Ьх„)н„, находится (не очень точно) пнп ЩЬ) по Ь. После определения Ь новое приближение вычисляется по формуле х„:= х„+ Ь Ьх„, и = О, 1, ..., Ф.