Семейство шестиэтапных методов шестого порядка Текст научной статьи по специальности «Математика»

УДК 519.62/.642 Вестник СПбГУ. Прикладная математика. Информатика... 2018. Т. 14. Вып. 3 МБС 65Ь06

Семейство шестиэтапных методов шестого порядка

И. В. Олемской, Н. А. Коврижных

Санкт-Петербургский государственный университет, Российская Федерация, 199034, Санкт-Петербург, Университетская наб., 7—9

Для цитирования: Олемской И. В., Коврижных Н. А. Семейство шестиэтапных методов шестого порядка // Вестник Санкт-Петербургского университета. Прикладная математика. Информатика. Процессы управления. 2018. Т. 14. Вып. 3. С. 215-229. https://doi.org/ 10.21638/11702^рЬи10.2018.303

Работа посвящена построению экономичного явного метода шестого порядка численного интегрирования систем структурно разделенных обыкновенных дифференциальных уравнений. Приведена общая схема метода, алгоритмически учитывающая выделенные структурные особенности описываемой канонической формы систем структурно разделенных уравнений. Выписаны условия порядка и упрощающие условия для рассматриваемого явного одношагового шестиэтапного метода и найдены условия их совместности. Дан алгоритм построения семипараметрического семейства явных шестиэтапных методов шестого порядка, и построена расчетная схема метода для фиксированных значений параметров. На решении двух тестовых задач проведено ее сравнительное тестирование с тремя явными одношаговыми методами шестого порядка. Ключевые слова: порядок, условия порядка, упрощающие условия.

В работах [1,2] выделен класс систем структурно разделенных обыкновенных дифференциальных уравнений. В [3] дан алгоритм приведения произвольной системы обыкновенных дифференциальных уравнений первого порядка к канонической форме систем структурно разделенных обыкновенных дифференциальных уравнений

у'г = Мх,У11 •••,Уг—1 ,У1+1,...,Уп), г = 1, 2,...,1, (1)

у'з = (х,У\,...,Уз-\), 3 = 1 + 1,...,п1 п> 2, (2)

П

р = ти ,х е \Х0,Хк] С М, У8 : [Хо ,Хк] , в = 1,...,п,

I

и :[Хо,Хк] X М-"* —^ М^, Р = ти, г = 1, 2,...,1,

и=г п

и :[Хо,Хк] X М— Мг;, Г3 = ^ т„, з = I + 1,...,п.

V=3

Группы уравнений (1), (2) структурно тождественны. Каждое уравнение одной из них занимает определенное место в последовательности уравнений своей группы. Его правая часть не зависит от искомых функций, поведение которых описывается этим и всеми последующими уравнениями той же группы. Для численного интегрирования систем (1), (2) в [1, 2] предложен явный одношаговый метод, в рамках которого построены расчетные схемы д-го порядка с наименьшим числом этапов т(д) = д — 1, 9 <5.

© Санкт-Петербургский государственный университет, 2018

В [4] для систем (1), (2) (при I = 1, п = 2) с перекрестной структурой связей были построены экономичные шестиэтапные методы шестого порядка. Здесь в развитие структурного подхода для систем (1), (2) со снятыми ограничениями на размерность групп уравнений (п > 2) строится шестиэтапный метод шестого порядка. Необходимость в интегрировании систем такого типа возникает, например, в задачах небесной механики, физики высоких энергий.

Метод интегрирования. Считаем, что известно точное решение ув(х), в = 1,...,п, системы (1), (2) в точке х € [Хо, Хь]. Не умаляя общности рассуждений, для простоты вывода будем полагать т8 = 1, в = 1,...,п.

Для численного интегрирования систем (1), (2) рассматривается явный одноша-говый метод. В предположении достаточной гладкости правой части этой системы приближение ха к точному решению уа(х + к), в = 1,...,п, в точке х + к € [Хо, Хь] ищем в виде

Уг (х + к) х ^ = уг(х) + к Ь^к^и, (к), г = 1,...,1,

ш = 1

(3)

Уз (х + к) х г^ = Уз (х) + к ^ (к), ] = I +1,...,',

ш = 1

причем кв,и = (к) вычисляются в строгой последовательности

к1,1,..., кп,1, к1,2,..., кп,2, к1,з, к2,з,... (4)

по формулам

кг,и = /г { х + С1,шк,у1 + нУ^ (н),У1+1 +

V и=1 и=1

Ш — 1 Ш — 1 \

+ к а1,2,и,и кь+1,и (к), ...,УП + Н^2 а1,2,и,и ки,и (км , г = 1, 2,...,1, (5)

и=1 и=1 )

(Ш Ш

х + С2,иН, У1 + Н^2

(н),У1+1 +

и=1 и=1

Ш Ш \

+ Н^2а2,2,и,-и кь+1,„ (к),..., Уз-1 + Н^2а2,2,и,-и к^-1,и (к) I , ] = I +1,...,п.

и=1 и=1 )

Здесь Ъчи, сд,и, ад,р,и,„, д € {1, 2}, р € {1, 2}, — параметры метода (3)-(5), д — номер блока, р — номер группы, т — число этапов (стадий), к — шаг метода. Для явного метода (3)-(5) с1,и = а1,1,1,1 = а1,2,1,1 = 0. Количество параметров метода равно 2т2 + 5т — 3.

Выбор параметров метода (3)-(5) д-го порядка осуществляется таким образом, чтобы методическая погрешность

\у3(х + к) — х3\ = 0(кд+1), в = 1,...,п.

