Введение в вычислительную физику - Федоренко Р.П.
ISBN 5-7417-0002-0
Скачать (прямая ссылка):
ночного соотношения х = а у + 0, получаем уравнение для у:
у= bay+ fib+ у, у(1)=0.
Проинтегрируем задачу справа налево, попутно определяя x(t) =
= а(*)у(<) + Р(0*
Перейдем к анализу метода прогонки, рассматривая для общности краевые условия вида Лх(0) + By(O) = C. В этом случае
а(0) = -В/А, 0(0) = С/А. Разбе-
ремся в том, действительно ли «большой» параметр (который так и остался в задаче) уже не страшен и процесс вычислений устойчив. Нам нужны некоторые оценки для a(t). Ограничимся физически наиболее естественными условиями при 1 = 0:
А>0, 5^0, т.е. а(0) > 0.
Рассмотрим поле направлений a. Ha плоскости ({, а) введем кривую (рис. 9)
а = 0, т.е. a(t) = y/b(t)/a(t) = 0(1).
При а = 0, очевидно, а > 0. Выделим области а > 0 (ниже кривой
а(()) и а< 0 (выше кривой a(t)). Несложный анализ показывает,
что
О S a(t) ^ max a (I) = 0(1). е
Этого нам достаточно для дальнейшего.
Посмотрим, что дает теория численного интегрирования, примененная к уравнению а = —Ьа2 + а. Если бы мы оценивали устойчивость численного интегрирования для этой системы по самой общей теореме, оперирующей с оценкой погрешности типа Zect (где С — константа Липшица правой части, є — локальная погрешность), то картина была бы пессимистической. В самом деле,
92
ОСНОВЫ ВЫЧИСЛИТЕЛЬНОЙ МАТЕМАТИКИ
[Ч.І
и все трудности были бы такими же, как и при методе фундаментальных решений. Ho
(—Ьа2 + а) = —2Ьа < О,
т.е. мы получили устойчивое решение, для которого специальная теорема о точности численного интегрирования не содержит экспоненциального множителя, и для шага т достаточно иметь только соотношение
(ba2) 1, т.е. 80т«1.
Итак, нас выручает устойчивость искомого решения a(t). В задаче для р та же самая ситуация:
^jj(—бар + а<р + /) = —Ьа « —40,
т.е. мы имеем дело с интегрированием устойчивой задачи.
Наконец, обратная прогонка. Ее уравнение имеет вид
у = bay+ (fib + <р), (bay)y = Ьа « 40 > 0.
Эта задача неустойчива вправо и, соответственно, устойчива влево. Ho ведь нам нужно интегрировать ее именно справа-налево! И здесь все в порядке, несмотря на присутствие большого параметра.
Заметим, что прогонку можно осуществить в обратном направлении — решая уравнение для а справа-налево. В этом случае (см. рис. 9) траектория а (г) «притягивается» к кривой а = = —y/b(t)/a(t) (a(t) < 0) и интегрируется устойчивая влево задача.
§ 10. Прогонка в разностной задаче Штурма-Лиувилля
Рассмотрим классическую краевую задачу Штурма—Лиувилля:
+ + r(t)xd) = /(с).
dt dt с краевыми условиями общего вида:
ах + рх =7, t = 0, O1X+ P1X = Y1, t = Т.
Начнем с построения разностной схемы, т.е. разностной аппроксимации задачи. Введем сетку, для простоты равномерную:
W*=o’ = пх> X = T/N, TV»I,
и счетные величины х„, п = 0, 1, ..., N.
§ 10] ПРОГОНКА В РАЗНОСТНОЙ ЗАДАЧЕ ШТУРМА-ЛИУВИЛЛЯ 93
Построим разностное уравнение, аппроксимирующее дифференциальное:
Pn+1/2 т Pn-1/2'
-L- Q І!±!_____Х"~\ -(- г х = /
' Чп 2т гпхп Jn>
где рп+1/2 = /?(/„ + т/2), Qn = ?(/„) и т.д. Это уравнение можно написать только при я *= I, 2, ...,Ar - 1 (во внутренних узлах сетки). Для дальнейшего уравнениям удобно придать стандартную форму:
апхп —I ^пхп + ^nxn+l ^п,
где ап, An, cn, dn — так называемые локальные коэффициенты схемы. Имеем для них выражения
— _L _I — ± J.I
2 Рп-1/2 х ^n' т2 ^п+1/2 т
ап + гп» Ar
Примем довольно естественные с физической точки зрения условия
р > 0, г < 0 и отметим важное соотношение An ^ а„ + сп. Аппроксимируем левое краевое условие:
дс, —ДCn а -7- + Px0 = 7,
и запишем его в стандартной форме:
^охо + coXi = — — Р, C0 = —, d0 = у.
Ради простоты ограничимся физически наиболее естественными условиями а > О, P < 0; следовательно A0 > с0.
Аппроксимируем правое краевое условие:
XN-XN-l , „
«1 т + PiXw —
и запишем его в стандартной форме:
aNxN-I ~ ^NxN = ^N-
Итак, мы получили специальную, но очень распространенную в приложениях систему линейных алгебраических уравнений:
94
ОСНОВЫ ВЫЧИСЛИТЕЛЬНОЙ МАТЕМАТИКИ
[Ч. I
Матрица системы имеет так называемую трехдиагональную (якобиеву) форму:
—b0 . C0 О
«і ~Ь\ cI
а2 ^2 С2
ап -ьп Сп
0 aN-I aN cN -bN
Такие матрицы часто появляются при аппроксимации дифференциальных уравнений разностными. Их специфика — большой порядок (N=Tfx) и огромное число нулей (так как операторы дифференцирования являются операторами локального типа: значение производной функции в какой-то точке зависит только от значений функции в сколь угодно малой окрестности этой точки).
Большую роль в вычислительной математике играют так называемые экономные методы решения подобных систем уравнений. Это такие методы, в которых количество операций пропорционально первой степени числа неизвестных, т.е. в данном случае O(N). Напомним, что если бы мы просто сослались на то, что получена система линейных алгебраических уравнений, которую можно решать любой стандартной программой, дело было бы довольно скверным. В общем случае решение системы N алгебраических уравнений с N неизвестными требует O(Ni) операций и 0(N2) ячеек памяти.
