Амосов А.А., Дубинский В.А., Копченова Н.В. Вычислительные методы для инженеров (1994) (1095853), страница 85
Текст из файла (страница 85)
Одним из наиболее популярных среди них является л]етод Нулерова четвертого порядка точности: У вЂ” 12 У(~п-1 Уп-1) + 10 1(~п Уе) + У(~и+1 Уп+1)). у 14.11. Жесткие задачи 472 1. Понятие о жестких задачах. В последние годы при решении задачи Коши явными методами Рунге — Кутты и Адамса значительное число исследователей сталкивается с весьма неожиданным и неприятным явлением. Несмотря на медленное изменение искомых функций расчет приходится вести, казалось бы, с неоправданно мелким шагом Ь. Все попытки увеличить шаг и тем самым уменьшить время решения задачи приводят лишь к катастрофически большому росту погрешности.
Обладающие таким свойством задачи получили название жест]сих. Сразу же подчеркнем, что жесткость является свойством задачи Коши (а ~е используемых численных методов). Жесткие задачи встречаются в самых различных областях науки и техники. Традиционными источниками появления таких задач явля— ются химическая кинетика, теория ядерных реакторов, теория автоматического управления, электротехника, электроника и т.д. Жесткие задачи возникают также при аппроксимации начально-краевых задач для уравнений в частных производных с помощью полудискретных методов (методов прямых).
Трудности, возникающие при численном решении жестких задач, продемонстрируем на таком примере. 11ример 14.22. Вычислим значения приближенного решения задачи Коши д (т) = — 25у (1) + соИ+ 25ап1, у (О) = 1, (14.126) используя метод Эйлера. Решением этой задачи является функция у (т) = а1п1+ е Ж Как нетрудно видеть (рис. 14.16), на начальном переходном участке 0 ~~ 1 ( 0.2 решение быстро меняегся. Однако уже через небольшой интервал времени переходная часть решения е ~5~ практически исчезает и решение становится медленно меняющимся. Рис. Ц.1о Естественно, что приближенное решение на переходном участке приходится вычислять, используя достаточно мелкий шаг. Однако при 1 ) 0.2, когда переходная часть решения, казалось бы, уже практически отсутствует, возникает желание перейти к вычислению со сравнительно крупным шагом.
Предположим, что при 1 = 0.2 найдено значение решения у (0.2) ~ 0.205 с точностью с = 10 з. Возьмем шаг Ь = 0.1 и будем вычислять решение при 1 1~ 0.2, используя метод Эйлера: 473 уп+1 = уп + Ь (-25уа + совтп + 25в1птп). Полученные значения приближений и точные значения решения приведены в первых трех столбцах табл. 14.6. Соответствующая ломаная Эйлера изображена на рис. 14.16.
Как нетрудно видеть, метод ведет себя неустойчивым образом и оказывается непригодным для решения рассматриваемой задачи при значении шага Ь = 0.1. На конкретном примере мы убедились в том, что условие абсолютной устойчивости (14.95) нарушать нельзя. В рассматриваемом случае Л = — 25 и это условие равносильно требованию Ь 4 0.08. Взяв удовлетворяющий этому условию шаг Ь = 0.025, получим вполне приемлемое решение (см. табл.
14.6). Т а б л и ц а 14.6 Неявный Метод Эйлера; Ь = 0.025 Точное решение Неявный Метод Эйлера; Ь = 0,1 метод Эйлера; Ь = 0.3 метод Эйлера; Ь= 0.1 0.205 0.478 0.714 0.886 0.980 Попробуем теперь воспользоваться для решения задачи (14.128) неявным методом Эйлера: уп+1 = уп + Ь (-25уп+1 + сойти+1 + 25Ыд,1). Напомним, что он А-устойчив и, следовательно, абсолютно устойчив для всех значений Ь.
Используя для нахождения решения формулу 474 0.2 0.3 0.4 0.5 0,6 0.7 О,В 0.9 1,0 1.1 1,2 1.3 1.4 1.5 0.205 0.296 0.389 0.479 0.565 0.644 0,717 0.783 0,841 0.891 0,932 0.964 0,985 0.997 0.205 0.287 0.404 0.460 0.597 0.600 0.787 0.683 0.996 0.664 1.277 0.451 1.759 -0.158 0.205 0.296 0.390 0.480 0.565 0.644 0.718 0.784 0.842 0.892 0.932 0.964 0.986 0.996 0.205 0.304 0.391 0.479 0.564 0.643 0.716 0.782 0.840 0.890 0.930 0,962 0.984 0.996 ун + Ь (соэ1н+1 + 25аш1н+1) Уп+1 1+ 25Ь Э при Ь = 0.1 получим значения, совпадающие с соответствующими значениями решения с точностью 2 10 Э (см.
табл. 14.6). Если же нас устраивает точность е = 6 10 з, то шаг Ь можно утроить. Вычисленные с шагом Ь = 0.3 значения также указаны в табл. 14.6. 11'(1) = Ау(т), (14.129) Предположим, что А — матрица простой структуры и все собственные числа этой матрицы имеют отрицательные вещественные части (КеЛ„( ( КеЛ,„1~ ... ь КеЛ1 < О). Как отмечалось в предыдущем параграфе, в этом случае решение задачи асимптотически устойчиво и представляется в виде Если среди собственных чисел Л, имеются числа с сильно различающимися значениями вещественных частей, то возникает проблема, связанная с наличием в решении у компонент, имеющих существенно 1 различные временные постоянные г; = —.
Через довольно корот- ~ КеЛ,~ кий интервал времени поведение решения будет определяться наиболее слабо меняющейся (медленной) компонентой решения. Так, если 1 ВеЛа ( ... ( НеЛг < НеЛ1 < О, то при 1 — 1о > т~ — — ~ — Л-т справедливо ~ ВеЛ~ ~ приближенное равенство у (1) м с1е еь В то же время применял 1(1-ц) 476 Приведенная в примере ситуация типична для жестких задач.'При использовании классических явных методов наличие в решении быстро меняющийся хесп~кой ко.нвонснщм даже на том участке, где ее значение пренебрежимо мало, заставляет выбирать шаг Ь из условия абсолютной устойчивости. Для жестких задач это ограничение приводит к неприемлемо малому значению шага Ь. Поэтому численное решение таких задач требует применения специальных неявных методов.
Простейшим из них является неявный метод Эйлера. На рис. 14.17 и 14.18 проиллюстрировано принципиальное отличие результатов вычислений, осуществляемых с помощью явного и неявного методов Эйлера. 2. Жесткие задачи для систем дифференциальных уравнений Рассмотрим задачу Коши для системы линейных дифференциальных уравнений с постоянными коэффициентами емый для решения задачи Коши метод должен обладать свойствами устойчивости, позволяющими подавлять наиболее быстро меняющуюся Л (1-го) (жесткую) компоненту с,„е ~ . е„погрешности Е(1) = а1Е ( С1+ атЕ ( )аз+ ...
+ вас См. Отметим, что для медленной компоненты решения временной постоянной является величина ( ш(п (КеЛ~)) ', а для жесткой компоненты 1~Ха погрешности — величина ( п1ах ~КеЛь!) 1. Их отношение и определяет 1< Мн степень жесткости задачи. Приведем применительно к системе (14.129) одно из определений жесткости. Пусть КеЛ1, < 0 для всех к = 1, ..., т.
Определим число жесткости системы (14.129) с помощью формулы П1ах ~ ЯаЛ1С~ 1~ Ма (14.130) ь ТьлТ 1<Хи Систему уравнений (14.129) назовем жесткой, если для нее 8 Э 1. Пример 14.23. Рассмотрим систему у' = -7501у — 2499», »' = -7497у - 2503»; у (0) = О, » (0) = 4.
Собственные значения матрицы козффициентов А= 7501 24991 -7497 -2503! Л»1 таковы: Л1 = -4, Л2 — — -10000. Здесь число жесткости з = = 2500 много Л1 больше единицы и позтому систему можно квалифицировать как жесткую. Общее решение системы имеет внд у (с) = с1е М + сто 1оооос» (~) = с1е 41 — Зс~е 1оооос Начальным значениям у (0) = О, » (0) = 4 отвечает решение у (1) — -41 + е-100001» (1) — 3е-41 + е-100001 477 Жесткая компонента решения здесь быстро затухает и через очень небольшой интервал времени решение будет практически совпадать с у (М) = — е 4~, г (1) = = Зе ~т.
Решение этой задачи с помощью явных методов очень неэффективно. Например, для метода Рунге — Кутты четвертого порядка точности условие 2.8 устойчивости 6 < в данном случае приводит к необходимости шах ВеЛ~~ использования очень мелкого шага Ь < 2.8 10 ". 3 а и е ч а н и е. В приведенном выше наиболее простом определении жесткости требовалось, чтобы ВеЛь < О для всех Е Более общие современные определения жесткости [9], [26] не исключают наличия собственных чисел Л, с положительными вещественными частями при условии, что для них В.еЛ; < шах ~ КеЛ~~. 1< Хт В случае, когда матрица А зависит от 1 (т.е.
решается система линейных уравнений с переменными коэффициентами), число жесткости также зависит от 1 и определяется по формуле (14.130), в которой уже Л~ = Л~(М). Так как э = э (Х), то система у' = А (М)у может оказаться жесткой на одном интервале времени 1 и нежесткой — на другом.
В настоящее время нет общепринятого математически корректного определения жесткости задачи Коши для системы нелинейных уравне- ний (14.131) ю' = Х(~, ю) Чаще всего эту задачу называют жесткой в окрестности точкй] (т, у (1)), если жесткой является соответствующая линеаризованная:~' система у'(~) = А» (1), А = ~„'(1, у ®). (14.132); В подтверждение правомерности такого определения жесткости можно ~ сослаться на то, что (см. предыдущий параграф) погрешность я ®' решения нелинейной системы удовлетворяет приближенному равенству.!;.
я (1) в Аа (1) и поэтому в малой окрестности точки (т', у (1))погрешности решений линеаризованной системы (14.132) и системы (14.131)" должны вести себя примерно одинаковым образом. К счастью, для того чтобы распознать жесткую задачу на практике, часто совсем не обязательно проводить математически строгое ее исследование. Если система уравнений (14,131) правильно моделирует 479 реальное физическое явление, включающее процессы с существенно различными временными постоянными, то соответствующая задача Коши должна быть жесткой.