Для компактного представления рассматриваемого метода используем табличную форму блочного (покомпонентного) представления шестиэтапного метода (3)-(5) интегрирования системы (1), (2) (табл. 1).

т

т

■Ш

ш

Таблица 1. Шестиэтапный метод (3)—(5) интегрирования системы (1), (2)

С1,то «1,1 «1,2 ,ги,д ¿>1,-ш

0 С1,2 С1,з С1,6 0 «1,1,2,1 «1,1,2,2 «1,1,3,1 «1,1,3,2 «1,1,3,3 «1,1,6,1 «1,1,6,2 «1,1,6,3 ...... «1,1,6,6 0 «1,2,2,1 «1,2,3,1 «1,2,3,2 «1,2,6,1 «1,2,6,2 ...... «1,2,6,6 Ь1,1 ¿>1,6

С2,-Ш «2,1,-ш,^ «2,2,-ш,^ ¿>2,-ш

С2Д С2,2 С2,3 С2,6 «2,1,1,1 «2,1,2,1 «2,1,2,2 «2,1,3,1 «2,1,3,2 «2,1,3,3 «2,1,6,1 «2,1,6,2 «2,1,6,3 ...... «2,1,6,6 «2,2,1,1 «2,2,2,1 «2,2,2,2 «2,2,3,1 «2,2,3,2 «2,2,3,3 «2,2,6,1 «2,2,6,2 ...... «2,2,6,6 Ь2,1 ¿>2,2 ¿>2,3 ¿>2,6

Условия порядка для шестиэтапного метода (3)-(5) шестого порядка с использованием предположений

ш ш —1 ш

^а1,1,™," = 53 а1,2,ш,и = С1,ш, = С2,ш, Р €{1, 2}, V = 1,...,6, (6)

V=1 V=1 V=1

представляют собой систему 314 нелинейных алгебраических уравнений:

= —ГТ' « = 0,1,...,5, Я,р,г,е,ге{ 1,2};

V ^ + 1

Е ^д.^д.^Е ая,р,^,цср,ц = 2 . (3 + /)' ^ = 1;

1

Е Е ад,р,1У,цср^ — + /) ' ^ ~~

Е Е Е ар,г,ц,£Сг,(; — + /) ' ^ ~~

1

Е ^ Е ад,р,1ЛМсг>,М^ ' ^ Е ад,г>1ЛМсг>,М^ — 4(5

Е V Е Е = ю(/+ „)(! + „)>

м «

(1,п) €{(0, 2), (0, 3), (1,1), (1, 2), (2,1)};

Е V Е Е Е = щ1 + 2и + 2ру

(¡,п,р) €{(0, 0,1), (1,0,1), (0,1,1)};

1

м

Е Е ад,р,1ЛМср>М = ТТТГ;

V м 30

(7)

( Ху ар,г,л,е сг,е ) ' ( Ху а4,р,и,^ср,^

1

72'

24(1 + I)'

^^ cq,v^2 ар,т,л,е^2 аг,е,е,ф се,ф

1

144'

л е ф

V л е Ф 360

^д,!' Ху Ху ар,г,ц,(; Ху аг,е,£,'ф Ху ач,Ь,'ф,фсЬ,ф „„

V л е ф Ф 72

Количество неизвестных параметров метода (при т = 6), подлежащих определению, равно 100. Для поиска решения системы (6), (7) были использованы упрощающие ограничения

с2,1 =0, Ъ1,2 = 0, Ъ2,2 =0; (8)

Еш 1 2 ой

1/=2 а1,1,и),иС1,1< = 2С1,«" ^ = 2,..., 6,

Т,1=-Ш = ^-уЧ«. (! - с1^)> ^ = 2,..., 6, / = 0,1;

1

Ч^^ЯМ,^ = -&1,2 (1 - С1,23); 1 2

ЕШ-1 1 2 ой

„=1 «1,2,№,1/62,)/ — - С.\ и) , IV — . . . , О,

Ъl,v а1,2^,и = Ъ2,и (1 — с2,и), ад =1,...,6; (9)

YZ=■2 «2,1 — ад — 2,..., 6;

