CC BY
15
УДК 519.67
ОСНОВНЫЕ ЭТАПЫ РАЗРАБОТКИ И ТЕСТИРОВАНИЯ ВЫСОКОТОЧНЫХ СПЕЦИАЛИЗИРОВАННЫХ ИНТЕГРАТОРОВ
Абдульмянов Т.Р., КГЭУ, канд. физ.-мат. наук, доцент, taguir2007@mail.ru Петрова Н.К., КГЭУ, канд. физм.-мат. наук, доцент, пк_petrova@mail.ru
В данной работе анализируются этапы построения ядра интегратора высокого прядка для численного интегрирования систем обыкновенных дифференциальных уравнений второго порядка. Получены разностные уравнения, аппроксимирующие исходные уравнения с точностью 0(Ь10) и исследована устойчивость решения разностных уравнений системы. Разработана компьютерная программа численного интегрирования. Тестирование работы интегратора выполнено при помощи точных решений. Ключевые слова: численные методы интегрирования систем, разностные уравнения, аппроксимация, устойчивость, сходимость к точному решению.
ВВЕДЕНИЕ
Моделирование физических процессов и явлений состоит из трёх основных этапов: экспериментальной части, этапа построения модели и этапа проверки адекватности модели. Такая классическая схема для этапов моделирования может быть существенно изменена в зависимости от исследуемого объекта. Например, исследуемый объект может быть доступен только для наблюдений. То есть, проведение прямых экспериментов может оказаться вовсе невозможным. В таких случаях экспериментальная часть моделирования будет отсутствовать. Второй этап, этап построения модели, может быть менее сложным или же, напротив, более сложным в зависимости от исследуемого объекта и методов моделирования. Этот этап моделирования и этап проверки адекватности модели с применением численных методов рассматривается в данной работе.
26
Теоретической основой применения численных методов интегрирования обыкновенных дифференциальных уравнений, уравнений в частных производных и систем дифференциальных уравнений является теорема о сходимости решения разностного уравнения к точному решению исходного уравнения [1]. Согласно этой теореме, необходимыми и достаточными условиями сходимости решения разностного уравнения к точному решению исходного уравнения являются условие аппроксимации и условие устойчивой зависимости разностного решения от начальных условий [2]. Причём, выполнение каждого из этих двух условий является одинаково важным для сходимости решения разностного уравнения к точному решению исходного уравнения. В ходе проверки аппроксимации разностного уравнения определяется погрешность аппроксимации, определяющая главным образом окончательную точность решения задачи. Однако, в теореме о сходимости необходимым условием является также и условие устойчивости решения разностного уравнения. Исследование устойчивости решения разностного уравнения почти всегда является сложной задачей. По этой причине часто этим исследованием пренебрегают, считая его не настолько необходимым, как условие аппроксимации. Однако в теоремах не бывает более или менее важных условий теоремы. Если не выполняется хотя бы одно из необходимых условий теоремы о сходимости, то нет никаких оснований считать, что получено решение исходного уравнения. Следовательно, проверка адекватности модели также оказывается формальной, а по сути, не является такой проверкой.
РАЗНОСТНЫЕ УРАВНЕНИЯ ДЛЯ ЯДРА ИНТЕГРАТОРА
Рассмотрим пример построения разностных уравнений десятого порядка для системы из четырёх неоднородных обыкновенных дифференциальных уравнений второго порядка:
27
<12 Р 2 <2 Я 2
—— = Р + V Д оояШ/ + р.), —— = -®1 Я + V К ео8(ю,/ + р.),
<2 2 V г г ^ 3 V
(1)
< 2£ 2 _ < 2 _ — = -ю22£ + Vнг 81П(^ + (), — = -Юз2Г + VЕ, 81П(^ + рг),
где Юг, CDi, р, Д, ^г, Нг, Кг - произвольные константы. Общее решение системы (1) имеет следующий вид:
2 Д 2 к
Р = Х 2 -2 С08(^ + ( ), Я = Х 2 2 С08(^ + ( ),
,=1 ^2 -аг ^з (2)
2 Н 2 Е
Б = V 2 -2 81П(^ + Р X Т = V 8Ш(^ + ( ).
Двухшаговые разностные уравнения для системы (1) с точностью 0(к10) имеют следующий вид:
иг (/+1) = 2и, (/) - и, (/ - 1) + /г (/)к2 + /г (2)(/)к4/12 +
+ /г (4)(/)кб/360 + /г (б)(/)к8/20160, г =1, 2, 3, 4, (3)
где/ (/) - значения правых частей уравнений (1) в точках ^ = к/, а /^О - значения производных функций / по времени I порядка к в точках ^ = к/. Начальные значения иг(0) и иг(1) определяются при помощи решения (2). Разностные уравнения (3) имеют явный вид. Поэтому их решения можно найти при помощи простых итераций. В случае линейных дифференциальных уравнений построение разностных уравнений, аппроксимирующих исходные уравнения с необходимой точностью, не представляет большой проблемы. Разностные уравнения (3) аппроксимируют уравнения системы (1)
28
с точностью 0(^10). Исследуем устойчивость решений разностных уравнений (3) к малым изменениям начальных условий и к малым изменениям правых частей уравнений (1). При помощи разностных уравнений (3) получим:
К (2) < 2 и (1) - иг (0) + С, и, (3) < 2и (2) - и, (1) + с = зи (1) - 2и (0) + 3с, и, (4) < 2 и (2) - и, (1) + с = 4 и (1) - 3 и, (0) + 6с,
и, (] + 1) < 2и, (]) - и, (] - 1) + с = (] + 1)и (1) - и (0) + дГс, (4)
где с - абсолютный максимум правых частей разностных уравнений (3), ф] - определяются как] - ый элемент числовой последовательности: 1, 3, 6, 10, 15, 21, 28, ... [3]. То есть, ф - является частным решением следующего рекуррентного уравнения третьего порядка:
ф+з - 2ф+2 + ф+1 - Ф = 0, (5)
с начальными условиями ф0 = 1, ф1 = 3, ф2 = 6. Характеристическое уравнение рекуррентного уравнения (5) имеет один действительный корень и два комплексно сопряжённых корня. Предположим теперь, что в начальных условиях и(1) и и(0), а также в правых частях уравнений произошли малые отклонения от их исходных значений, то есть, изменились и стали равны соответственно и,(1) + 8,- , и,(0) + 5,- и с + у,. Подставляя эти новые значения вместо прежних начальных условий и правых частей в неравенство (4) получим, что левая часть этого неравенства при таких малых изменениях не может иметь большие изменения, связанные с этими малыми отклонениями. Следовательно, решения разностных уравнений (4) будут устойчивыми.
29
СРАВНЕНИЕ РЕЗУЛЬТАТОВ ЧИСЛЕННОГО ИНТЕГРИРОВАНИЯ С ТОЧНЫМ РЕШЕНИЕМ
При разработке компьютерных программ для численного интегрирования дифференциальных уравнений и систем дифференциальных уравнений возникают проблемы, которые аналогичны проблемам поиска точных решений. В частности, невозможно построить интегратор, при помощи которого можно было бы выполнить численное интегрирование всех дифференциальных уравнений или систем дифференциальных уравнений. Также как невозможно найти универсальный метод точного решения всех дифференциальных уравнений. Однако можно построить ядро интегратора, которое можно применить в численном интегрировании целого класса уравнений или систем уравнений.
30
: ...........и—
-1М2 '-
1 1001 2001 3<Ю1 4001 5 (¡01 «001 7001 3001
Щ-12
1 1001 200 1 3001 4001 5001 6001 7001 8001
Рис. 1. Величина погрешностей АХ?;) - вертикальная ось системы координат (/ = 1, 2, 3, 4) в порядке их следования; по горизонтальной оси - номер шага интегрирования
для системы уравнений (1)
Рассмотрим интегратор, разработанный в рамках данной работы. Основу этого интегратора составляют разностные уравнения (3). Проверим работу этого интегратора при помощи применения интегратора к системе дифференциальных уравнений (1) с заранее известным точным решением (2). Сравним решение (2) системы (1), взятое в узловых точках, с численным решением этой системы. На рис. 1 представлены графики погрешностей:
Д1 (Гг) = Р(Г,) - и1(Гг), Д2(Гг) = Я(Г.) - ^(Гг), Дз(Гг) = ВД - №(0), Д4(Гг) = Т(Г,) - и4(Гг),
где Г = к г. Значения произвольных констант равны: ф 1 = 0,33; ф2 =
31
=0,5; Ю2 = 0,83; Ю3 = 0,53; щ = 0,63; щ = 0,37; Б1 = 0,22; В2 = =0,46; = 0,26; ^2 = 0,39; Н1 = 0,78; Н2 = 0,3; Е1 =0,29; Е2 = 0,44. Шаг интегрирования к = 0,1. Результаты вычисления погрешностей АД^) показывают, что точность аппроксимации O(k10) гарантирована на интервале десятков тысяч шагов интегрирования при шаге интегрирования к = 0,1. При помощи численных экспериментов уменьшением шага интегрирования можно убедиться в сходимости решения разностного решения к точному решению. В случае отсутствия точного решения такая проверка может служить подтверждением сходимости решения разностного уравнения к точному решению, несмотря на то, что не является доказательством сходимости.
Рассмотрим теперь возможность применения построенного интегратора к системе уравнений линейной теории физической либрации Луны, которые близки к системе (1):
Л2 Р 2П ^ 1 _ Л2 Я 2п & 2 , —2 Р+п ■ а, —- = -щ2 я+к ■ п ■ а, Л2 2 Л сИ2 3 Л 1 1
(6)
72 о Л2Т
= -щ25 + (1 + ¿1 )01 + п(1 -кхк2)■ Я, ~Т = -а1Т + (1 + к2а -п(1 -к1к2)■ Р,
Ж2 м
где Ql = А1 БтЩ I + ф1) + А2 sin( щ I + ф2), Q2 = Bl•COS( щ I + ф1) + £2-СО$(®2 ^ + ф2),
к1 = 0,00063190362; £2 = 0,00040398188; п = 0,229970835; А1 = = -0.08925507; А2 = 0,002480931; В1 = -0,00048902224; £2 =
= -0.0049433714; щ = 0,2308957; щ = -0.00262428; Ю2 = п; Ю3 = =(к1к2)1/2п; ф1 = 1,627905082; ф2 = 0,727650661.
Все константы выражены в радианах. Сравнивая системы (1) и (6) заметим, что точное решение системы (6) также можно получить в
32
виде решения (2). Константы Бг, Гг, Иг, Ег в этом случае будут определяться по следующим формулам:
Бх = пБх + Лг щ; Б = пБг + Лг^ щ; Г = кхпЛх - Б\ ШХ; Гг =
=к\ п Л2 - Бг^ щ;
Их = Гхп(1 - кхкг)/(щ2 - щ2) + Л: (1 + к\); Иг = = Ггп(1 - кхк2)/(щ2 - щ2) + Лг(1 + кх);
Ех = - Бхп(1 - кхк2)/(®22 - щ2) + Бх(1 + кг);
Ех = - Бгп(1 - кхкг)/(щ2 - щ2) + Бг(1 + кг).
Следовательно, к таким случаям построенный интегратор может быть применён без изменений. По общему виду уравнений (6) нетрудно также определить и другие уравнения, к которым интегратор может быть применён без изменений. Однако, в случае нелинейных уравнений в программу необходимо будет внести существенные изменения. На рис. 2 представлены графики погрешностей Л/(^) для случая системы уравнений (6). В случае системы (6) уравнения связаны между собой не только основными частотами, но более существенно: функции Р и Я из первого и второго уравнений (6) входят соответственно в четвёртое и третье уравнение. В результате, появляются соответствующие изменения в графиках погрешностей Л/(^).
-4Е-10 -1-
1 Д001 8001 12001 16001 20001 24001 2Ю01
33
Рис. 2. Величина погрешностей Ау(Ь) - вертикальная ось системы координат (/ = 1, 2, 3, 4) в порядке их следования; по горизонтальной оси - номер шага интегрирования
для системы уравнений (6).
ЗАКЛЮЧЕНИЕ
Моделирование большинства современных динамических систем не представляется возможным без применения методов численного интегрирования дифференциальных уравнений и систем дифференциальных уравнений. На этапе предварительного анализа уравнений и систем уравнений часто удаётся выделить основ-
34
ные возмущения и основные частоты, которые влияют более существенно по сравнению с остальными. В таких случаях общую задачу интегрирования можно разделить на два этапа. На первом этапе - решить задачу для основных частот. На втором этапе - для возмущений высшего порядка. Такое разделение может исключить появление неустранимых погрешностей: погрешности, допущенные на первом этапе, могут оказаться неустранимыми на втором этапе интегрирования. То есть, интегратор для основных частот (ядро полного интегратора) должен быть высокоточным. Только в этом случае можно приступать к учёту возмущений высшего порядка. Разработка такого ядра интегратора для промежуточного численного решения систем дифференциальных уравнений была основной целью данной работы. В результате выполнения работы:
1) разработана компьютерная программа, применимая для численного интегрирования достаточно широкого класса систем дифференциальных уравнений второго порядка;
2) разностные уравнения, на основе которых разработан интегратор, аппроксимируют основные уравнения с точностью 0(^10);
3) решение разностных уравнение устойчиво к малым изменениям начальных условий и правых частей уравнения.
Увеличение числа шагов интегрирования всегда приводит к увеличению погрешности. То есть, какова бы не была точность аппроксимации, всегда будут накапливаться погрешности вычислений и, в конечном счёте, погрешность станет недопустимой. Для компьютерных вычислений была использована двойная точность вычислений встроенных функций. Начальные значения Цг(0) и Цг(1) также были определены с двойной точностью. В противном случае необходимая точность 0(^10) интегрирования не может быть получена.
35
Источники
1. Самарский А.А. Теория разностных схем. М.: Наука, 1983. 616 с.
2. Самарский А.А., Гулин А.В. Устойчивость разностных схем. М.: Наука, 1973. 305 с.
3. Абдульмянов Т.Р. Алгоритмы и методы решения задач дискретной математики с применением компьютерных вычислений. Казань: изд. КГЭУ, 2011. 156 с.
References
1. Samarskii A.A. Teoriya raznostnykh skhem. M., Nauka, 616 (1983).
2. Samarskii A.A., Gulin A.V. Ustoichivost' raznostnykh skhem. M., Nauka, 305 (1973).
3. Abdul'myanov T.R. Algoritmy i metody resheniya zadach diskretnoi matematiki s prime-neniem komp'yuternykh vychislenii. Kazan': izd. KGEU, 156 (2011).
Information
Abdulmyanov T.R., Petrova N.K.
THE MAIN STEPS OF THE DEVELOPMENT AND TESTING OF HIGH-PRECISION SPECIALIZED INTEGRATORS
This paper analyzes the stages of building the core integrator of high strand for the numerical integration of systems of ordinary second order differential equations. The difference equations approximating the original equations with accuracy O(h10) are obtained and investigated the stability of solutions of difference equations system. The computer program for numerical integration is considered. The testing of program is realized using exact solutions. Keywords: numerical methods for the integration of systems, the difference equations, the approximation, the stability, convergence to the exact solution.
Дата поступления 22.01.2015.
36
