Явный метод типа Рунге - Кутты пятого порядка Текст научной статьи по специальности «Математика»

Вычислительные технологии

Том 10, № 2, 2005

ЯВНЫЙ МЕТОД ТИПА РУНГЕ - КУТТЫ ПЯТОГО ПОРЯДКА

И.В. ОЛЕМСКОЙ Санкт-Петербургский государственный университет, Россия

e-mail: OLIV@pobox.spbu.ru

An explicit one-step method is suggested for numerical integration of the ordinary differential equation system of a special type. Economical numerical schemes of the fifth order are constructed for the integration of both systems and second-order differential equations of this special type. The algorithm of allocation of structural features of the special type is described.

1. Метод интегрирования

В [1-3] выделен класс систем обыкновенных дифференциальных уравнений:

УО = /о(х,Уо,Ш,.. .,Уп);

у' = Мх,Уо,... ,Уг-1 ,уг+1,... ,Уп), г =!,...,I;

Vj

fj(x,Vo,---,Vj-i), j = l + 1,

n,

(1) (2) (3)

x e [Xa,Xk] С R, Vs :[Xa,Xk ] Rrs, s = 0 ,...,n fo :[Xo,Xk] x Rm Rr0, fi : [Xo,Xk] x Rm-£S=iRri, fj :[Xo,Xk] x Rm-£rs Rrj,

i = 1,...,/,

j = l + 1,... ,n,

E rs = m.

s=0

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

m0(q) = q, ms(q) = q - 1, q < 4,

1 , . . . , n.

© Институт вычислительных технологий Сибирского отделения Российской академии наук, 2005.

n

s

Очевидно, что эффективность от применения этого метода при интегрировании систем

т I / т \

(1)-(3) тем выше, чем больше отношение У] / ( I. Здесь т3 — весовые коэф-

s=ro+1 ' \^=1 /

фициенты. Отсутствие общей группы уравнений (1) (го = 0) в рассматриваемой системе может еще больше повысить эффективность метода. Это было продемонстрировано в [4,5] для систем

У1 = /l(x,У2), (4)

У2 = /2(х,У1)

(Уг, /г — функции размерности гг, г = 1, 2). Например, в [4] для систем (4) построены четы-рехэтапные расчетные схемы пятого порядка, которые на треть экономичнее формального метода Рунге — Кутты.

В данной работе для систем более общего вида

Уг = /i(x,Уl,■■■,Уi-1,Уl+1,■■■,Уn), г =1,...,/, (5)

у'з = /3 (x,Уl,■■■,Уj-l), 3 = 1 + 1,■■■,n, ()

строится четырехэтапный метод пятого порядка. То есть так же, как ив [4], общая группа уравнений отсутствует (г0 = 0), но снято ограничение на число уравнений в группах

(2), (3), п > 3-

Приближенное решение ищем в виде

4

Уs (х + К) ^ Zs = Уs(ж) + ^ Ъзр кяр(Н), в = 1, ■ ■ ■ , п, (6)

р=1

где функции кар — к^(к) вычисляем в строгой последовательности

к11, к21, ■ ■ ■ , kn1, к12, к22, ■ ■ ■

по схеме

кгр — К/г(х + ciph, У1(х) + ^ ^ агр1пк1п, ■ ■ ■ , Уг-1(х) + ^ ^ аг,р,г-1,пкг-1,п,

П=1 п=1

р- 1 р- 1

У1+1(х) + X аг,р,1+1,пк1+1,п, ■ ■ ■, Уп(х) + ^ агРппкпЯ); (7)

П=1 П=1

рр

кю = К/(х + У1(х) + X] азр1пк1п, ■ ■ ■ , Уз-1(х) + X аз,р,з-1,пкз-1,п), (8)

П=1 П=1

сгр = 0, г = 1, ■■■,/, Су1 = 0, 3 = I + 1, ■ ■ ■ При построении расчетных схем в [3] было получено несколько важных соотношений,

связывающих параметры метода , Ъяр, csp■

с1р = ■ • • ■ = С1р — сгр1 с1+1,р = ■ • • ■ спр = Сjр ,

Ъ1р = ■ ■ ■ = Ъ1р = Ъгр, Ъ1+1,р = ■ • • ■ Ъпр — ЪjP ,

агр1п ■ • • ■ аг,р,г—1,п а, гргП' аг,р,1+1,п ■ • ■ ■ агрпп — aipjn'.

ajpln ■ • ■ ■ ajpln — ajpщ, aj,p,l+1,n ■ • ' ■ = aj,P,j-1,^n — а " зрзп

р

р

Использование ограничений (9) существенно упрощает (не в ущерб экономичности) построение условий порядка, вывод и реализацию расчетных схем. Руководствуясь этими соображениями, здесь будем изначально использовать равенства (9). Это позволяет представить метод (6)-(8) в компактной табличной форме (табл. 1).

Условия порядка для метода (6)-(8) пятого порядка даже с использованием равенств (9) и стандартных предположений

4=1

-■вр

(10)

представляют собой систему 102 нелинейных алгебраических уравнений с 50 неизвестными. Для поиска решения этой системы помимо упомянутых выше ограничений (9), (10) были использованы следующие упрощающие предположения:

1

С»V 2 Ciq,

1С2

2 '

Е

ajq'jv СIV

Е

v=2

q

ajqiv

v=2 4

V=v

4

Е*

V=V+1 V

Л 1 2

aivjv ^ (1

1

q— 1

^ ] aiqjv cjv — ^ '

V=1

1С2

2 с»'

Е'

v=2

3 ^'

Ьте (1 )

ajviv (1 ) ? ^ 2, . . . , 4,

(11)

'ш aivjv

bjv (1 - с^),

ajvjv

bjv (1 - с^) , V — 1,..., 4,

v=v+1

1,..., 4.

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

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

По характеристикам (число этапов, порядок точности на шаге) полученная в рамках структурного подхода расчетная схема (см. табл. 2), как и для системы (4), на треть экономичнее лучших существующих [6, 7].

Необходимость в интегрировании систем такого типа возникает, например, в задачах небесной механики, физики высоких энергий и т.д. Как следствие полученного результата (теорема 1.1) выпишем экономичную расчетную схему прямого интегрирования систем дифференциальных уравнений второго порядка вида

г" = дСМ^у'), у'' = / (М,^).

Приведенная стандартным образом х1 — г, х2 — г', х3 — у, х4 — у' к нормальной форме систем первого порядка (после переобозначения — х2, £2 — х3, £3 — х1, £4 — х4) система попадает в класс рассматриваемых (5).

(12)

р

2

V

Таблица1