Т<1=ги 4^4^0.2,1,^«, = ь1и, (1 - с'-^), т = 1,... ,6, / = 0,1; Е!!=2 Ь2,1,^,2 = ^1,2 (1 - С12);

ЕШ 1 2 ой

„=1 а2,2,и,,»С2,„ = -С т = 2,. .. ,0,

4^4,,/«2,2,,/,«,= (1 - ад = 1,---,6, / = 0,1.

Методика построения вычислительных схем (с использованием упрощающих ограничений (8) и (9) в рамках структурного подхода) достаточно подробно рассматривалась в [1, 2]. Поэтому здесь, опуская как представление условий порядка, так и их исследование, приведем лишь окончательный результат — систему-следствие:

1

=ад = о, 1,..., 5; (10)

ХГ=1 ам,™,^,, = ^-ус^, / = 0,1, го = 1,..., 6; (И)

Т,1=ги = (1 - с^1), / = 0,1, IV = 2, ...,6;

Е®=2 „«1,1,^,2 = ^1,2 (1 -

Е6 у 2 2 _ ^

^=3 — уз'

Е^/а1,2,™^ = С1,ш, V =1,...,6; 1

Е™-1 2 Ой

„=1 а112,ги,1уС2,1, — 2С!.«" ^ о, . . . , О,

EV=w а1,2^,ш = Ь2,™ (1 - С2,ш), V = 1,..., 6;

Е6 , V ^ V — 1 1

„=3 0.1,2,^2^ - -;

Е6 г \-^v—1 2 1

„=3 2^д=2 а1,2,1/,мс2,и - Т7'

Ео т 2 — 1 2 _

1У=3 2^и= 2 а1,2,1/,мс2,и — 77^

Е6 7_ 1 3 1

м 24' ад = 0,1,..., 5;

EV=1 а2,1,™^ = С2,ш, V = 1,..., 6;

Е!У=2 «гд.-ш^С!^ = ад = 2,..., 6;

4^4,¡^2,1,^,«, = -р-^Ъ^ (1 - с^1), / = 0,1,

Е!!=2 Ь2,1,^,2 = ^1,2 (1 -

Е6 Г V—1 2

v=3 с1,^м=2 а1,2,^м 2^и=3 а2,1,м,«с1,-

Е6 Г V—1 2

v=4 С1,^ 2^м=3 а1,2,^м 2^и=3 а2,1,м,ис2

л 6 *"> *-о 1

ад = 1,..., 6;

1

72' 1

72'

Е6 ь, 2 V 2

ЕШ £

v=1 а2,2,™^ с2

/ = 0,1, ад = 2,..., 6;

Е„=ш = ^-руЬг,«, (1 - 4^)' ¿ = ад = 2,..., 6;

Е6 7_ 2 V1 V 2 1

¡/=3 Ь2,ис2,^ 2^м=2 "2,2-

18'

Е6 7 \ ^V 1 2

v=3 с1,^2^м=2 а'1,2^,м1^и=2 а2,2,м,ис2,и =

72'

(12)

(13)

(14)

(15)

(16)

(17)

(18)

(19)

(20) (21)

(22)

(23)

(24)

(25)

(26)

(27)

(28)

(29)

(30)

(31)

(32)

(33)

Система-следствие (10)-(33) состоит из двух взаимосвязанных (параметрами , с2^,Ь1^, Ь2^) блоков <1> и <2>, где <1> — уравнения (10)-(21), <2> — (22)-(33). Каждый блок, в свою очередь, содержит по две группы структурно разделенных по параметрам а1}1^д, а1,2^,д, а2,1^,д, а2}2^}д ограничений. Две группы: <1,1> (уравнения (11)—(14)) и <1,2> (равенства (15)-(21)) — входят в первый блок ограничений.

1

Второй блок содержит группы <2,1> равенств (23)-(29) и группы <2,2> ограничений (30)-(33).

Задача о поиске решения исходной системы условий порядка сведена к поиску частного решения системы-следствия (10)-(33) с более простой структурой связей, состоящей из 105 уравнений с 95 неизвестными. Последовательно в соответствии с ограничениями для всех четырех групп уравнений находим решение системы-следствия. Компактное изложение поиска решения представим в форме алгоритма его построения.

Алгоритм построения явного метода шестого порядка. Используя основные предположения (8) и три пары упрощающих ограничений: (12), (25), (31) при и> = 6, I = 0,1, определяем значения

С1,в = 1, С2,6 = 1, а1,1,в,в = 0, а2,1,б,б = 0, а2,2,б,б = 0.

Параметры метода С1,2 = а1, 01,3 = а2, 01,4 = аз, 02,3 = а4, 02,4 = а5, 0.1,1,3,2 = ав, а2,2,4,з = а^ могут быть выбраны в качестве свободных параметров, подчиненных лишь некоторым очевидным ограничениям.

Разрешая две линейные относительно Ьд<и, ц = 1,2, системы (10), (22), определяем параметры сд^, обеспечивающие их совместность и весовые коэффициенты

5 Сд,4 Сд,3 - 3 Сд,4 - 3 Сд 3 + 2

С-д 5 = ---------,

10 Сд,4 Сд,3 - 5 Сд,4 - 5 Сд,3 + 3 5 С2д3(10 С2д 4 - 8 Сд,4 + 1) - 5 Сд,3(8 С2д4 - 7 Сд,4 +1) + 5С24 - 5 Сд,4 + 1

Ьд,1

60 Сд,4 Сд,3(5 Сд,4 Сд,3 - 3 Сд,4 - 3 Сд,3 + 2)

5 С2 4 - 5 Сд,4 + 1

(у ^ == ___

Ъ 60 Сд,3 (1 - Сд,3 )(Сд,4 - Сд^)(3 Сд,4 + 6 Сд,3 - 2 - 10 Сд,4 Сд,3 + 10 Сд,4 С^ - 5 С^) '

ъ __5 с\ з - 5 съз + 1_

д' 60 Сд 4 (сд>4 - 1) (сд>4 - сд>3) (6 Сд 4 + 3 съз - 2 - 10 Сд 4 съз + 10 с2 4 съ3 - 5 с2 4)'

_ 5 4з(10 <?чЛ - 12 Сд 4 + 3) - 5 с,,3(12 <?чЛ - 15 од,4 + 4) + 15 - 20 од,4 + 6

9'6 60(5СдАСд<3-2СдА-2Сд<3+1)(1-СдА)(1~Сд<з)

Ьд,5 = 1 - (Ьд, 1 + Ьд,3 + Ьд,4 + Ьд,в).

Затем находим параметры, обеспечивающие совместность уравнений:

1) первой группы <1,1> первого блока (11)-(14):

_ -61,4(601,2(60^4 - 501,4 + 1)в1,1,з,2 + 01,3(501,301,4 - 01,4 - 01,з)) _ а1Д'4'3 6 (5с1,з - 5с1,з + 1) С1,з

2) второй группы <1,2> первого блока (15)-(21):

