циальных уравнений методом Рунге -Кутта с автоматическим выбором шага. Параметры:

X - .начальное значение иезависимой .переменной х\ t/[k] - .начальные значения искомых (функций Ук{х); п-порядок системы уравнений;

f{x, у, п, z)-процедура, представляющая .интегрируемую систему (совокупность функций fh), т. е. процедура, осуществляющая счет правых частей уравнений и выдающая в качестве результата .массив z значений производных yh в любой заданной «точке» (х, у и т......Уп), т. е. z{k]=fh{x, г/[1]...., у[п]);

eps, eta - значения, определяющие точность численного интегрирования;

xfin - конец интервала интегрирования;

j in -выходной па1раметр, представляющий решение в точке x=xfin;

prim - булевский параметр, значение которого .при первом обращении в процедуре runige должно быть true. Если же нужно получить значения функций у в нескольких промежуточных точках Хо, Xi,..., Хп, то процедуру необходимо вызывать повторно п раз с х==хи, xfin-Xh+i для ft=0,1,...", п-1. В последнем случае повторные обращения к range для экономии машинного времени можно осуществлять с prim=ialse;

сотр - нелокализованный (глобальный) идентификатор вещественной працедуры-фун.кции сотр {а, Ь, с), значением которой является абсолютное значение разности мантисс avib после того как порядки этих величин выравнены до наибольшего порядка параметров а, b п с;

rklstep{x, у, h, xh, yh) -процедура, интегрирующая систему уравнений на одном шаге h методом Рунге - Кутта с начальными значениями х, у[к]. Выходными параметрами являются xh=x+h и yh[k], причем последний вектор есть решение системы в точке xh. Параметры п, /, Z входят в rklstep как нелокализованные объекты;

SS, hs - глобальные параметры, значения которых не должны изменяться между повторными обращениями к процедуре range. При первом обращении к range значения этих параметров могут быть не заданы. Переменная hs должна иметь тип real, а ss - тип integer.*

* Эта процедура содержит некоторые новые идеи, связанные с идеями Гилла [lOi] и Фрёберга Iilli). Однако должно быть ясно,

procedure runge(x,y,n,f,eps,eta,prim,xfin) result: (yfin); value x,y; real x,eps,eta,xfin; integer n; Boolean prim; array y,yfin; procedure f; begin real xl,x2,x3,h; integer k,j; Boolean out; -array z,yl,y2,y3[l:nl; procedure rklstep(x,y,h)result:(xh,yh); value x,h; real x,h,xh; array y,yh; begin integer j,k; array wll:nla[l:5I; af 1]: = aI2]: = a[5]:=0.5 Xh; a[3: = a[4]: = h; xh: = x;

for k:= 1 step 1 until n do yhfk]:=wlk]:=ylk]; for j: = 1,2,3,4 do begin f(xih,w,n,z); xh:=x+a{jl;

for k: = l step 1 until n do begin yblkj:=уЩк] + aO +1]Xzlkj/S;

wfk]:=yM+afj]XzIk] end к end j end rklstep; start: if prim then

)egin h:=xfin-x; ss:=0 end else h:=hs; out:=false; labl: if (x+2.01Xh-xfin>0) = (h>0) then begin hs: = h; out:=true;

h:=(xfin-x)/2 end;

rklstep(x,y,2Xh,xl,yl); lab2: rklstep(x,y,h,x2,y2); - • "

rklstep (x2,y2,h,x3,y3); for k: = 1 step 1 until n do if comp(ylIk],y3lk],eta)>eps then go to labS; x:=x3;

if out then go to fin; for k:=l step 1 untiLn do ylk]:=y3[k]; if ss=5 then begin ss:=0; h:=2Xh end; ss: = ss + l; go to labl; labS: h:=0.5xli; out:=false; xl:=x2;

for k:=l step 1 until n do yl[k]:=y2[kl; go to lab2;

что в отношении времени вычислений и ошибок округления она, возможно, не будет оптимальной. {Прим. П. Наура к алгоритму 9

fin: for k: = l step 1 until n do у11п1к]:=уЗ[к]

end runge;

Свидетельство к алгоритму 96

Процедура runge алгоритма 96 отличается от процедуры runge алгоритма 9а только тем, что в заголовке локализованной в ней подпроцедуры rklstep добавлен список значений value х, h\. Это требуется сокращенным АЛГОЛом-60 [19]. В пояснительном тексте к алгоритму 9а были сделаны некоторые уточнения. В частности, в соответствии с «Подтверждением» Г. Тачера (см. ниже) был изменен тип переменной ss.

Алгоритм 96 был транслирован на машине М-2120 в системе ТА-1М, и с его помощью были успешно решены следующие задачи.

1. Интегрирование уравнения у=-2ху для входных данных х-1; j/[l]=0.5; п=1; prim=true; xfin=2 (первый из примеров в «Свидетельстве к алгоритму 218а»). Точное решение у=0.2.

При eps=eta = 10- был получен результат г ш[1] = = 0.20000139 с абсолютной погрешностью 0.139x10-5. Алгоритм 218а при этих же условиях дает результат 0.199999623 с абсолютной погрешностью 0.377x10-е.

При eps=eta=lO- было получено точное значение j/fin[1] =0.20000000 с выходными значениями глобальных величин ss = l и As = 0.0625.

Процедура / задавалась описанием

procedure f (x,y,n,z);

value x,n; real x; integer n; array y,z; zl[l]:=-2XxXyIl]t2;

2. Интегрирование системы уравнений у1=Уг, у2=-Уи имеющей решение yi=sinx и y2=cosx (первый пример в нижеследующем «Подтверждении» Г. Тачера), для входных данных

х=0, у[1]=0, у12]=:1, п = 2, prim=true, xfin=7.

Для eps=eta-l(y- было получено решение yi= =0.65697307, /2=0.75391108, ss = 5 и /is=0.109375. Это решение отличается от результатов Г. Тачера на Ау1 =0.13x10-4 и Ау2=0.84 X. 10-.

Для eps = eG= 10-8 было получено 1=0.65698657, «/2=0.75390227, ss = l, As=0.0546875. Отличия от результатов Г. Тачера Ai/i==0.55X Ш"» и А2=0.43Х10-е.