CC BY
31
УДК 519.622.2
Экспериментальное исследование явных методов решения обыкновенных дифференциальных уравнений с расширенной областью устойчивости
Чадов С.Н., асп.
Излагается экспериментальное исследование явных методов численного решения обыкновенных дифференциальных уравнений применительно к жестким задачам. Приводится краткое описание явных методов решения обыкновенных дифференциальных уравнений с расширенной областью устойчивости и сравнение точности методов и вычислительных затрат.
Ключевые слова: численные методы, жесткие системы ОДУ, методы Рунге-Кутта-Чебышёва.
Experimental study of explicit ODE solvers with extended stability domain
Chadov S.N., PhD student
Рresents an experimental study of explicit ODE solvers applied to stiff problems. A short description of the explicit methods with extended stability regions is given; these methods as well as some other well-known explicit methods are applied to the sample stiff problems with comparison of accuracy and computational costs.
Keywords: numerical methods, stiff ODE systems, Runge-Kutta-Chebyshev methods.
Введение. При решении жестких систем обыкновенных дифференциальных уравнений возникает проблема устойчивости выбранного численного метода. По этой причине для таких задач приходится выбирать либо неоправданно малый шаг интегрирования, либо использовать специальные устойчивые методы, которые, как правило, являются неявными. Неявные методы, в свою очередь, требуют больших вычислительных затрат на решение на каждом шаге нелинейной системы алгебраических уравнений. Однако для достаточно большого класса жестких задач все же возможно использование явных методов, специально сконструированных таким образом, чтобы иметь более широкий интервал устойчивости, чем классические методы. Расширение интервала устойчивости, как правило, достигается за счет уменьшения точности вычислений, однако во многих случаях удается найти компромисс между точностью и устойчивостью. Ниже проводится сравнение нескольких методов с расширенной областью устойчивости с классическими явными методами на модельных жестких задачах.
Методы. Существует несколько подходов к созданию стабилизированных явных методов. Один их них [1-3] основан на явном конструировании полинома устойчивости таким образом, чтобы он имел заданные свойства. Обычно эти полиномы задают область устойчивости, вытянутую вдоль вещественной оси, и строятся на основе полиномов Чебышё-ва. Поэтому методы типа Рунге-Кутты с расширенной таким образом областью устойчивости получили название методов Рунге-Кутты-Чебышёва. Принципы их построения можно найти в работе [1]. Мы будем рассматривать методы порядков с 3-го по 8-й и в дальнейшем
метод Рунге-Кутты-Чебышёва порядка э будем обозначать РКОэ.
Другой подход, рассмотренный в работах [4, 5], предполагает оценку наибольших собственных значений матрицы Якоби и стабилизацию расчетной схемы в точках, в которых эти значения достаточно велики. В [5] такие методы названы адаптивными. Из них будем использовать методы 1-го и 2-го порядков, рассмотренные в [5], а также стабилизированную версию классического метода Рунге-Кутты 4-го порядка, также рассмотренную в [5]. Обозначим эти методы соответственно ЭкУ0|180У1, Бку01180у2 и РК4-э1аЬ.
Из классических методов рассмотрим следующие [6]: Эйлера 2-го порядка (Еи1ег2); Хойна 2-го порядка (Иеип2); Рунге-Кутты 4-го порядка (РК4); а также методы Богацки-Шампайна 3-го порядка (ББ3) [7], Дорманда-Принса 5-го порядка (РР5) [8].
Система уравнений. Для тестирования вышеназванных методов будем использовать несколько хорошо известных жестких задач.
Задача 1 (ОРЕСО) представляет собой модель, описывающую реакцию Белоусова-Жаботинского:
^ = 77,27 (у2 + У1 (1 - 8,375 • 10~5У1 - у2))
1Т = 7727 (-у + у1>у 2)
±у3 = 0,161 (у, - У3)
У1(0) = 1; У2(0) = 2; У3(0) = 3; ( = 0..30.
Задача 2 (УйРОЬ) представляет собой уравнение Ван-дер-Поля:
а ~ У2
= E ((1 - У12 )у 2
Уі
¿У 2 &
Е = 8500; t = 0..1.
Результаты. Приведем результаты решения тестовых задач различными методами. В качестве показателя точности использовалось следующее выражение:
(
3 = RMS
Z,
Уі - У І
\
У;
где у и yi - полученное и точное решения соответственно; функция RMS(X) - rooted mean square:
RMS( X) =
ШІ
N
В качестве меры вычислительных затрат было выбрано количество вычислений правой части уравнения (Р). Данная характеристика, на наш взгляд, более показательна, чем время выполнения, поскольку не зависит от производительности вычислительной системы, а также является показательной не только для относительно простого модельного уравнения, но и для достаточно широкого класса задач, имеющих спектр, качественно подобный спектру этого модельного уравнения. Эксперимент будем проводить следующим образом: выберем некоторый порог точности &. Будем считать, что если &<&, то метод обеспечивает достаточную точность. Найдем для каждого метода максимальное значение шага по времени, при котором метод обеспечивает достаточную точность.
Результаты решения задачи 1 (ОРБСО) приведены табл. 1.
Таблица 1. Результаты решения задачи OREGO
Метод Точность, 3 Вычислительные затраты,¥
RKC3 2,410-5 3,0105
RKC4 4,1 і 0-5 2,5105
RKC5 9,610-5 1,8105
RKC6 1,810-4 1,5105
RKC7 3,6-10-4 1,2105
RKC8 6,0-10-4 1,0105
Skvortsovl 8,1 10-3 5,3103
Skvortsov2 9,210-3 1,5104
RK4-stab 5,510-3 2,2104
Euler2 2,1 10-3 3,0105
Heun2 4,210-6 9,0-105
BS3 5,510-8 4,5-105
RK4 5,1 10-10 6,0105
DP5 4,510-12 9,0-105
Анализ полученных результатов (табл. 1) показывает, что обе группы специализированных методов обеспечивают значительно большую производительность, чем классические явные методы. Так, количество вычислений правой части уравнения для метода РКО8 почти в 6, а для метода Эку01180у1 - в 120 раз
меньше, чем для метода РК4, при том, что оба этих метода обеспечивают вполне приемлемую для многих приложений точность. Стоит также отметить, что для классических методов высокого порядка полученное решение имеет значительно избыточную точность (4,5-10-12 для метода РР5, что на 10 порядков выше заданного нами значения 9), и соответственно, вычислительные затраты на решение данными методами неоправданно высоки.
Рассмотрим задачу 2. Численные эксперименты показали, что для этой задачи методы Скворцова в виде, приведенном в работе [5], не являются оптимальными. Так, изменив некоторые константы, нам удалось добиться более качественного решения данной задачи. В частности, параметр с, представляющий собой коэффициент, вычисляемый в зависимости от оценки наибольшего по модулю собственного значения, т. е. в конечном итоге управляющий жесткостью вычислительной схемы, будем вычислять по-другому. В методе первого порядка с для нежестких компонент в нашем варианте алгоритма будет вычисляться
1 ь
по формуле с = — + —. В методе 2-го порядка
1
b
24
и для же-
для нежестких компонент с =
стких с = 1,13а^-+1. Результаты решения зада-
а -1
чи 2 приведены в табл. 2 с учетом этих изменений.
Таблица 2. Результаты решения задачи VDPOL
Метод Точность, 3 Вычислительные затраты,¥
RKC3 0,03 8,0105
RKC4 0,05 5,.5105
RKC5 0,03 5,0105
RKC6 0,05 5,8105
RKC7 0,05 6,7105
RKC8 0,06 5,3105
Skvortsov1 0,06 3,0105
Skvortsov2 0,05 4,0-105
RK4-stab 0,02 3,3105
Euler2 0,05 6,7105
Heun2 0,04 1,0106
BS3 0,02 4,3105
RK4 0,03 6,7105
DP5 0,02 6,7105
Анализ полученных результатов (табл. 2) показывает, что данная задача гораздо труднее поддается решению стабилизированными методами; их преимущество перед классическими методами довольно незначительно (менее чем в 1,5 раза). Это можно объяснить тем, что данная задача требовательна не только к устойчивости, но также и к точности методов. Нужно отметить, что полученные результаты противоречат результатам измерений, приведенных в [5]. Это объясняется тем, что в [5] измерение точности проводилось только в ко-
нечной точке интервала интегрирования, тогда как мы измеряем точность на всем интервале.
Заключение
Проведенные эксперименты показали, что стабилизированные явные методы могут с успехом применяться хотя бы для некоторых жестких задач. Обе группы методов (Рунге-Кутты-Чебышёва и адаптивные) показывают достаточно неплохие результаты по сравнению с классическими методами, однако на обеих задачах эффективность адаптивных методов оказалась значительно выше.
Список литературы
1. Van der Houwen P.J., Sommeijer B.P. On the
Internal Stability of Exlplicit m-stage Runge-Kutta methods for large m-values. // Z. Angew. Math. Mech 60. - 1980.
2. Medovikov A.A. Third Order Explicit Method for the Stiff Ordinary Differential Equations // Numerical analysis and Its Applications , First International Workshop, WNAA'96 Rousse, Bulgaria, June 24-26, 1996.
3. Abdulle A. Fourth Order Chebyshev Methods with Recurrence Relation // SIAM J. Sci. Comput. - Vol. 23. - № 6. - Р. 2041-2054.
4. Скворцов Л.М. Простые явные методы численного решения жестких обыкновенных дифференциальных уравнений // Вычислительные методы и программирование. -2008. - № 9. - С. 154-162.
5. Скворцов Л.М. Явные многошаговые методы численного решения жестких обыкновенных дифференциальных уравнений // Вычислительные методы и программирование. - 2008. - № 9. - С. 409-418.
6. Амосов А.А, Дубянский Ю.А., Копченов Н.В. Вычислительные методы для инженеров: Учеб. пособие. -М.: Высш. шк., 1994.
7. Bogacki P., Shampine L.F. A 3(2) pair of Runge-Kutta formulas // Applied Mathematics Letters. - (1989). - 2 (4). -Р. 321-325.
8. Dormand J.R., Prince P.J. A family of embedded Runge-Kutta formulae // Journal of Computational and Applied Mathematics. - 1980. - 6 (1). - Р. 19-26.
Чадов Сергей Николаевич,
Ивановский государственный энергетический университет, аспирант кафедры высокопроизводительных вычислительных систем, e-mail: [email protected]
