Федоренко Р.П. Введение в вычислительную физику (1185915), страница 13
Текст из файла (страница 13)
Но значит ли это, что его решение (х„) в какойто мере аппроксимнрует решение дифференциального уравнения? Так мы приходим к другому фундаментальному понятию — к сходи- мости разностной схемы. Правильнее и аккуратнее было бы сказать: сходимость решения разностного уравнения к решению дифференциального при т - О. Но мы будем использовать «жаргонное» выражение «сходимостьь. Перейдем к аккуратному оформлению этих понятий. Рассмотрим ' следующие математические объекты: 1) дифференциальное уравнение — = Д2', Т), 2'(0) = х„ и его решение 2'(г); 2) ограничение функции 2'(г) на сетку (г„), т.е. сеточную функцию (2'„)„", где 2'„= 2'(г„); 3) разностное уравнение (например, явную схему Эйлера или еще какую-нибудь) и его решение (х„).
Мы будем иметь дело с последовательностью сеток, соответствующих уменьшающимся шагам т - О, так что при стремлении к полной аккуратности следовало бы каждый объект пометить еще индексом т: 2"„, х„* н т.д. Этот индекс т неявно везде следует иметь в виду. 0 п р е д е л е н и е. Говорят, что разностное решение сходится, если Вш ~~2"„— х'„~~ = О, и = О, 1, ..., М=?Ут. «О Если установлена оценка ~!2'„' — х'„~~ ц Стг, и = 1, 2, ..., М, 5 51 6! иитяп иювлиии зядлчи коши для систем одт где С не зависит от т и и, то говорят, что установлен р-й порядок сходимости, а схема имеет р-й порядок точности.
Установление сходимостн, а лучше, и порядка сходимости (еще лучше, с хорошей оценкой константы С) есть основная цель теоретическою обоснования метода приблилсенного решения. Этому вопросу посвящен б 7. Схема Адамса. Опишем общую конструкцию схем численного интегрирования, достоинством которой является ее экономичность. Каждый шаг интегрирования требует только одного вычисления правой части с, в то же время порядок точности метода может быть (формально) любым желаемым. В методах Рунге-Кутты (они описаны в б 7) число вычислений / на шаг равно порядку точности метода.
Итак, пусть задача решается на равномерной сетке, значения х„(и все предшествующие х„„х„з, ..., хр) уже найдены. По значениям ус = у(х„ +,) для с = О, 1,..., р (р определяет порядок точности метода) в узлах сс с„,я+с построим интерполяционный полипом степени рс с-о Ею можно применить и для экстраполяции функции у [х(с) [ на интервале [с„, с„+ т = с„с[. Теперь используем очевидное тождеспю с+с х(с„ + т) = х(с„) + $ у[х(с)[ асс. с (б) Заменяя у[х(с)) интерполяционным полнномом и вычисляя интеграл, получаем формулу х„+, = х„ + т(а с + асУ~ с + ...
+ а У ), (7) где ао, а,, ..., ар — некоторые универсальные (не зависящие от шага т) числа. Они очевиднмм образом вычисляются через интегралы от базисных ннтерполяцнонных полиномов 1'(с). При вычислении а, делается замена переменных с = с„+ т,т н рассматривается стандартный шаг по ~, равный единице.
Оценим погрешность аппроксимации, предполагая /(х), а следовательно, н решение х(с) достаточно гладкими. Погрешность экстраполяции [[у[х(с) [ — с". (с)[[ = О(тя+с) (см. з 3). При интегриро- (ч. т ОснОВы Вычисли3тльной мАтемАтики б2 ванин по (1„, ( 1 в оценке погрешности вычисления интеграла появляется еще множитель т. Переписывая (б), (7) в форме, дающей в пределе дифференциальное уравнение, получаем х((„ + ) — к((„> аобу(х„) + а, Ях„ ~ ) + ... + арф( ха р ) + 0( тР + ' ) . Итак, порядок погрешности аппроксимации равен числу используемых в (7 ) точек р + 1 .
Неудобством метода является необходимость помнить некоторое число прошлых значений г„ 3 Это, конечно, мелочь, если не считать самого начала процесса интегрирования, когда прошлого просто нет. Приходится несколько первых шагов выполнять нестандартно, например методом Рунге-Куттм (см. з 7). Метод Адамса является характерным примером схемы, формальный порядок которой превмшает порядок дифференциального уравнения.
Стандартный алгоритм начинает работать лишь при задании, кроме начальных данных хса еше и значений х„х, ..., х . Таким образом, общее решение разностного уравнения содержит больше, чем нужно, произвольных постоянных и, следовательно, какие-то лишние «решения». Полезно иметь представление о том, во что переходят лишние решения в пределе при т О. Рассмотрим простейшее уравнение х *= ах и две схемы второго порядка — примитивную схему (4) и квалифицированную схему Адамса второго порядка: Ха Ы вЂ” х„ " = — У(х„) — — У(х„,), У(х) Вя ах. В этом простом случае можно вычислить и проанализировать общие решения разностнмх уравнений.
Они ищутся в стандартной для однородных разностных уравнений с постоянными коэффициентами фОрМЕ Х„= Сита, + Сздз, ГдЕ рн аз — КОРНИ ХараКтсрнетнЧЕСКОГО уравнения, С„С2 — постоянные, определяющиеся в данном случае начальными данными хо и хи Характеристические уравнения получаем, подставляя ра в уравнение.
Для простейшей схемы (4) имеем а' — 2 д — 1 О. .. а,, а Л + ( )'. Первый корень (при (ат~ м 1); ( )2„1. р( 3) ат „1. П(ЕЗ) аа еаа.(1 .1 О(тз)) и Р(УЕ) бз 3 51 ИНТЕГРИРОВАНИЕ ЗАДАЧИ КОШИ ДЛЯ СИСТЕМ ОДУ Это решение в пределе дает. решение дифференциального уравнения. Второй корень: д ~ 1+От дь ( 1)АЕ ААЧ т.е. сеточная функция никакого разумного предела не имеет (из-за множителя ( — 1)" = ( — 1)'~'). Это есть паразитическое решение, появившееся из-за превышения порядка разностного уравнения над порядком дифференциального. В общем решении х„= С,д", + С дз, и дая того чтобы схема имела второй порядок точности, нужно обеспечить соотношения С = х„+ О(тз), С О(т~). Аналогичные выкладки для схемы Адамса приводят к характеристическому уравнению дз = д + — д сд — — ат.
з г г Его корни: з дьз з+ 4 От~ Предоставим читателю убедиться, что д — ат/2, д, = 1+ ат+ азтз/2+ О(т') = е" (1+ О(гз)) Таким образом, дп' е"(1+ О(тз)), а паразитическое решение дз ~ (ат)" Очень быстро стремится к нулю (мы, конечно, считаем ~ ат ~ ж 1, например ат м 0.1). Итак, выбор х, в схеме Адамса должен обеспечить соотношение С, = хв+ О(тз). Полезно проверить, что выбор х~ = хо + т/(хо) обеспечивает требуемые соотношения. Более высокий интеллектуальнмй уровень схемы Адамса (по сравнению с примитивной схемой с центральной разностью) сказался в том, что паразитическое решение этого метода очень быстро убывает — как (ат)п'. Численное интегрирование на ЭВМ. Представленный выше анализ погрешностей приводит к выводу, что точность численного интегрирования тем выше, чем меньше шаг т.
Это верно только до известного предела — до тех пор, пока погрешности округления, связанные с конечной разрядностью машинной арифметики, остаются пренебрежимо малыми, Реальная расчетная формула (метода Эйлера, для определенности) в действительности прн реализации на ЭВМ имеет вид х„> = х„+ т/(х„) + с„. 64 основы вычислитввьиой мАтвмАтики 1ч.
1 Величина в„обычно определяется погрешностью округления щш машинном представлении х„, т.е. имеет величину 10 '~х~ (г = 12 на БЭСМ-б, г= б-;- 7 на ЕС, г =16 на ЕС при двойной точности). Погрешность округления т/, если г вычисляется с машинной точностью, несущественна, так как обычно ~ тД м ~ х~ (х мало меняется за шаг). Но если функция г вычисляется сложным алгоритмом, ее погрешность может определять величину в„. Таким образом, машинная вычислительная формула имеет внд *."„-*."+ ос)(1.~*а~-ф. где Ь вЂ” относительная погрешность вычисления /, т~ — относительная погрешность представления х. Можно трактовать эту формулу как точную с погрешностью вычисления г.
И теперь ясно видно, что, начиная с некоторых малых величин, дальнейшее уменьшение т приводит к падению точности. Интегрирование уравнений высокого порядка. Пусть требуется проинтегрировать уравнение 14к — 4 — — У(х), к(0), ..., х(0) заданы. Не составляет труда построить разностное уравнение: хм+а — 4кп+1+ах„-4х„1+хк в -- -Их.). Данные Коши позволяют вычислить значения хо, хи хм хз, по ним находим х4, затем хз и т,д. Однако здесь конечная разрядность машинной арифметики имеет еще худшие последствия. Машинная расчетная формула имеет вид х„+, —— 4х„+, — бх„+ 4х„, — х„з+ т~У(х„) + т1! х~. Таким образом, погрешность округления т~ ~ х ~ эквивалентна погрешности вычисления У порядка т1 х!/т4. К счастью, есть простой выход — переход от уравнения четвертого порядка к системе уравнений первого порядка: х'=хз, хз кз, хз=х4, хх=г(х').
Именно по этой причине теория численного интегрирования строится для систем уравнений первого порядка. Замечание. Выше без определений были использованы некоторые понятия (аппроксимация, точность и т.п.). Смысл их достаточно прозрачен. Он уточняется в следующем параграфе. Б 61 АБСТРАКТНАЯ ФОРМА НРНБЛНЖБННОГО МЕТОЛА й 6. Абстрактная форма прнбпнигвиного метода Приближенное интегрирование задачи Коши послужит нам удобным примером, на котором можно будет ввести основные обьекты общего приближенного метода и установить связи между ними. Настоящий параграф носит, так сказать, идеологический характер, в нем появляются фундаментальные понятия теории приближенных методов вычислений. Итак, мы исходим из задачи, записанной в общей форме: А,(к') = Р'".