Введение в вычислительную физику - Федоренко Р.П.
ISBN 5-7417-0002-0
Скачать (прямая ссылка):
y(t, h) = x(t) + h bx(t),
bx = fx[x{t)\ bx(l) — r(t).
86
ОСНОВЫ ВЫЧИСЛИТЕЛЬНОЙ МАТЕМАТИКИ
[Ч. I
Линеаризация краевых условий строится очевидным образом. Пояснения заслуживает способ вычисления матрицы fx[t\ = fx(x(t)). В расчетной схеме в качестве этой матрицы использовалась кусочно-постоянная матрица
Решение линейной краевой задачи осуществляется методом ортогональной прогонки, описанным в § 18. He следует думать, что шаг численного интегрирования совпадает с шагом сетки т, он в несколько раз меньше. В процессе интегрирования запоминается ограничение функции bx(t) на сетку {tn} — сеточная функция {6х„}.
Определение шага. Определяется функция R(K) — норма невязки сеточной функции [хп + h находится (не очень точно)
min R(h) по А. После определения h новое приближение вычисляется по формуле
хп •¦= хп + h 6хп, /1 = 0, 1, ..., N.
Выше описан стандартный шаг итерационного процесса. Итерации продолжаются до получения достаточно малой невязки. Отметим, что во многих прикладных задачах одной из наиболее трудоемких операций является линеаризация, т.е. вычисление f х. Конечно, такая операция трудна при достаточно сложной форме правой части. В рассматриваемом примере это не так. Трудоемкой операцией, вообще говоря, является и вычисление шага h, требующее нескольких вычислений невязки. Естественно, число операций пропорционально числу интервалов сетки N, и следующее ниже усовершенствование имеет целью повысить точность, не увеличивая существенно N или даже совсем его не меняя. Это достигается изменением определения невязки.
Уточненная формула вычисления невязки. Имея сеточную функцию {*„}, восполняем ее до непрерывной функции X(t) с помощью кусочно-гладкого интерполяционного аппарата (с помощью сплайна, например). Теперь, в принципе, можно говорить
о непрерывной функции r(t) = х — f(x(t)) и вычислять ее норму
(J г2 dt)112. Разумеется, норма вычисляется по какой-нибудь хорошей квадратурной формуле, так что на самом деле такая «непрерывная» невязка вычисляется в дискретном наборе узлов квадратуры (в описываемых ниже расчетах на каждом интервале (tn, tn+j) интеграл вычислялся по квадратуре Гаусса с тремя узлами). Остальные элементы вычислительной схемы остаются без изменений.
§8]
РЕШЕНИЕ КРАЕВЫХ ЗАДАЧ ДЛЯ СИСТЕМ ОДУ
87
Приведем результаты, характеризующие эффективность алгоритма и достигнутую точность.
В табл. 7 представлены следующие данные: і — номер итерации, R — норма невязки, h — шаг на ї-й итерации. Начальное приближение выбрано таким, что краевые условия выполнены; они линейны, поэтому не нарушаются в процессе итераций. Невязка со-
Таблица 7
/ 0 I 2 3 4 5
R 15.60 14.15 13.73 12.40 11.36 10.63
h - 0.155 0.053 0.153 0.130 0.115
і 6 7 8 9 10 11
R 9.31 6.84 1.43 0.20 0.007 0.0002
h 0.197 0.413 1.04 1.03 1.0 1.0
держит только компоненту х — /; невязки в краевых условиях остаются нулевыми. Расчет проводился при шаге основной сетки х = г/50 = 0.4.
Таблица 8 характеризует точность расчета. В ней приведены значения компонент вектора х в нескольких характерных точках
Таблица 8
дс1(10.4) л1 (20) х2(5. 6) дс2(12) X3(O) л3(10.4)
а) -1.0808 -1.1808 0.1244 -0.0455 -0.9621 -0.0438
б) -1.091 -1.190 0.1214 -0.0444 -0.9863 -0.04308
а) -1.091 -1.190 0.1213 -0.0443 -0.9863 -0.04312
дс4(5.6) х4(12) дс5(0) дс5<5.6) дс5(12) Jt5 (20)
а) 1.1838 0.9867 0.6442 —0.09101 0.02496 0.00915
б) 1.1817 0.9867 0.6529 -0.08878 0.02434 0.008806
а) 1.1820 0.9866 0.6529 -0.08874 0.02433 0.008622
времени. Каждая величина представлена тремя числами: а) результат, полученный за одиннадцать итераций при грубом вычислении невязки; б) результат, полученный после одной дополнительной итерации с более точным вычислением невязки,
в) «точное» решение.
88
ОСНОВЫ ВЫЧИСЛИТЕЛЬНОЙ МАТЕМАТИКИ
[Ч.І
§ 9. Метод дифференциальной прогонки
Начнем изучение самого, вероятно, популярного в вычислительной математике алгоритма. С ним связаны существенные успехи в решении сложных задач, и его изобретение считается по праву ярким событием в истории развития современной вычислительной математики. Любопытно, что метод прогонки предназначен для решения задачи, на первый взгляд не содержащей никаких проблем.
Пусть требуется решить краевую задачу для системы двух уравнений
X=^ay + /, у—Ьх + у, 0<К1, а> О, Ь > О,
с краевыми условиями х(0) = T0, у(1) = 0. Коэффициенты а и Ъ ради простоты считаем постоянными. Очень важны их числовые значения. Пусть а & Ь 40. Эти значения взяты не случайно, в них суть дела. Как мы увидим дальше, а и Ь — «большие параметры». Именно то, что они большие, определяет специфику задачи и требует разработки принципиально новых алгоритмов.
Объясним причины, по которым описанный в § 8 метод сведения краевой задачи к задачам Коши не работает. Проанализируем известный («школьный») метод. Ради простоты положим / = 9 = 0 (не в них дело). Найдя решение двух задач Коши:
а коэффициенты O1 и а2 определим, подставив его в краевые условия. Что же из этого получится? Посмотрим, какие последствия имеют большие значения коэффициентов а и Ь (и в самом ли деле они большие).