С1,4

«1,2,4,3 — -?-^- x

6 С2,3 (5 С213 - 5 С1,3 + 1)|(15 С1,3 - с23 + 4с2,3 - 3 С1Л (01,3 - С1,4)

х | (20 с1,3 - 30 с213 + 150 с1,4 с213 - 100 с1,4 с1,3 - 6 + 20 с1А) с23 -- 45 С1,4 с2,3 + 5 С13 + 20 С1,4 01,3 - 4 01,4 ^;

3) второй группы <1,2> первого блока и первой группы <2,1> второго блока (15)-(21), (23)-(29):

3 С1,3 - 2 С2,3 с2,2 = -

02,1,3,2

15 С1,3 02,3 - 10 02,3 + 2

02,3(3 01,3 - 2 С23)(15 С1,3С2,3 - 10 02,3 + 2)

02,1,4,4 =

2 012(45 с2,3 - 30 01,3 + 4)

с2,4(15 с1;з с2,4 с2,3 + 2 с2,3 - 10 с2,4 с2,3 ~ 3 с1;з + 2 с2,4)(с2,3 - с2,4) . 6 С1А(сЛз - с1,4)(5 с|3 - 5 с2,з + 1)

4) второй группы <1,2> первого блока, первой группы <2,1> второго блока и второй группы <2,2> второго блока (15)-(21), (23)-(33):

15 с2,302,3 - 10 01,302,3 + 2 С2,3 - 01,3

«2,2,3,3 — -:-тг~о--п-'

45 с2 3 - 30 с2,3 + 4

а2,2,4,4 = {(600 с2,3 - 270 с2,3 - 2025 с2,3с2,3 - 6750 4,3 с2,3 - 3900 с3,3 с2,3 -

- 520 с2,3 - 9000 с213 с2,3 - 1200 с2,3 с2,3 + 9000 с2,3 с2,3 + 6750 с2,3 с313 + 120 с13 +

+ 160 С2,3 + 3600 г^з С1,3 + 3150 с2,3 02,3 - 16)4 4 + + (630 с2,3 С2 3 - 180 с2 3 - 120 С1,3 С2,3 - 900 с2 3 с2 3 + 2025 с2 3 с2 3 - 120 ¿\\ 3 01,3 + + 405 с31 3 - 1350 с23 с31 3 - 320 с13 - 1350 022,3 с213 + 36 с13 -- 24 С2,3 + 160 с2,3 + 1200 с33 с^с2^ + + (90 с2,3 С2,3 - 135 С2,3 с3 3 - 60 с3,3 С13 + 24 01,3 02,3 + 180 с2,3 с2,3 - 120 с2 3 С13 + + 40 с32 3 - 18 с213 - 8 с2,3)с2,4 + (5400 с33 с13 - 720 с33 + 96 с22 3 + 4050 с22 3 с3,3 -

- 1620 с2,3 с2,3 + 1680 с\ 3 - 27 000 с5, 3 с21 3 - 1200 ¿2,3 - 360 с2 3 с13 - 72 с13 с23 -- 810 с23 с31 3 + 540 с213 с2 3 - 20 250 с\ 3 с3 3 + 20 250 с5, 3 с313 - 8100 с13 с213 -

- 14 400 с4,3 с1 3 + 32 400 с4> 3 с21 3 + 10 800 ¿2,3 с13)а2,2,4,3} х

-1

6 с2,4((15 С13 с2,з + 2-10 с2,з) с2,4 + 2 с2,з - 3 01,3)

1

(5 с23 + 5 С2,3 + 1)(45 С2 3 - 30 01,3 + 4)'

Из рассмотрения уравнений, связывающих параметры первой группы

<1,1> первого блока, находим их, последовательно разрешая линейные системы:

5) (11) (при ад = 2, I = 0,1):

1

«1,1,2,1 — «1,1,2,2 — 2С1'2'

6) (12) (при т = 5, I = 0,1):

1(л \ 1Ь1>5п \

«1,1,5,5 — т711 — с1,5«1,1,6,5 — 7г 7-ц — с1,5),

2 2 Ь1,6

X

7) (11) (при т = 3, I = 0, 1):

1 С1,з - С1,2 1 С1,2 «1,1,3,1 — оС1,3 — «1,1,3,2-, «1,1,3,3 — оС1,3 — «1,1,3,2-;

2 С1,з 2 С1,з

8) (12) (т = 2, I = 0,1) и (13):

»1,з (С1,5 - С1,з)(1 - С1,з)

«1,1,4,2 — 7- • т-77-тт«1,1,3,2,

»1,4 (С1,5 - С1,4)(С1,4 - 1)

»1,з (С1,4 - С1,з)(1 - С1,з)

«1,1,5,2 — — - ' 7-гр-7Т«1,1,3,2,

»1,5 (С1,5 - С1,4)(С1,5 - 1)

»1,з (С1,4 - С1,з)(С1,з - С15) «1,1,6,2 — 7- ' -77-77|-ч-«1,1,3,2;

»1,6 (1 - С1,4)(1 - С15)

9) (11) (I = 0,1, т = 4, 5, 6) и (12) (I = 0,1, т = 4):

(5 с2 4 - 5 С1,4 + 1)(С1,2 - С1,з) С1,4

аМЛ1 " '(543-5^3 + 1)43-аМА2 +