Четырехэтапный метод (6)—(8) интегрирования системы (5)

cis aisiv aisjv bis

0 0 0 bii

Ci2 ai2ii ai2i2 ai2 j i bi2

Ci3 ai3ii ai3i2 ai3i3 ai3ji ai3j2 bi2

Ci4 ai4ii ai4i2 ai4i3 ai4i4 ai4j i ai4j2 ai4j3 bi4

cjs ajsiv ajsjv bjs

Cji ajiii aji5i bji

Cj2 aj2ii aj2i2 aj2ji aj2j2 bj2

Cj3 aj3ii aj3i2 ai3i3 aj3ji aj3j2 aj3j3 bj3

Cj4 aj4ii aj4i2 aj4i3 aj4i4 aj4ji aj4j2 aj4j3 aj4jj4 bj4

Применение для ее интегрирования метода (6)-(8) (см. табл. 1) и последующее сведение расчетной схемы к алгоритмически простому виду дают вычислительную схему четырех-этапного метода пятого порядка прямого интегрирования системы (12) вида

z(t + h) = z(t) + hz' (t) + hJ2Bipkip + O(h6), z' (t + h) = z' (t) + ^ КК + O(h6),

p=i p=i

y(t + h) = y(t) + hy' (t) + h^Bjpkjp + O(h6), y' (t + h) = y' (t) + X bjpkjp + O(h6),

p=i p=i

(13)

где функции kip, kjp из соображений явности вычисляются в строгом порядке kii, kj1, ki2, ... , ki4, kj4 по правилам

/ p—1 p—i p— i

kip = h g I t + Ciphj z + Ciph z + h 53 Aipidkidj y + Ciphy' + h 5] A ipjdkjd, y + aipjdkjd

V d=i d=i d=i

p p — i p

kjp = h f ( t + Cjph, z + Cjph z' + h Ajpidkid, y + Cjphy' + h Ajpjdkjd, z' + 53 ajpidkid I j

V d=i d=i d=i /

4 4 p—i

Bip 53 bjdajdipj Bjp 53 bidaidjpj Aipid 53 aipjsajsidj

d=p d=p+i s=d

p p p

Aipjd = 53 aipisaisjd, Ajpid =53 ajpjsajsid, Ajpjd =53 ajpisaisjd-

s=d+i s=d s=d+i

Значения параметров aipjd, ajpid, bip, bjp — в табл. 2, а параметров Aipid, Aipjd, Ajpid, Ajpjd, Ajpid j Ajpjd Bipj Bjpj Cipj Cjp в табл. 3.

4

4

4

4

Таблица2

Метод (6)—(8) пятого порядка интегрирования систем (5)

С-гв агзЪ aгsjv Ъгв

0 0 0 82 + 77У6 285+ 1140

4 У6 15 15 2 У6 15 30 2 У6 15 30 4 15 У6 15 297 351-У 1337 764

1 У6 2 8 19 19У6 160 640 9 9У6 32 128 1 У6 10 40 9 32" 9У6 128 7 7У6 32 128 2432 + 64У6 2415+ 345

7 + у6 10 + 20 19971 , 142933У6 29375 1 940000 64143 772839У6 41125 1316000 263168 , 110052У6 205625 1 205625 3 У6 10 20 4977 9400 4419У6 18800 2213 , 9809У6 9400 + 112800 61 , 4469У6 940 + 22560 18184 , 51676У6 250401 1 250401

ajsгv ajsjv

2 У6 15 30 2 \/6 15 30 2 \/6 15 30 0

2 У6 5 10 1 У6 10 40 3 3У6 10 40 3 3У6 10 40 1 У6 10 40 4 У6 9 36

2 + у6 5 + 10 1337 , 1947У6 1250 + 5000 4551 1083У6 1750 1000 8448 , 496У6 4375 + 625 6 + 3У6 25 + 200 17 + 27У6 50 + 200 3 10 У6 20 4 + у6 9 + 36

1 103 83У6 38 76 2901 , 11721У6 382 + 5348 72 272У6 23 161 62874 , 49236У6 83467 1 83467 3 + 3У6 8 + 8 1 У6 4 4 9 8 8 0 1 9

ТаблицаЗ

Метод пятого порядка (13) интегрирования систем (12)

Сгв А • • Aгsгv А • • Aгsjv Вг.з

0 4 У6 15 15 1 У6 2 8 7 + у6 10 + 20 0 11 4У6 225 225 209 19У6 231 21У6 2560 640 2560 640 4902303 | 1320067У6 14304369 5743641У6 840941 , 340324У6 9400000 1 9400000 13160000 13160000 1028125 1 1028125 11 4У6 225 225 363 33У6 77 7У6 2560 640 2560 640 34311 8373У6 55983 , 23087У6 5933 , 4713У6 188000 94000 376000 1 376000 75200 1 75200 82 + 77\/6 285 + 1140 927 1881У6 2674 5348 1552 , 176У6 2415 + 805 6986 | 16412У6 83467 1 250401

В

Й В

¡¿¡с

Е Й

н

Р

Р

Е

й

Й

Ы Т

ч о Й о

ы

о

К

А

С.

2 у/б 15 30

2_ У6

5 10

2 + У6

5 + 10

А

11 2У6

450 225

11 Ув

160 40

679 + 537л/6 800+ 1600

33 зУб 800 200

18281 + 3917л/6 57963 12591\/6 7464 + 3096\/6 100000+ 50000 140000 70000 21875+ 21875

2157 2211л/6 1120 2240

276 + 114У6 175 + 175

А

11 у/6 100 25

_13_7У6 _81_ + 17У6

250 250 500 + 250 _ 1 + 5У6 4 + 16

л/6

3 3\/6

4 16

В

3±_

1 + 4 + 36

1_ У6

4 36 0

со

0

1

0

2. Приведение систем к специальному виду

Структурный метод интегрирования систем т обыкновенных дифференциальных уравнений (СОДУ) самого общего вида

4 = к ,гт), к = 1,..., т, (14)

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

Задача состоит в том, чтобы указать порядок следования уравнений и нумерации переменных исходной системы (14), при котором для преобразованной системы

пг' = п^(ж,пг) (15)

с выделенными особенностями (1)-(3) эффект от применения структурного подхода будет

т т

максимальным, т.е. максимально отношение У] ип(8)/ ^ ип(8).

«=Г0 + 1 8=1

