Введение в вычислительную физику - Федоренко Р.П.
ISBN 5-7417-0002-0
Скачать (прямая ссылка):
Как известно, точное общее решение рассматриваемой системы уравнений с постоянными коэффициентами имеет вид
где Xi, Yi — некоторые числа (которые легко ВЫЧИСЛЯЮТСЯ, HO они нам не нужны), C1, C2 — произвольные постоянные, a X1, X2 — корни характеристического уравнения
а) х‘ = аУ, X1(O) = I, б) х2 = ау2, X2(O)=O,
у1 = Ьх1, у1 (0) =0; у2 = Ьх2, уг(0) = 1,
ищем решение в виде линейной комбинации:
—X а
§9]
МЕТОД ДИФФЕРЕНЦИАЛЬНОЙ ПРОГОНКИ
89
Таким образом, почти любое частное решение (в том числе построенные нами два линейно независимых решения) есть сумма (с примерно равными коэффициентами C1, C2) двух экспонент: одной —
сильно растущей (типа e40t), второй — сильно убывающей (типа е~ш).
Теперь обратимся к исходной задаче. Прежде всего подчеркнем, что выбор больших значений коэффициентов а и b не был произвольным, он продиктован практикой. Наиболее близким содержательным примером задачи, качественные черты которой хорошо передает разбираемая модельная, является прохождение излучения через слой большой оптической толщины, например прохождение потока нейтронов, источником которого является ядерный реактор, через слой защиты. В этом случае одно краевое условие х0 = Sf0 задает поток нейтронов, падающий на внутреннюю поверхность защиты, второе условие у(1)=0 означает отсутствие потока нейтронов, падающих на внешнюю границу защиты. Интересующая же нас величина х(1) имеет смысл потока нейтронов, выходящего со стороны внешней границы защиты.
Искомое решение есть функция типа e~40tS?0 (в этом, собственно, и состоит назначение защиты: ослабить поток нейтронов примерно в е40«1016 раз). Мы же пытаемся получить его в виде линейной комбинации двух линейно независимых решений, в каждом из которых решающую роль играет именно растущая экспонента. Получить функцию типа е~4Ы в виде линейной комбинации решений, в которых главную роль играют компоненты типа e40t (они должны взаимно погаситься), — очень трудная вычислительная задача, сопровождающаяся резким падением точности.
Следует принять во внимание и накопление погрешностей вычислений. В рассматриваемом случае, допустим, при интегрировании задач Коши методом к-ю порядка погрешность аппроксимации у левого конца траектории имеет величину порядка тк3?0 (х — шаг численного интегрирования). Ее последствия у правого конца траектории достигнут величины е40т*,Г0, и нужно, чтобы они были заметно меньше искомого решения, т.е. величины порядка e~40Sf0.
Итак, имеем ориентировочное соотношение для шага численного интегрирования:
e40xk<SCe~40, Т.е. Tt-tKlO-30.
Даже при к = 5 получаем х < IO-6, т.е. нужно брать сетку с миллионом узлов. Это чудовищно. Забегая вперед, укажем, что
90
ОСНОВЫ ВЫЧИСЛИТЕЛЬНОЙ МАТЕМАТИКИ
[Ч. I
фактически такие задачи решаются методом прогонки при условиях IX1 2x|<scl. В нашем случае х« 1/40, т.е. вполне приемлем, например, расчет при т « 1/200 и даже при х « 1/100.
Убедившись в провале «школьного» метода решения простой краевой задачи, приступим к описанию метода прогонки. Отметим только в качестве «морали»: в вычислительной математике важна не столько внешняя форма задачи, сколько качественные свойства искомых решений. Абсолютно одинаковые внешне задачи часто требуют существенно разных методов. Ta же самая задача при a fa Ъ s» 5 (или на интервале 0 ^ г =S 0.1) без всяких затруднений может быть решена только что скомпрометированным методом фундаментальных решений.
Рассмотрим детально метод дифференциальной прогонки. Будем искать связь между компонентами решения вида
х(() = а(Є)у(0 + Р(0.
где а(/), (1(1) — неизвестные пока «прогоночные коэффициенты». Получим уравнения для них. Дифференцируя прогоночное соотношение:
х = ау + а у + Э,
и учитывая, что х = ау •+¦ /, у = Ъх + ip, имеем
ay + f = ay -1- a(bx + <р) + р.
Заменяя х на ау 0:
ay + f = а + Ьа2у + Ь$а + а<р 4- (3,
и приводя подобные члены (коэффициенты при у и единице), полу-
чаем
у [а + Ьа2 — а] + IP + bfia — /] = 0.
Приравнивая нулю коэффициенты при у и единице, приходим к
уравнениям для а и fk
а + Ьа2 — а — 0, (3 + а$Ь + а<р — / = 0.
Эти уравнения дополним начальными данными, используя стандартный прием метода прогонки. Левое краевое условие х(0) = Sf0
запишем в виде того же самого прогоночного соотношения:
х(0) = а(0)у(0) -I- 0(0). Очевидно, следует положить а(0) = 0,
P(O)¦ = Sf0. Итак, получены задачи Коши для a(t) и Р(г)- Они могут быть проинтегрированы (например, численно), и можно считать, что функции a(t) и f}(t) у нас уже есть.
Перейдем к следующему характерному элементу прогонки — разрешению правого краевого условия. Имея условие у( 1) = 0 и
9]
МЕТОД ДИФФЕРЕНЦИАЛЬНОЙ ПРОГОНКИ
91
Рис. 9
прогоночное соотношение при t = 1: х(1) = а(1)>’(1) + 0(1), легко находим значение лг(I) = р( 1).
Наконец, рассмотрим заключительный этап прогонки. Опять-таки отклоним напрашивающийся рецепт: раз мы знаем х(1) и >>(1), можно (формально) интегрировать задачу Коши справа налево. Ho эта задача так же неустойчива, как и задача Коши, решаемая слева направо. Мы воспользуемся уравнением у = Ьх 4- <р. Заменяя х из прого-