с1,з - с1,4 + с1,з С1,4(15 4 з - 20 01,3 + 5 см + 3) 6 (5с1,3 _ 5с1,з + 1)с1,з 2 С1,4 - 10 С1,4 С1,з - С1,з + 15 С 1 з С1,4

«1,1,4,4 — -„ , 2

«1,1,5,4

6 (5 4з - 5 С1,з + 1) »1,4 (С1,4 - 1)(2 «1,1,4,4 + С1,4 - 1)

«1,1,6,4 =

2 »1,5 (1 - С1,5 )

»1,4 (2 «1,1,4,4 (С1,4 - С1,5) + (2С1,5 - С1,4 - 1)(1 - С1,4))

«1,1,5,з =

2 »1,6 (С1,5 - 1)

с1,5 («1,1,5,2 с1,2 + «1,1,5,4 с1,4 + «1,1,5,5

«1,1,6,з

2 С1,з С1,з

1 «1,1,6,2 С1,2 + «1,1,6,4 С1,4 + «1,1,6,5 С1,5

2 С1,з С1,з

ш

«1,1,ш,1 = С1,ш - ^^ «1,1,ш,^, т = 5, 6.

и=2

Параметры «1,2,ш,„ второй группы <1,2> первого блока находим

10) из равенства (17) (при т = 5):

»2,5 (1 - С2,5)

«1,2,6,5 —-;-;

»1,6

11) из системы уравнений (15), (16) (при т = 3):

С1,з (2С2,2 ~ С1;3) _ 4з

«1,2,3,1 — -~-1 «1,2,3,2 —

2 С2,2 , , , 2 С2,2

Разрешая линейную систему уравнений (18), (19), (21) относительно вспомогательных переменных = ^ 1=^+1 b1,vc1,va1,2,v,^ = 2, 3,4), определяем их значения

& = 26=3 bl,vci,va1,2,v,2 =

(4223 2 3 3 \

_ (с2,5 + с2,3 с2,5 + с2,4 - С2,4 - С2,з С2,4 С2,5+ С2,3 С2,4 - С2,з - с2,б) &

С2,2 (С2,4-С2,2)(С2,3-С2,2)

5+15 с2,4 С2,з ~ 8 с2,3 - 8 с2,4 ^^

120 С2,2 (С2,4 - С2,2)(С2,3 - С2,2) '

& = 6=4 bl,vCl,va1,2,v,3 =

2 2 3 4 3 2 3

C2, 4 C2 5 C2, 2 + C2, 4 c2 5 - C2, 4 C^3 5 + 5 - C3 5 - C2 , 4 C2 , 5 C2 , 2 + C2 5 C2 , 2 - 5 C2, 2

-:-:-02,5 -

(C2,3 - C2,2 ) C2,3 (C2,4 - C23)

5+15 C2,4 C2,2 - 8 C2,4 - 8 C2,2

(35)

120 (02,3 - 02,2) 02,3 (С2,4 - 02,3) ' С4 = 26=5 Ь^01,^01,2,„,4 =

2 2 3 4 3 2 3

_ с2,5 с2,3 С2,2 + С2,3 с2,5 ~ с2,3 <^5 + с| 5 - - С2,5 С2,з С2,2 + С2,2 - С2,2 ^

(С2,4 - С2,з) (С2,4 - С2,2) С2,4

5+15 С2,з с2,2 - 8 С2,з ~ 8 с2,2 ^^

120 (с2,4 - С2,з) (С2,4 - С2,г) С2,4 '

12) из трех линейных систем: (17) (при т = 4) и (36); (17) (при т = 3) и (35); (17) (при т = 2) и (34) последовательно определяем параметры

С4 - Ь2,4(1 - 02,4)

«1,2,5,4 — -Г-/-ТТ->

Ь1,5(С1,5 - 1)

Ь2,4 01,5(1 - 02,4) - С4

«1,2,6,4 - -Г-/-ТТ->

Ь1,5(С1,5 - 1)

(1 ~ С1;4) «1,2,4,3 Ь 1,4 + £з + (с2,3 ~ 1)^2,3 &1,5 (С1,5 - 1)

(С1,5 ~ 4)«1,2,4,3^1,4 + С1;5(с2,з ~ 1)&2,3 + £з 6*1,6 (1 — с1,Б)

а1,2,5,3 = а1,2,6,3 =

а _ (ci,4 ~ l)bl,4C2,3«l,2,4,3__Ь__

Ь1,5С2,2(С1,5 - 1) &1,5(С1,5 - 1)

(1 - с1,з)с1,3 &1,3 + (1 - Cl,4)Ci,4 &1,4 2Ь1,5С2,2(С1,5 - 1)

(C1,5 - C1,4)b1,4C2,3a1,2,4,3 Ь

а1,2,6,2

b1,6C2,2(c1,5 - 1) Ь1,б(в1,5 - 1)

(С1,5 - <=1,3)43 &1,3 + (С1,5 - С114)С^4 »1,4 2»1,6С2,2(С1,5 - 1)

( С2 , з - С2 , 2 ) « 1 , 2 , 4 , з С21 , 4

«1,2,4,1 —--т;--Ь С14,

С2,2 2 С2,2

2

С2,з«1,2,4,з С1,4

«1,2,4,2 —---Ь ~-;

С2,2 2 С2,2

13) параметры «1,2,5,1, «1,2,6,1 из уравнений (15) (при т = 5, 6):

ш — 1