Здесь и в дальнейшем пг = (г^),^ ),..., ¿Л(т))Т, п^ = (^п(1),^п(2),..., <£п(т))Т, п = (п(1),п(2),... ,п(т)) — перестановка элементов множества 1т = {1, 2,... , т}, определяющая порядок следования уравнений системы (14). Например, п(г) = j означает, что j-е уравнение системы (14) будет занимать г-е место в системе (15). Далее, и8 — весовые коэффициенты, либо задаваемые пользователем, либо вычисляемые как относительные доли затрат (число арифметических операций, время) вычисления правой части от общего числа вычислений всех компонент вектор-функции ^ = (^1, ... , <^т), э = 1,... , т. Обозначим через Р = {п} множество всех перестановок из т элементов 1т = {1, 2,... , т}.

В [8] рассматривался частный случай решения этой задачи — алгоритм приведения к виду (1)-(3) в предположении выполнения равенств и1 = и2 = ... = ит. Здесь будет дано обобщение алгоритма [8] на случай неравноправных правых частей и1 = и^ УЬ, d £ /т, Ь = d.

2.1. Основные понятия

Определение 2.1. Матрицу А = {ач}т'^=1 будем называть структурной матрицей системы обыкновенных дифференциальных уравнений (14) и обозначать А(^), если

ч

0, если ^г не зависит от ,

1, если ^г зависит от .

Так, структурная матрица системы (1)-(3) имеет блочный вид:

(

А(/) =

* * * * * *

* °Г1 хгх 0Г1хП * * *

* * * О г\хг\ * * *

* * * * 0гг+1Хгг+1 °П+1Хгп

* * * * * * 0г„хг„

Здесь и в дальнейшем OriXrj — нулевая матрица размерности гг х г]. Символом * обозначены блоки, которые могут быть и ненулевыми.

Определение 2.2. Множество, элементами которого являются номера компонент искомой вектор-функции г = {¿1, ...,..., гт}, от которых не зависит правая часть 1-го уравнения СОДУ (14), будем называть 1-м горизонтальным структурным множеством СОДУ (14) и обозначать Ег(ф), т. е.

Ег(ф) = О € 1т

о а(ф) =

Определение 2.3. Множество, элементами которого являются номера уравнений СОДУ (14) с правыми частями, не зависящими от ]-й компоненты вектор-функции г = {г1,... , гт}, будем называть ]-м вертикальным структурным множеством СОДУ (14) и обозначать Н] (ф) :

Н] (ф) = {г € /„

ац = 0, А(ф) = }.

Пусть г = (го, г1,... , гп) € ({0} и N х Введем в рассмотрение

к-1

5(г, &) = X Гг, 5(г, 0) = 0, 5(г, п +1)

т.

г=0

Определение 2.4. Будем говорить, что СОДУ (14) имеет нуль-структуру с параметрами / € {0} и N п € N / < п, г = (г0, г1,... , гп) € ({0} и N х Ж", и обозначать ее если для горизонтальных множеств справедливы включения

{5(г, г) + 1,..., 5(г, / + 1)} С Ег(г,г)+М(г)(ф), г = 1,...,/, п>/> 0 ,

{5(г,^ ) + 1,..., 5(г, п +1)} С Ег(г,])+М(])(ф), г = / + 1,...,п, п>/ > 0 ,

где = 1,... , гв; 5 = 1,... , п.

Можно ввести это же понятие, используя определение структурной матрицы. Определение 2.5. Будем говорить, что СОДУ (14) имеет нуль-структуру с параметрами / € {0} и Ж, п € Ж, / < п, г = (г0 ,г1, ..., гп) € ({0} и N х Ып, и обозначать ее ^[/,п,г], если для элементов структурной матрицы А(ф) СОДУ (14) справедливы равенства

а^(з)и(в) = 0, 5 = 1,... ,п,

где

£(в) = 5(г,в) + 1,...,6(г,5 + 1) V (5) =

5(г, 5) + 1,... ,5(г, / + 1), если 0 </<п 5 < /,

5(г, 5) + 1,..., 5(г, п + 1), если 0 < /<п 5 >/.

Так, СОДУ (1)-(3) имеет нуль-структуру ^$т[/,п, г].

Структурная матрица СОДУ (14) с нуль-структурой ^5т[0,п,г], п € N, имеет вид

(

А(Ф)

* * *

* °Г1 хгх Ог1 ХГп

* * * Огпхгп

\

Определение 2.6. Будем говорить, что нуль-структура п, г] СОДУ (14)

предельная, и обозначать ее если справедливы невключения

{¿(г, 1),...,£(г,/ + 1)} СЕ^Ы, /> 0, {£(г,/ + 1),...,5(г,п +1)} С Ег(г>г+1)(^), / > 0.

Определение 2.7. Под объемом нуль-структуры п, г] СОДУ (14) будем пони-

т

мать величину ^ и и обозначать ее ^б^/, п, г]|, где и — весовые коэффициенты.

^=го+1

Можно показазать, что п, г] | > [/,п,г]|, п € Р.

Сформулируем с использованием введенных понятий две рассматриваемые здесь задачи.

Задача 1. Найти перестановку п* € Р такую, для которой справедливо неравенство

[0,п,г]| > |^5т^[0,п,г]| Уп € Р. Задача 2. Найти перестановку п* € Р такую, для которой справедливо неравенство

|^5тп^[/,п,г]| > Уп € Р.

Определение 2.8. Множество ..., С 1т будем называть элементным

нуль-структурным основанием СОДУ (14), если для В^ = и^=1 П и3 = 0, э = д,

справедливы включения С Ек(<^) Ук € ^, j = 1,... , d.

Определение 2.9. Для любого множества С С 1т величину ^ иг будем называть

г&О

весом множества С и обозначать V[G]:

V [С] = Е

= иг. гее

Пусть множество Ж С /т. Перенумеруем элементы этого множества с использованием чисел натурального ряда Ш = {г1, г2,... , гд}. Здесь и в дальнейшем |Ш| = д — мощность множества Ш. Введем в рассмотрение

8

Л,(г, э) = У гг, Л,(г, 0) = 0, ^(г^) = |Ш |, г = (г1 ,г2,...,г^) € N.

г=1

Теоpема 2.1. Для того чтобы множество В = {г1, г2,... , гд} было элементным нуль-структурным основанием ... , = г8, 5 = 1,... СОДУ (14), необходимо

и достаточно, чтобы

}С Егр Ы,

1 и/ ГММУ

Р = 1,..., д(г, d) — —2—

{Щr,S-1) + 1,Щr,S-1)+2, . . . , 4(^,8-1)+к} С ЕгНг,э-1) + к к = 1,..., г8, 5 = 1,..., d. Здесь [а] — целая часть числа а € Л.

