Н.Н. Калиткин - Численные методы (1133437), страница 54
Текст из файла (страница 54)
Разлагая репгение и(х) по формуле Тейлора на интервале сетки х„.:х~х„ег и обозначая и(х„) =и„, получим 1 и„, =ил+А„и„'+-й lг,',и„'+ ..., Ь„=х„ег — х„. (13) Стоящие в правой части производные можно найти, дифференцируя уравнение (7) требуемое число раз: и'=7'(х, и), и"=„— 7'(х, и) =7„+Цл (14) и т. д.
В принципе, если 7(х, и) имеет гг-е непрерывные производные по совокупности аргументов, то в разложении (13) можно удержать члены вплоть до О (иезг). Однако использовать для расчетов формулу (13) с большим числом членов невыгодно. Во-первых, даже при сравнительно простой правой части выражения для производных могут оказаться громоздкими. Во-вторых, если правая часть известна лишь приближенно, то находить ее производные нежелательно. В простейшем случае, подставляя (14) в (13) и ограничиваясь только первым членом разложения, получим схехгу ломаных"): у„„=у„+й„~(х„, у„), гг„=х„„— х„.
(15) 11оскольку при такой замене можно найти только приближенные значения искомой функции в узлах, то будем обозначать эти значения через у„в отличие от точных значений и„=-и(х„). Для численного расчета по схеме ломаных достаточно задать начальное значение гго=тг. Затем по формуле (15) последовательно вычисляем величины у„у„..., уи. Геометрическая интерпретация этой схемы дана на рис. 41, где изображено поле интегральных кривых.
Использование только первого члена формулы Тейлора означает движение не по интегральной кривой, а по касательной к ней. На каждом шаге мы ') Она была предложена Эйлером н называется также схемой Эйлера. 244 ОБЫКНОВЕННЫЕ ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ !ГЛ. 71П заново находим касательную; следовательно, траектория движежения будет ломаной линией. Исследуем сходимость метода ломаных, предполагая правую часть 1(х, и) непрерывной и ограниченной вместе со своими первыми производными: !1 ( =- М„~ 1„! ( М„~ 1л , '~ М, (отсюда следуе~, что )ил~~М,=М,+ + М,МО). Рассмотрим погрешность приближенного решения гл =у„— ил. Вычитая (15) из (13), получим соотношение, связывающее погрешности в соседних узлах сетки: Ук г„„=г„+Ь„[((х„, ул)— ! ! — ! (хл, и„)) — -- ллси„" = =г.
(1+А!').— -2 йси." (16) г'г Рис. 4! (члены более высокого порядка малости здесь спущены). Последовательно применяя рекуррентное соотношение (!6), выразим погрешность на произвольном шаге через погрешность начальных данных т — ! кл — ! т — 2 =ЕО П (1+Ю -2 — ~! )2„'и" и (1 1 А!л) (1у) л=О А=-л Ф ! Отсюда нетрудно дать асимптотическую оценку погрешности. Заметим, что при малых шагах сетки т — ! П (1+)21О)„= П ехр (й|,)„= л О л О 1-т — ! Гт-1 =ехр~.У', (!!(„)„~=ехр~ ~ 1„(т, и(т))!(т л=О кк причем в качестве верхнего предела интеграла можно взять х, ибо ошибка при атом остается в пределах общей точности преобразований.
Аналогично преобразуя второй член (17), получим !к ! к !кт *. = ., -! () !. к ) - -,' 5 к к ! ! ! ! -~()' !. кк~ а ! кк К к Здесь й(х) — непрерывная функция, дающая в каждом узле хл величину шага л„; в качестве такой функции можно выбрать линейный сплайн. 245 ЗАДАЧА КОШИ ~е„',(М(х ) шах й„=О(шахй„), О<а<а (19а) где "а рл ррр р--; ( р "р р»р(1 'р.,рр)~р„'Рр' 'р. ррррр Таким образом, при й — 0 ггрибргиокенное решение сходится к точному равномерно (на ограниченном отрезке! х — ха! ==а) с первым порядком точиостгг. 3 а меч а н не !.
Оценка погрешности (19) является мажорантной. Для функций со знакопеременными производными эта оценка может быть сильно завышена по сравнению с асимптотической оценкой (18). Замечание 2. Экспоненцнальный член в оценке (18) характеризует расхождение интегральных кривых (см. рис.
41); если он очень велик, то исходная задача 1(оши плохо обусловлена. П р и м е р. Проинтегрируем по схеме Эйлера задачу Коши для уравнения (3): и'(х)=ха+и', О.==х=-1, и(0)=0. В таблипе 18 даны численные решения у(х), полученные на сетках с шагами й=-1, '!а и гр4, столбик у(х) будет пояснен в п. 10. Приведено также точное решение и (х), вычисленное методом Пикара (см. пример в п.
3). Видно„что схема Эйлера для получения удовлетворительной точности требует гораздо более малого шага, чем использонанный здесь. Таблица 18 Рассмотрим структуру погрешности (18). Первое слагаемое справа связано с погрешностью начального значения е,=у,— и„ которая умножается на ограниченную (благодаря ограниченности производных) величину. Начальное значение можно задать точно и считать, что еа — -О. Остановимся на втором слагаемом. Оно обусловлено тем членом формулы Тейлора (13), который был отброшен при выводе схемы ломаных (1б). Оценим это слагаемое сверху; заменяя все функции под интегралами их модулями и вынося шахй(х) за знак интеграла, получим ЗАДАЧА КОШИ Ограничимся только написанными членами, так как уже они обеспечивают четвертый порядок точности. Для вычисления решения в следующей точке запишем дифференциальное уравнение в интегральной форме "л., и илп,=ил+ ~ 7(х, и(х))ггх=и„+ ~ Г(х)г(х (29) кп кл и подставим в него интерполяционный многочлен (28).
Получим формулу Адамса для переменного шага ! у =у.+й.Г(х.)+--й",Г(х., х.,)+ + — й'„(2йл+Зй„,) Г(хл, хл. „хл,) + + 1!2 ййп (Зйл+ 8йлйп -!+ 4йлйл-г+ 6Еп г+ 6йп — гйк-г) х хГ(х„, хл „хл „хл,), где йл=хл„— хл. (30) Эта формула имеет четвертый порядок точности. Если отбросить последнее слагаемое, получим формулу третьего порядка точности. Аналогично получаются формулы низших порядков. Формула первого порядка совпадает со схемой ломаных. Чаще пользуются менее громоздким вариантом формулы (30), рассчитанным на постояиный шаг интегрирования. Вместо разделенных разностей вводят конечные разности ААГл пл = р(Г (хл, хл „..,, хл р), приблизительно равные р-й производной в точке (хл+х„р)/2, и получают у г=у +йГл+ и й'Ь'Г„+ТййлЛ'Г„+-ц-)гА'Г„.
(31) Остаточный член этой формулы равен (251/750) йлГгч(х). Метод без изменений переносится на системы уравнений первого порядка типа (25). Чтобы начать расчет методом Адамса, недостаточно знать у(х,). Для начала расчета по формуле (30) надо знать величину решения в четырех точках х„, х„ х„ х, (а при формуле р-го порядка точности — в р точках). Поэтому надо вычислить недостающие значения у„каким-либо другим методом — методом Рунге — Кутта, или разложением по формуле Тейлора (13) — (14) с достаточно большим числом членов. При работе на ЭВМ это вдвое увеличивает обьем программы.
Кроме того, формулы (30) громоздки, а несложные формулы (31) рассчитаны только па постоянный шаг и требуют нестандартных действий при смене шага: надо перейти к формулам (30), сделать по ним четыре шага и снова вернуться ЗАДАЧА КОШИ 4 ]] (21), равномерно сходится к точному решению с погрешностью 0(шахй»), т. е. двучленная схема Рунге — Кутта имеет второй порядок точноспш. Формула (21) имеет неплохую точность и нередко используется в численных расчетах.
При этом обычно полагают либо упн апн хп п«Г» хт» Рис. 42. Рис. 43. а=1, либо а='/». В первом случае получается схема особенно простого вида 1 ! у„„=у„+]1~(х„+ — и, у„+ — й(„). (22) Ее смысл поясняется рис. 42. Сначала делаем половинный шаг 1 по схеме ломаных, находя у„+ н»=у„+--п(„. Затем в найденл' ной точке определяем наклон интегральной кривой у,', + и» = =] (х„+]ль у«+ ы»). По этому наклону определяем приращение функции на целом шаге у„и„=у„+йу„'+ и», Геометрическая интерпретация второго случая У„„,=У,+ — Д(х„, У„)+1(х„+й, У„+й(„)) 6 изображена на рис. 43.
Здесь мы сначала грубо вычисляем по схеме ломаных значение функции у„.,=у„+Й1„и наклбн интегральной кривой у„'+]=~(х„»„у„»») в новой точке. Затем находим сРеДний наклон на шаге У„' си»=(У«'+Уи'4.1)2 и по немУ Уточняем значение у„,, Схемы подобного типа нередко называют «предиктор — корректор». Табл и и а 19 пи и (х) ь=с,а 0,000 0,250 0,00 0,50 1,00 0,000 0,031 0,317 0,000 0,042 0,339 0,000 0,042 0,350 248 ОБЫКНОВЕННЫЕ ДИФФЕРЕНЦИАЛЬНЫЕ УРАБНЕНИЯ !ГЛ УИ! В таблице 19 приведен численный расчет по схеме (22) того же примера, который рассмотрен в таблице 18 (п.5). Из таблицы видно, что схема второго порядка точности дает существенно лучшие результаты, чем схема ломаных; уже расчет на грубой сетке с й = 0,5 можно считать удовлетворительным.