«1,2,ш,1 = С1,ш - ^^ «1,2,-ш,^, т = 5, 6. и=2

Параметры «2,1,ш,„ первой группы <2,1> второго блока определяем из уравнений (23)—(29), последовательно разрешая линейные системы:

14) (23), (24) (при т = 2, 3):

С2,2 (2 С1,2 - С2,2) с2,2

«2,1,2,1 — -, «2,1,2,2 — -,

С1,2 С1,2

(С1,2 - С1,з) С2,з(2 С1,з - С2,з) с2,з С1,2 «2,1,3,1 — -«2,1,3,2 Н--~-, «2,1,3,3 — т---«2,1,3,2!

С1,з 2 С1,з 2 С1,з С1,з

15) (25) (при т = 5, I = 0,1):

»1,5 (1 - С1,5)2 »1,5 (1 - С1,5) (1 + С1,5 - 2 С2,5) «2,1,5,5 — 777-77-7; «2,1,6,5 — -777-77-ч-5

2 »2,5 (1 - С2,5) 2 »2,6 (1 - С2,5)

16) (25) (при т = 2, I = 0,1) и (26):

_ »2,3«2,1,3,2 (С2,3 ~ 1) (С2,5 ~ С2,з) ' ' ' »2,4(С2,4 - 1) (С2,5 - С2,4)

»2,з«2,1,з,2 (С2,з - 1) (С2,4 - С2,з)

«2,1,5,2 —

«2,1,6,2 =

»2,5(С2,5 - 1) (С2,5 - С2,4) »2,3«2,1,3,2 (С2,5 ~ С2,з) (с2,4 ~ С2,з) .

»2,6(1 - С2,5)(1 - С2,4)

17) (25) (при т = 4, I = 0,1):

1 \(л ч , (1 ~ С1,4)2 ^ \

«2,1,5,4 — --7Т \ и — С2,4;«2,1,4,4»2,4--~-»1,4 ? ,

»2,5(С2,5 - 1) 2

1 Г/ ч , , (1 -С1,4)(1 + С1,4 -2С2,5)^ \

«2,1,6,4 — -77-п (С2,5 — с2,4/)«2,1,4,4»2,4 Н--71-»1,4 ? ,

»2,6(1 - С25)1 2 )

18) (24) (при т = 4, 5, 6, I = 0,1):

2

С2,4 С1,2 С1,4 «2,1,4,3 — ~---«2,1,4,2--«2,1,4,4,

2 С1,з С1,з С1,з

«2,1,-ш,3 = -| - ^ а2,1,ги,1/С1^ 1, № = 5, 6,

С1,3 ^ и=2,и=3 '

ш

а2,1,ш,1 = С2,ш - ^ , ^ = 4, 5, 6.

и=2

Параметры а2,2,ш,„ группы <2,2> второго блока находим, разрешая линейные системы уравнений: (30) (при ад = 2, 3, 4, 5, 6) и (31) (при ад = 5, 4, 3):

С2,2

«2,2,2,1 — «2,2,2,2 — ,

»2,3^2,3 - 1)(2 а2,2,3,3 + С2,3 - 1) , »2,4^2,4 - 1)а2,2,4,3

«2,2,5,3 - -777-77-^--1--Г-Тл-^-'

2 »2,5(1 - С2,5) »2,5(1 - С2,5)

»2,4(С2,5 - С2,4) , »2,3^2,5 - С2,3) , »2,3(1 - С2,3)(1+ С2,3 - 2 С2,5)

«2,2,6,3 — —;-77-«2,2,4,3 + —;-77-Г~ «2,2,3,ЗН--7Т,-Тл-\->

»2,6(1 - С2,5) »2,6(1 - С2,5) 2 »2,6(1 - С2,5)

»2,4(с2,4 ~ 1) ,

«2,2,5,4 — ТП-Тл-7 1 1 «2,2,4,4 + &1А - ^ [,

2 »2,5(1 - С2,5) I >

«2,2,6,4 = 7Т,-ТТ--Г {2 (с2 5 - С2 4)(22,2,4,4 + (1 - С2,4)(1 + С2,4 ~ 2 С2,б)

2 »2,5(1 - С2,5) I- ■>

«2,2,5,5 = тг(1 — с2,Б), «2,2,6,5 = 0 (1 — с2,б), 2 2 »2,6

1 ( С2 м

1 I С2,м

«2,2,м,1 = с2,м--< ~тт- + У^(С2,2 - с2}1/)а2}2}1Л}1/ 1, /X = 3, 4, 5, 6,

С2,2 2

' 4 и=3 '

— \~тг ~ ^«2,2,^Л, Я- = 3,4, 5, 6.

2,21 2 ±3 J

«2,2,м,2 1 г,

С2,2 2

' 4 и=3

Представленный выше алгоритм построения решения системы-следствия условий порядка метода (3)—(5) шестого порядка численного интегрирования системы (1), (2) является доказательством следующей теоремы.

Теорема. В рамках шестиэтапного метода (3)-(5) при ограничениях (6), (8), (9) существует метод шестого порядка (семипараметрическое семейство расчетных схем) численного интегрирования системы (1), (2).

Преимущество предложенного метода (3)-(5) численного интегрирования системы (1), (2) перед существующими в классе явных одношаговых методов типа Рунге— Кутты по характеристикам: порядок, число этапов — очевидно, так как известные методы Рунге—Кутты для получения приближения шестого порядка требуют как минимум семи стадий.