Теоpема 2.2. Для того чтобы на множестве перестановок Р = {п} существовала такая п* € Р, что СОДУ (15) имела бы нуль-структуру ZSm?[0,п, г], необходимо и достаточно, чтобы существовало элементное нуль-структурное основание В?(^1,... , ип), для которого = Вт\В? (^1,...,^п), а |^в| = гв, 5 = 0,...,п. Причем |ZSm2[0, п, г]| =

V [В?] .

Теорема 2.3. Для того чтобы на множестве перестановок Р = {п} существовала такая перестановка п* € Р, что СОДУ (15) имела бы нуль-структуру ZSm2[/,n, г], необходимо и достаточно, чтобы существовали два взаимно непересекающиеся множества — элементныее нуль-структурные основания В? и В? (^1+1,... , ^п), для которых = 1т\ [В?(о;1,... , и В?(^1+1,... , ^п)] , |^в| = гв, 5 = 0,... , п. Причем ?[/,п,г]| = V [В?] + V [В?] .

Доказательство теорем 2.1 — 2.3 конструктивно [8]. На их базе и построен алгоритм выделения структурных особенностей и приведения произвольной системы к виду (1)-(3).

2.2. Алгоритм решения задач 1, 2

Строим множество П(1'0) = {г € 1т : г € Ег(ф)}. Если множество П(1'0) = 0, то Уп € Р, | ZS[/, п, г] | = 0. Это значит, что никакое изменение порядка следования уравнений исходной системы не может обеспечить выделения групп уравнений, имеющих структурные особенности. Представляет интерес случай, когда:

I. П(1'0) = 0. Полагаем значение номера шага алгоритма ^ = 0 и вводим в рассмотрение пару множеств В1 и В2, для определенности будем именовать их рекордными. В них будут храниться лучшие по суммарному весу на момент прохождения по дереву перебора упорядоченные множества, построенные в результате работы алгоритма.

II. При рассмотрении множества необходимо выделить два случая.

П.Л. Если множество П(1'м) = 0, то:

П.Л.1. Строим множества , г € П(1'м), з € П(1'м), г = з, являющиеся пересечением г-го горизонтального структурного множества Ег(ф) и з-го вертикального структурного множества Н] (ф) с множеством П(1'^):

Д((1,]/))(ф) = Ег(ф) п Н](ф) П г € з € .

11.Л.2. Затем строим множество возможных продолжений S(1'м) на ^-м шаге алгоритма при выбранных элементах в(0), в(1), ... , в(м-1), называемых в дальнейшем узловыми:

S= {7 = (г,з) € х = з : {г,з} С Я^}.

Отдельно следует оговорить несколько специальных случаев:

— так, при П(1'0) = {г1}, г1 € 1т, и |S(1'0)| = 0 упорядоченные множества определяются единственным образом: В1 = 0, В?2 = {г1};

— если П(1,0) = {г1,...,гк}, {г1,...,гк} € 1т, к > 1, и |(1,0)| = 0, то упорядоченные множества определяются: В?1 = {г*}, В?2 = {з*}, где г* € П(1'0) : тг* = шахгеП(1,о) тг; з* €

П(1'0)\{г*} : = шахг€П(1,0)\{г»} Ыг;

— для П(1'0) = {г1, г2}, {г1, г2} € 1т и S(1'0) = {(г1,г2)} упорядоченные множества определяются единственным образом: В?1 = 0, В?2 = {г1, г2};

— если П(1'0) = {г1, г2}, {г1, г2} € 1т, а S(1'0) = {(г1,г2), (г2,г1)}, то из двух равновесовых пар упорядоченных множеств выбираются, например, В?1 = 0, В?2 = {г1, г2}.

Все перечисленные случаи заканчиваются переходом на п. V алгоритма.

11.А.3. При рассмотрении построенного множества 5(1'м) возможны два случая.

П.А.З.а. Если 5(1'м) = 0, то проход по этой ветви закончен. Среди еще не рассмотренных элементов П(1'м) выбираем элемент г* (для определенности центральный), обладающий максимальным весом иг*. Прежде чем перейти к следующему этапу, необходимо сделать проверку на перспективность дальнейшего рассмотрения этого варианта

м-1

иг* + ^ {в^в^}

5=0

>1 (V В1 + V В2

2

г* € : и

тах иг

г€П(1>м)

(16)

В случае, если справедливо неравенство (16), строим упорядоченное множество В1. Количество элементов множества В1 = {г1, г2,... , ¿к} будет нечетно (к = 2^ + 1) и сформировано по правилу

¿5+1 = в(5), ¿к-5 = в(5), £ = 0,1,...,^ — 1, гм+1 = г* €

?(5)

иг* = тах иг.

г€П(1,м)

После чего поднимемся вверх по дереву перебора для построения второго упорядоченного множества В2 на п. III алгоритма.

В случае бесперспективности продвижения по этой ветви (неравенство (16) не выполняется) или отсутствия претендента на роль центрального элемента, опускаемся на один уровень вниз по дереву перебора. Для чего полагаем ^ := ^ — 1 и переходим на п. 11.А.3 алгоритма при ^ > 0.

П.А.З.Ъ. Если 5(1'м) = 0, то в качестве узлового элемента в(м) на ^-м шаге может быть выбран элемент в * € 5 (1'^), удовлетворяющий двум требованиям. Во-первых, среди не рассматривавшихся ранее в качестве узловых на ^-м шаге ему отвечает множество Б(1*'м), имеющее наибольший вес:

в * = (в*,в2*) € 5(1'^ : V

п(1'^)

> V

п(1'^)

Ув € 5(1'^).

Во-вторых, при проверке на целесообразность дальнейшего продвижения по ветви в(0), , в* справедливо неравенство

V

м-1

и{в1,в2}

5=0

+ V

Б

(1'М)

> 2 ^

В1

+ V

В2

(17)

Если неравенство (17) справедливо, то элемент в * становится узловым на ^-м шаге (уровне) алгоритма в(м) = в*. Формируем множество

п(1',+1) = Б^ив^^}

и, увеличивая номер уровня ^ := ^ + 1, переходим на п. 11.А.

В случае невыполнения любого из перечисленных выше требований, предъявляемых к претенденту на роль узлового элемента, возможны два варианта:

— при ^ > 0 необходимо вернуться на предыдущий уровень (для этого уменьшаем его номер на единицу (^ := ^ — 1) и переходим на п. 11.А.3 данного алгоритма);

— при ^ = 0 работа алгоритма по перебору закончена. Рекордные упорядоченные множества В1 и В2, определенные на этот момент, и являются искомыми. Переходим на п. V алгоритма.

*

*

11.В. Если множество П(1'м) = 0, то построение возможного варианта упорядоченного множества В1 закончено. Критерием целесообразности использования В1 и дальнейшей работы алгоритма по построению множества В2 является выполнение неравенства

V и{в!£,,в°}

5=0

> 2' v

B1

+ V

B2

(18)

При выполнении неравенства (18) формируется упорядоченное множество B1 = {i1, ... , ik}• Оно будет содержать k = 2ß элементов, которые определяются по правилу: ij+1 =

, ik—5 = , С = 0,1,... , ß — 1. Затем для построения второго упорядоченного множества B2 переходим на п. III рассматриваемого алгоритма.

В случае невыполнения неравенства (18) опускаемся по дереву перебора на один уровень ниже, полагая ß := ß — 1 и переходя на п. II.A.3.

III. Строим множество П(2>0) = n(1'0)\B1 при сформированном упорядоченном множестве B1, полагая значение номера шага алгоритма v = 0. Следующий п. IV алгоритма практически повторяет п. II с той лишь разницей, что имеет другие критерии отсечения и передачи управления.

IV. При рассмотрении множества необходимо учитывать два варианта перебора.

IV.A. Если множество Q(2'v, = 0, то:

IV.A.1. Строим множества D^», i G Q(2'v), j G Q(2'v), i = j, являющиеся пересечением i-го горизонтального структурного множества ЕДф) и j-го вертикального структурного множества Hj (ф) с множеством

D((j)^) = ад n Hj(ф) П Q(2'v), i G Q(2'v), j G Q(2'v).

IV.A.2. Строим множество возможных продолжений S(2'v) на v-м шаге алгоритма при выбранных узловых элементах а(0), а(1), ... , a(v-1):

S(2'v) = {7 = (i, j) G Ü(2'v) x H(2'v,,i = j : {i, j} С D2'v>}.

Отдельно следует оговорить несколько специальных случаев.

Если множество П(2'0) = {j} состоит из одного элемента, то дальнейший перебор не имеет смысла. Суммарный вес построенных упорядоченных множеств B2 = {j} и B1 имеет максимальное значение. Полагаем B1 = B1, B2 = B2 и переходим на п. V.

При П(2>0) = {j

1,... , jk}, k > 1, и

S(2>0) = 0 из элементов множества выбираем j * G

B1

+ Wj*

>V

B1

+ V

B2

, производим

П(2'0) : и]» = шах]еП(2,о) и]. В том случае, если V

замену рекордных множеств В?1 = В1, В?2 = {з*}. После проверки на максимум переходим на п. II.А.3 для построения нового упорядоченного множества В1. Причем, если число элементов В1 четно, необходимо понизить уровень первого дерева (^ := ^ — 1).

При П(2'0) = {з'1,з'2} и S(2'0) = {(з1,з2)} полагаем В?1 = В1, В2 = {з1,з2} и переходим на п. V.

При П(2>0) = {з\,з2} и ^2>0) = {(зЧ,з'2), (з'2,з'1)} полагаем В1 = В1, В2 = {з^} и переходим на п. V.

IV.A.3. При рассмотрении построенного множества необходимо выделять два

случая.

IV.A.3.a. Если S(2'v)

0, то продвижение по этой ветви закончено. Прежде чем

перейти к следующему шагу алгоритма, необходимо сравнить вес полученных упорядоченных множеств с весом рекордных

V- 1

V

В1

+ Wi

+ V и {а?

(?) Л?)

а

> V

В1

+ V

В

(19)

?=0

где г* € : и>г* = шах^Если справедливо неравенство (19), то необходимо про-

извести замену рекордных упорядоченных множеств В1 = В1 и В2. Количество элементов множества В2 = {г1, г2,... , ¿к} будет нечетно (к = + 1) и сформировано по правилу

г?+1

а

(?)

гк—?

а

(?)

£ = 0,1,

V - 1,

(■V+1

г* <Е :

:

шах wi

геп(2>^)

После проверки на успешность полученного варианта (19) независимо от результата сравнения осуществляем возврат по дереву перебора, для чего полагаем V := V — 1 при V = 0 с дальнейшим переходом на п. 1У.А.3 алгоритма. В случае же, если V = 0, переходим на п. 11.А.3 для построения нового упорядоченного множества В1. Причем, если число элементов В1 четно, необходимо понизить уровень первого дерева (^ := ^ — 1).

1У.Л.3.Ъ. Если Б(2'^ = 0, то в качестве узлового элемента а(^ на v-м шаге алгоритма может быть выбран элемент а* € Б(2^), удовлетворяющий двум требованиям. Во-первых, среди не рассматривавшихся ранее в качестве узловых на V-м шаге ему отвечает множество Ва'^, имеющее наибольший вес:

а* = (а*, а*) € Б(2^) : V

В

(V)'

> V [В<2^] У7 € Б(2'

V)

Во-вторых, при проверке на целесообразность дальнейшего продвижения по ветви а

(0)

а

(1)

а

(V-!)

а*

V

В1

В

(V)'

> V

В1

+ V

В2

(20)

справедливо неравенство

V— 1

"+ ^ {а1?),а2?)} + V ?=0

При выполнении этих требований элемент а* становится узловым: а(^ = а* на v-м шаге (уровне) алгоритма.

Для перехода на следующий верхний уровень дерева перебора формируем множество ^(2^+1) = вО)2('1Уг))\{а1гУ), а2 } и, увеличивая номер уровня (V := V + 1), переходим на п. 1У.А.

В случае невыполнения неравенства (20) при V > 0 необходимо вернуться на предыдущий (нижний) уровень. Для этого уменьшаем его номер на единицу (V := V — 1) и переходим на п. 1У.А.3 данного алгоритма.

В случае невыполнения неравенства (20) при V = 0 работа алгоритма по поиску наилучшего упорядоченного множества В2 при фиксированном В1 закончена. Необходимо перейти на п. 11.А.3 для построения нового упорядоченного множества В1. Причем, если число элементов В1 четно, необходимо понизить уровень первого дерева (^ := ^ — 1).

Если множество = 0, то необходимо иметь в виду, что:

При V = 0 на множестве 1т не существует упорядоченного непустого множества В2 с заданными свойствами. При В1 = 0 рекордными становятся В1 = В1, В2 = 0, и переходим на п. У.

]^.В.2. При V = 0 проход по этой ветви закончен. При выполнении неравенства

V [В Ч + V

11^1?

?=0

а

)}

> V

В1

+ V

В2

(21)

*

1

2

*

а

V

необходимо произвести обновление рекордных упорядоченных множеств В?1 = В1 и В?2 = {¿1,... , ¿к}. Причем В?2 будет содержать к = элементов, которые определяются по правилу г^+1 = ¿к- = £ = 0,1,...,^ — 1. После сравнения весов множеств и обновления рекордных упорядоченных множеств (если это потребуется) опускаемся по дереву перебора на один уровень ниже, полагая V := V — 1 и переходя на п. ^.А.3 алгоритма.

Работа алгоритма по пп. 1-^ будет закончена естественным образом на нулевом уровне дерева перебора при бесперспективности продолжения перебора с построенными упорядоченными множествами В?1 и В?2. Только после этого выполняется следующий пункт алгоритма.

Замечание 3.1. Элементы любого упорядоченного множества В5 = {¿1 ,г2,...,гк}, в = 1, 2, построенного по пп. 1-^ данного алгоритма, удовлетворяют условию (2.5) теоремы 2.1. Покажем это. Действительно, по построению справедливы следующие включения:

{¿1,¿к} с Еп (ф) п Щк (ф) п >0) = Ь>0), 0(* ,1) = Ь, 0)\{^1 ^ },

{¿2, ¿к-1} с е12 (ф) п нгк_, (ф) п '1) = Ь'1),

П(5'2) = Ь^{¿2, ¿к-1};

........................ (22)

{¿4(5), ¿к+1-^)} С Е^.)(ф) П Нгк+1_((.)(ф) П П(^)-1) = Ь(^)-1),

0(в'4(в)) = Ь (5'4(5)-1)\{г4(в),гк+1-4(5)},

0(5,4(5)) = Г 0, к = 2*(в),

\ {¿4(^+1}, к = 2^(в) + 1.

И в силу того, что Ь(8>4(8)-1) с Ь(8>4(в)-2) с ... Ь(5,1) С Ь(5'0), справедливы включения {¿р, ¿р+1,... , ¿к+1-р} С Егр(ф) П Нгк+1-р(ф), ¿(1) = ¿(2) = V, р = 1,... ,ф). Что и требовалось доказать. □

Замечание 3.2. Для любого упорядоченного множества В5 = {¿1,¿2,... , ¿к}, в = 1, 2, построенного по пп. 1-^, можно провести разбиение на классы.

Для этого полагаем, что ¿1 € Ш1 Элемент ¿2 будет принадлежать множеству Ш1, если ¿1 € Ег2 (ф). Если же ¿1 £ Ег2 (ф), то = {¿1} и ¿2 € ш2.

Предположим, что {¿1,г2} С Ш1. Тогда, если {¿1,г2} С Егз(ф), то ¿3 € Ш1, иначе (если {¿1, ¿2} С Егз (ф)) ¿3 € Ш2

Допустим, что уже произведено разбиение множества В5 = {¿1, ¿2,... , ¿к} до элемента ¿и^-1+а включительно, т. е.

Ш1 = {¿1, ¿2, . . . , ¿и1}, ^2 = {¿и1+1, ¿2, . . . , ¿и2 }, ..., Ш«-1 = {¿и^+Ъ...^ },

и установлено, что {¿и?-1+1,... , ¿и5-1+^} С ш«.

Элемент ¿и5_1+^+1 упорядоченного множества будет принадлежать этому же классу ¿и^-1+й+1 € ш«, если справедливо включение {¿и5-1+1,... , ¿и^_1 +Л С Еги?_1+^+1 (ф). После чего естественным образом переходим к рассмотрению следующего элемента ¿и5_1+^+2 на принадлежность множеству ш«.

В случае же, если {¿и?_1+1,... , ¿и5-1+^} С Еги 1+^+1 (ф), формирование класса ш« закончено и ш« = {¿и?_1+1,... , ¿и?}, и« = и«-1 + ¿и?+1 € ш«+1. После этого проверяем на

принадлежность множеству ^?+1 элемент ги?+2. И так до тех пор, пока не исчерпаем все элементы множества В

Для проведенного разбиения упорядоченного множества В5 на классы

справедливо представление В5 = ^=1 причем ^ С (<^) Уд € й = 1,... ,р(з), ^ П = 0, £ = д, которое обеспечивает выполнение условия (2.6) теоремы 2.1. А это значит, что все построенные упорядоченные множества В5 являются (с учетом замечания 3.1 и утверждения теоремы 2.1) элементными нуль-структурными основаниями

V. Руководствуясь правилом замечания 3.2, проведем поочередно разбиение рекордных множеств на классы. Сначала сделаем это для упорядоченного множества В2 = В2 (^1,^2,... , а затем для В1 = В*(^+1,^+2,... ,

Замечание 3.3. Результатом работы алгоритма является построение двух взаимно непересекающихся элементных нуль-структурных оснований В1, В2 — пары, имеющей максимальный суммарный вес.

Действительно, были рассмотрены пары взаимно непересекающихся упорядоченных множеств В1, В2 на множестве 1т. Из рассмотрения были исключены лишь те, которые заведомо не могли дать максимальной суммы (ограничения (19)-(21)), а также те, которые не удовлетворяли условиям (16)-(18).

Условия (16)-(18) выражают симметричность построения В1 на множестве 1т и В2 на множестве 1т\В1.

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

V [ВЧ + V [В2] > V

V [В1 ] < 2 (V

в1

в1

V

+ V

в2

в2

(23)

Но в таком случае неравенство (23) справедливо только при

V [В2] > 2

V

в1

+ V

в2

А это значит, что упорядоченное множество В2 С 1т\Р1 С 1т уже рассматривалось при построении первого упорядоченного множества. Полученное противоречие и доказывает утверждение.

Замечание 3.4. Руководствуясь утверждениями теоремы 2.2 (для решения задачи 1) и теоремы2.3 (для решения задачи 2), для любой пары В^(^+1,^+2,... ,^га) и В2(^1, ... , ^г) элементных нуль-структурных оснований (полученной в результате работы пп. 1-У алгоритма) можно построить перестановку п € Р, на которой СОДУ (15) имеет нуль-структуру ^Б^/, п, г], / € N Для этого необходимо перенумеровать элементы классов

^ = {^й(Г'5) + 1, . . . ,^й(Г'5)+Г3 }, г5 = 5 =1,...,п,

элементных нуль-структурных оснований В2 (^1, ... , ) и В2 (^г+1, ^г+2,... , и мно-

, дГ0}, где г0 = |^о|. После чего искомая перестановка

жества = 1т\ (В1 и В2) = {^1,.. представима в виде

1 2 ... Г0 Г0 + 1 #1 #2 ... #го £<5(гМ)+1

п

¿(г,/ + 1) ¿(г,/ + 1) + 1 ... т

^й(Г'1+1) ^й(Г'1+1) + 1 . . . д<5(г-'П+1)

Покажем, что ZS^^n, r] для любой пары нуль-структурных оснований B1, B2, сформированной пп. I-IV алгоритма, является предельной.

Предположим, что это не так, т.е. справедливы включения

, • • • ^¿(r.i+i)} с EÖÄ(ril) (ф) , 1 > 0, {g¿(r>г+l), • • • с EöÄ(r,i+1) (ф) , 1 > о.

А это значит, что как для упорядоченного множества B1 = |^г(Г)г+1), g5(r,i+1)+1, • • • , g5(r,n+i)}, так и для B2 = {g5(r,1), g5(r,1)+1, • • • , g5(r,i+1)-1} выполняются включения алгоритма (22). И в силу того, что V[B 1j > VfB1] и V[B1 ]+V[B2] > V[B 1]+V[B2], пара B1, B2 вообще не могла быть построена в рамках приведенного алгоритма. Полученное противоречие и доказывает предельность нуль-структуры. Таким образом, нуль-структура ZS^/, n, r], / £ N, на любой перестановке п, построенной по приведенному выше правилу, является предельной.

VI. Руководствуясь замечанием 3.4, на базе рекордных элементных нуль-структурных оснований B1, B2 строим искомую перестановку п*, которая является решением задачи 1 или задачи 2.

Отдельно следует отметить три важных случая.

1. Суммарный вес построенной пары упорядоченных множеств B1 и B2 равен максимально возможному: V[B1] + V[B2] = V[Q(1'0)]• Дальнейший перебор улучшения дать не сможет. Рекордными делаем множества B1 = B1 и B2 = B2 с дальнейшим переходом на п. V алгоритма.

2. Вес упорядоченного множества B2 равен максимально возможному: V[B2] = ^^[О(2'0)] • Дальнейший перебор на втором дереве улучшения дать не сможет. Переходим на п. II.A.3 для построения нового упорядоченного множества B1^ Если число элементов B1 четно, то необходимо понизить уровень первого дерева (ß := ß — 1)

3. Для использования алгоритма при решении задачи 1 необходимо считать постоянно B1 = 0 и B1 = 0 . Начать работу алгоритма следует с п. III и ограничить только перебором по верхнему дереву, т. е. требование спуститься по дереву перебора на поиск нового упорядоченного множества B1 означает (при решении задачи 1) окончание перебора с переходом на п. V.

2.3. Пример применения алгоритма выделения структурных особенностей

Продемонстрируем работу алгоритма решения задачи 2 на примере выделения нуль-структуры ^* ^

[/,п, г] максимального объема СОДУ

У

F(x,y), y, F £ R1.

(24)

W5

Структурная матрица исходной системы ) с весовыми коэффициентами и2 = и3 = 2, = и6 = 3, и = и7 = 4, и4 = 5 имеет нуль-структуру [1, 2, г], г = (5,1,1). Ее объем

|ZS7f[1, 2, r] | = 7-

A(F)

0101010 1 0 0 1 0 0 0 0010111 1 0 1 0 0 0 0 1110 0 10 111110 1 1001010

A(n*F)

1001011 1 0 0 0 1 0 0 0100100 1110 111 0 1110 0 0 0101100 \1011100/

Перестановка (переобозначение) п* = ^з426175^ найдена за восемнадцать шагов алгоритма. Его работа схематично представлена в табл. 4. Структурная матрица преобразованной системы А(п*Р) имеет нуль-структуру ZS7П»F[3, 5, г], г = (1,1,1,1, 1, 2) с объемом ^[3, 5, г] | = 21. Эффективность применения структурного подхода при интегрировании исходной системы определяется отношением

т т

X ws^ws = 7/23,

«=Г0 + 1 8=1

а преобразованной —

m m

X wf(s)^№n(s) = 21/23.

s=ro+1 s=1

3. Тестирование

Для сравнения была взята задача из [7]. Во-первых, она подходит по наличию структурных особенностей рассматриваемого типа и позволяет продемонстрировать не только преимущество полученной схемы, но и работу предложенного здесь алгоритма. Во-вторых, на этой задаче уже проводилось тестирование, и среди испытуемых методов одним из лучших был метод Дормана — Принса. Эти обстоятельства и определили выбор как тестовой задачи

у1 = 2хУ2/5 У4 = /1(х,У2,У4), = 10х ехР (5(Уз - 1)) У4 = /2(х,Уз,У4), (25)

= 2ху4 = /а (Х,У4), = -2x ln У1 = /4(x,yi)

yi(0) = 1, x G [0,10], F =(/i,...,/4)T,

так и оппонента — метод Дормана — Принса. Точное решение задачи Коши (25) имеет вид

y1(x) = exp (sin x2), y2(x) = exp (5 sin x2), y3(x) = sin (x2) + 1, y4(x) = cos (x2).

Весовые коэффициенты задаем, полагая w1 = w2 = w4 = 10, = 1. В соответствии с определениями и классификацией, введенными выше, система (25) уже имеет нуль-структуру ZS 4 f[1, 2,r], r = (2,1,1) с объемом |ZS4F[1, 2,r]| = 11. Отметим, что наличие общей группы уравнений сразу не позволяет применить метод, предложенный здесь (см. табл. 2). Для интегрирования системы (25), имеющей такой вид, могут применяться экономичные расчетные схемы структурного метода [3].

Используя алгоритм выделения структурных особенностей исходной системы (25), можно построить четыре перестановки, на которых преобразованные системы имеют нуль-структуры с равным максимально возможным объемом. Рассмотрим две из них:

ZSVf[2, 3, r1 ], r1 = (0,1,1, 2) и ZS4n2f[1, 3, r2], r2 = (0, 2,1,1), |ZSVf[2, 3, r1]| = |ZSVf[1, 3, r2]| = 31

1 _ ,1234 \ 2 ( 1234

на перестановках п = '4213/ и = \ 3 1 4 2 соответственно.

Таблица4

Схема работы алгоритма при решении задачи 2

N |5 )| в )/г* V[Dв7«* ]М* Вр У[ВР] 2 £ У[В<*] ¿=1

1 ^ = 0 ^(1,0) = {1, 2, 4, 5(1>0) = {(4, 5), (4, 7), (5, 7), 13 в(0) = (4, 5) 14 0

5, 6, 7} (1, 5), (1, 7),...,(2,6)} 0

2 ^ = 1 П(1>1) = {2, 7} 5(1'1) = {(2, 7), (7, 2)} 2 в(1) = (2, 7) 6 0

3 ^ = 2 0(1>2) = 0 5(1>2) = 0 0 В1 = {4, 2, 14 0

7, 5}

4 V = 0 О(2'0) = {1,6} 5(2'0) = 0 0 г* = 1 4 В2 = {1} 4 18

5 ^ = 1 п(1'1) = {2, 7} 5(1'1) = {(2, 7), (7, 2)} 2 в(1) = (7, 2) 6 В1 = {4, 7, 2 51 14 18

6 V = 0 О(2'0) = {1,6} 5(2'0) = 0 0 г* = 1 4 2, 5} В2 = {1} 4 18

7 ^ = 0 О(1'0) = {1,2,4, 5(1>0) = {(4, 5), (4, 7), (5, 7), 13 в(0) = (4, 7) 14 18

5,6, 7} (1, 5), (1, 7),...,(2,6)}

8 ^ = 1 П(1'1) = {2, 5} £(М) = {(2, 5)} 1 в(1) = (2, 5) 6 18

9 ^ = 2 0(1>2) = 0 5(1'2) = 0 0 В1 = {4, 14 18

2, 5, 7}

10 V = 0 О(2'0) = {1,6} 5(2'0) = 0 0 г* = 1 4 В2 = {1} 4 18

11 ^ = 0 О(1'0) = {1, 2, 4, 5(1'0) = {(4, 5), (4, 7), (5, 7), 13 в(0) = (5, 7) 12 18

5, 6, 7} (1, 5), (1, 7),...,(2,6)}

12 ^ = 1 П(1>1) = {4} 5(1>1) = 0 0 г* = 4 5 В1 = {5, 4, 7} 12 18

13 V = 0 О(2'0) = {1, 2, 6} 5(2'0) = {(2,6)} 1 а(0) = (2,6) 5 18

14 V = 1 0(2'1) = 0 5(2'1) = 0 0 18

15 ^ = 0 О(1'0) = {1, 2, 4, 5(1'0) = {(4, 5), (4, 7), (5, 7), 13 в(0) = (1, 5) 11 18

5, 6, 7} (1, 5), (1, 7),...,(2,6)}

16 ^ = 1 П(1>1) = {7} 5(1>1) = 0 0 г* = 7 4 В1 = {1, 7, 5} 11 18

17 V = 0 О(2'0) = {2, 4, 6} 5(2'0) = {(4, 6), (4,2), (2,6)} 3 а(0) = (4,6) 10 18

18 V = 1 0(2'1) = {2} 5(2'1) = 0 0 г* = 2 2 В2 = {4, 2, 6} 10 21

"1 = {4}, ^2 = {2}, = {6}, = {1}, = {7, 5}, ^0 = {3} п* = (3, 4, 2, 6,1, 7, 5).

В

Й В

Е Й

н

Р

Р

Е

й

Й

р

Ы

Т

р

о Й о

Ы

о

К

А

о

со

Т а б л и ц а 5

Зависимость максимальной глобальной погрешности от величины шага

- №) - !ё У тахж€[0,10] 1 Уг(х) - У г | |

ПК 4 БОРЫ 5(4)7Р РС (п1) РС (п2)

2.000 -1.3229 -0.3749 -0.4042 -0.4593

2.500 1.5100 1.9966 2.0944 2.0439

3.000 2.9692 4.4891 4.5915 4.5412

3.500 4.8711 6.9900 7.0907 7.0407

4.000 6.8443 9.3635 9.6442 9.7863

4.500 8.8857 9.7154 9.9318 9.6638

5.000 9.6123 9.8461 9.6891 9.7538

Ниже приведены структурные матрицы исходной системы и преобразованных на соответствующих перестановках:

А(Р)

0 1 0 1

0 0 1 1

0 0 0 1

1 0 0 0

( о

1 1

V 1

0

1

о о

А(п2Р)

о о

0

1

0

1

о о

Изменение порядка следования уравнений (на перестановках п1 и п2), во-первых, позволило повысить эффективность структурного метода до максимально возможной, во-вторых, дало возможность использовать для интегрирования преобразованных систем полученную здесь расчетную схему (см. табл. 2). Интегрирование производилось на интервале х £ [0,10] с постоянным шагом к. В серии расчетов шаг интегрирования к изменялся от 10-2 до 10-5, причем каждое следующее значение было в \/10 раз меньше предыдущего. Помимо метода Дормана — Принса [7] (БР5(4)7) для демонстрации точностного превосходства полученного метода в серии вычислений использовался и четы-рехэтапный (равнозатратный полученному) классический метод Рунге — Кутты четвертого порядка (ИК4). Поскольку для исходной системы рассмотрены две нуль-структуры — [2, 3, г1] и р[1, 3, г2], и результатов реализации метода (расчетных схем) два —

РС(п1) и РС(п2) соответственно. Результаты тестирования — зависимость нормы максимальной глобальной погрешности от величины шага интегрирования — представлены в табл. 5. Тестирование подтвердило теоретическое ожидание результатов эксперимента. При равных порядках погрешности полученный метод (см. табл. 2) требует меньшего числа этапов, чем метод Дормана — Принса (БР5(4)7). При равных же затратах РС(п1) и РС(п2) на порядок точнее метода Рунге — Кутты четвертого порядка.

Список литературы

[1] ОлЕМСКОй И.В. Численный метод интегрирования систем обыкновенных дифференциальных уравнений // Мат. методы анализа управляемых процессов. Л., 1986. С. 157-160.

[2] ОлЕМСКОй И.В. Экономичная расчетная схема четвертого порядка точности численного интегрирования систем специального вида // Процессы управления и устойчивость: Тр. XXX научн. конф. СПб.: НИИ химии СПбГУ, 1999. С. 134-143.

[3] ОлЕМСКОй И.В. Структурный подход в задаче конструирования явных одношаговых методов // Журн. вычисл. математики и мат. физики. 2003. Т. 43, № 7. С. 961-974.

[4] ОлЕМСКОй И.В. Четырехэтапный метод пятого порядка точности численного интегрирования систем специального вида // Журн. вычисл. математики и мат. физики. 2002. Т. 42, № 8. С. 1179-1190.

[5] ОлЕМСКОй И.В. Методы типа Рунге — Кутты интегрирования систем и дифференциальных уравнений второго порядка специального вида // Вычисл. технологии. 2004. Т. 9, № 2. С. 67-81.

[6] Современные численные методы решения обыкновенных дифференциальных уравнений / Под ред. Дж. Холла, Дж. Уатта. М., 1979.

[7] Хайрер Э., Нерсетт С., Ваннер Г. Решение обыкновенных дифференциальных уравнений. Нежесткие задачи. М., 1990.

[8] ОлЕМСКОй И.В. Алгоритм выделения структурных особенностей // Николай Ефимович Кирин: Сб. стат. под ред. В.В. Жука, В.Ф. Кузютина, 2003. С. 234-251.

Поступила в редакцию 18 июня 2004 г.