Тестирование. Для демонстрации работы алгоритма построения экономичных методов шестого порядка при значениях свободных параметров а\ = а2 = ^

1 7 4 1 14

«з = «4 = уд, «5 = «6 = -Ер «7 = 5| построена шестиэтапная расчетная схема метода ИКБбб (табл. 2). При выборе свободных параметров аи, V = 1,..., 7, в рамках очевидных ограничений руководствовались обеспечением выполнения требований на весовые коэффициенты 0 ^ ^ 1.

Тестирование расчетной схемы КК866 проводилось на решении двух задач: движении космического аппарата (КА) по орбите Аренсторфа и движении КА в окрестности точки либрации Кроме предложенной схемы ИКВбб тестирование проводилось для схем Бутчера [5], Цитураса [6] и Хаммуда [7].

Таблица 2. Шестиэтапная расчетная схема шестого порядка ИКБвб

,1 ,w ,д 11,2,та,9 £>1,та

0 0 0 7 90

1 6 1 12 1 12 1 6 0

1 4 7 120 1 6 1 120 3 32 6 32 16 46

1 2 9 20 8 6 26 16 1 12 1 28 16 32 16 224 2 16

3 4 3 40 3 6 3 20 1 4 1 8 61 448 6 32 46 112 6 64 16 46

1 32 106 32 36 124 106 1 7 4 7 0 93 392 126 66 136 392 446 1064 360 133 7 90

c2,w «2,1,та,д «2,2, та, д £>2,та

0 0 0 17 336

1 Б 2 26 3 26 1 10 1 10 0

7 16 98 676 77 226 448 676 1 90 7 18 1 16 76 224

4 Б 4 26 12 26 16 26 4 6 19 90 8 136 14 26 119 1360 276 912

1 6 17 1080 11 72 103 1080 3 20 19 360 19 378 811 2692 31 3360 11 3240 6 12 24 96

1 166 436 33 29 612 146 328 146 644 436 0 1783 3664 863 1666 261 1016 40469 74386 960 661 0 29 480

Задача 1: орбита Аренсторфа. Плоское движение КА с координатами (х1,х2) в гравитационном поле, создаваемом Землей (0,0) и Луной (1,0), описывается [8] системой обыкновенных дифференциальных уравнений

Х1 + М Х1 - М

Х\ = Х\ 2X2 — ——--¡Л-

Di

х2 = Х2 - 2xi - м'тг - М77", Di D2

D2

где

Di = ((xi + м)2 + ¿2)3/2,

0.012277471,

D2 = ((xi - n')2 + x2)3/2

M = 1 — ц.

При начальных условиях

¿(0) = (0.994,0),

¿(0) = (0, 2.00158510637908252240537862224)

КА движется по орбите с периодом Tper = 17.0652165601579625588917206249. При переопределении параметров: ¿i = xi, ¿2 = ¿2, xз = ¿2, ¿4 = ¿i система приобретает структуру, позволяющую применить предложенный экономичный метод интегрирования. На рис. 1 приведен график зависимости глобальной погрешности Err от длины шага интегрирования tstep на интервале [0,Tper].

-log10(&r)

Рис. 1. Результаты для задачи орбиты Аренсторфа

Задача 2: точка либрации. Плоское движение неуправляемого КА во вращающейся системе координат в окрестности точки либрации Li = (1,1) системы Солнце— Земля описывается [9] системой обыкновенных дифференциальных уравнений

Х1 = Х2 + У1, Х2 = -Х1 + У2,

У1 =8(xi - 1) + (y2 - 1), У2 = -4x2 - У1.

При начальных условиях

ж(0) = + ——-——e, O^j, y(0) = (0,1 + 6)

КА движется по периодической орбите вокруг Li. При переопределении параметров Xi = xi, Х2 = У2, Хз = Х2, Х4 = yi система обыкновенных дифференциальных уравнений приобретает перекрестную структуру и также может быть проинтегрирована предложенным методом. Для тестирования было выбрано значение е = ущ. На рис. 2 приведен график зависимости глобальной погрешности Err от длины шага интегрирования tstep на интервале [0,Tper].

Наклоны ломаных на обоих графиках показывают, что для всех тестируемых методов зависимость глобальной погрешности от длины шага имеет шестой порядок. Это означает, что RKS66 действительно соответствует заявленному шестому порядку точности. При этом методы-оппоненты требуют семи вычислений правой части СОДУ на шаге (метод Цитураса требует восьми, однако на каждом шаге последнее вычисление правой части совпадает с первым вычислением на следующем шаге, поэтому фактически их требуется семь), тогда как RKS66 — только шести.

-log10(£rr)

12

■■■+■■ RKS66 ---RK67 Butcher

11

V

К

10

4

5

3

7

6

9

8

-1

-0.8

-0.6 -0.4 -0.2 0

Рис. 2. Результаты для задачи движения КА вокруг L1

Литература

1. Olemskoy I. V. Fifth-order four-stage method for numerical integration of special systems // Comput. Math. Math. Phys. 2002. Vol. 42. P. 1135-1145.

2. Olemskoy I. V. Structural approach to the design of explicit one-stage methods // Comput. Math. Math. Phys. 2003. Vol. 43. P. 918-931.

3. Олемской И. В. Модификация алгоритма выделения структурных особенностей // Вестн. С.-Петерб. ун-та. Сер. 10. Прикладная математика. Информатика. Процессы управления. 2006. Вып. 2. С. 55-64.

4. Olemskoy I. V., Eremin A. S., Ivanov A. P. Sixth order method with six stages for integrating special systems of ordinary differential equations // 2015 Intern. Conference on "Stability and Control Processes" in memory of V. I. Zubov. SCP-2015. Proceedings. 2015. P. 110-113.

5. Butcher J. C. On Runge—Kutta processes of high order // Journal of the Australian Mathematical Society. 1964. Vol. 4. P. 179-194.

6. Tsitouras Ch., Famelis I. Th. On phase-fitted modified Runge—Kutta pairs of order 6(5) // Intern. conference of Numerical Analysis and Applied Mathematics: Extended Abstracts. Crete, 2006. P. 1962-1965.

7. Chammud (Hammud) G. M. A three-dimensional family of seven-step Runge—Kutta methods of order 6. Numerical Methods and Programming // Advanced Computing. 2001. Vol. 2. P. 159-166.

8. Хайрер Э., Нерсетт С., Ваннер Г. Решение обыкновенных дифференциальных уравнений. Нежесткие задачи / пер. с англ. И. А. Кульчицкого, С. C. Филиппова; под ред. С. C. Филиппова. М.: Мир, 1990. 512 c. (Hairer E., Norsett S., Wanner G. Solving ordinary differential equations.)

9. Шиманчук Д. В., Шмыров А. С. Построение траектории возвращения в окрестность колли-неарной точки либрации системы Солнце—Земля // Вестн. С.-Петерб. ун-та. Сер. 10. Прикладная математика. Информатика. Процессы управления. 2013. Вып. 2. C. 76-85.

Статья поступила в редакцию 22 февраля 2018 г.; принята к печати 14 июня 2018 г. Контактная информация:

Олемской Игорь Владимирович — д-р физ.-мат. наук, проф.; [email protected] Коврижных Николай Александрович — аспирант; [email protected]

A family of sixth-order methods with six stages

I. V. Olemskoy, N. A. Kovrizhnykh

St. Petersburg State University, 7—9, Universitetskaya nab., St. Petersburg, 199034, Russian Federation

For citation: Olemskoy I. V., Kovrizhnykh N. A. A family of sixth-order methods with six stages. Vestnik of Saint Petersburg University. Applied Mathematics. Computer Science. Control Processes, 2018, vol. 14, iss. 3, pp. 215-229. https://doi.org/10.21638/11702/spbul0.2018.303

The paper deals with construction of an efficient explicit sixth order method for structurally partitioned systems of ordinary differential equations. The general scheme of the method is presented, which algorithmically uses structural properties of a system of differential equations. The order conditions and the simplifying conditions for the considered explicit six-stage method are written down and their consistency is determined. The general solution with seven free parameters is obtained and a computational scheme for certain values of free parameters is constructed. Its performance on test problems is compared to three other explicit sixth-order methods.

Keywords: order, order conditions, simplifying conditions. References

1. Olemskoy I. V. Fifth-order four-stage method for numerical integration of special systems. Comput. Math. Math. Phys., 2002, vol. 42, pp. 1135-1145.

2. Olemskoy I. V. Structural approach to the design of explicit one-stage methods. Comput. Math. Math. Phys., 2003, vol. 43, pp. 918-931.

3. Olemskoy I. V. Modifikatsiya algoritma vydeleniya strukturnykh osobennostei [Modification of structural properties detection algorithm]. Vestnik of Saint Petersburg University. Series 10. Applied Mathematics. Computer Sciences. Control Processes, 2006, iss. 2, pp. 55-64. (In Russian)

4. Olemskoy I. V., Eremin A. S., Ivanov A. P. Sixth order method with six stages for integrating special systems of ordinary differential equations. 2015 Intern. Conference on "Stability and Control Processes" in memory of V. I. Zubov. SCP-2015. Proceedings, 2015, pp. 110-113.

5. Butcher J. C. On Runge—Kutta processes of high order. Journal of the Australian Mathematical Society, 1964, vol. 4, pp. 179-194.

6. Tsitouras Ch., Famelis I. Th. On phase-fitted modified Runge—Kutta pairs of order 6(5). Intern. conference of Numerical Analysis and Applied Mathematics. Extended Abstracts. Crete, 2006, pp. 19621965.

7. Chammud (Hammud) G. M. A three-dimensional family of seven-step Runge—Kutta methods of order 6. Numerical Methods and Programming. Advanced Computing, 2001, vol. 2, pp. 159-166.

8. Hairer E., Nersett S. P., Wanner G. Solving ordinary differetial equation I: Non-stiff problems. 2nd ed. Berlin, Heidelberg, Springer-Verlag Publ., 2008, 528 p. (Springer Series in Computational mathematics.) (Russ. ed.: Hairer E., Nersett S. P., Wanner G. Reshenie obyknovennykh differen-tsial'nykh uravnenii. Nezhestkie zadachi. Moscow, Mir Publ., 1990, 512 p.)

9. Shymanchuk D. V., Shmyrov A. S. Postroenie traektorii vozvrascheniya v okretnost' kollinearnoy tochki libratsii sistemy Solntse—Zemlya [Construction of the return trajectory to the neighborhood of the collinear libration point of the Sun—Earth system]. Vestnik of Saint Petersburg University. Series 10. Applied Mathematics. Computer Sciences. Control Processes, 2013, iss. 2, pp. 75-84. (In Russian)

Author's information:

Igor V. Olemskoy — Dr. Sci. in Physics and Mathematics, Professor; [email protected] Nikolay A. Kovrizhnykh — postgraduate student; [email protected]