Текст подпрограммы и версий
de27r_p.zip , de27e_p.zip 		Тексты тестовых примеров
tde27r_p.zip , tde27e_p.zip

Подпрограмма:  DE27R (модуль DE27R_p)

Назначение

Вычисление решения задачи Коши для системы обыкновенных дифференциальных уравнений первого порядка в конце интервала интегрирования методом Адамса с постоянным шагом.

Математическое описание

Решается задача Коши для системы M обыкновенных дифференциальных уравнений первого порядка 

          Y ' = F (X, Y) ,

          Y = ( y1, ..., yM ) ,
          F = ( f1 (X, y1,..., yM),..., fM (X, y1,..., yM) )
 с начальными условиями, заданными в точке XN:
          Y(XN) = YN ,     YN = ( y10,..., yM0 ) ,

многошаговым предсказывающе - исправляющим методом пятого порядка точности. Суть метода состоит в том, что сначала по по явной формуле Адамса строятся предсказанные значения приближенного решения, которые затем исправляются по неявной формуле Адамса.

Решение вычисляется в одной точке XK, которая является концом интервала интегрирования. Интегрирование ведется с постоянным шагом, значение которого задается пользователем при обращении к подпрограмме.

Березин И.С., Жидков Н.П. Методы вычислений, т.2. Физматгиз, М., 1960.

Беленький В.З. Стандартная программа для интегрирования системы дифференциальных уравнений методом Адамса. Сб. "Вычислительные методы и программирование", вып. 3. Изд-во МГУ, 1965

Использование 

procedure DE27R(F :Proc_F70; var M :Integer; XN :Real;
                var YN :Array of Real; XK :Real; var H :Real;
                var Y :Array of Real; var DELTY :Array of Real;
                var DF :Array of Real; var RF :Array of Real;
                var YP :Array of Real; var RY :Array of Real;
                var RFN :Array of Real; var IERR :Integer);


   Параметры
F - имя подпрограммы вычисления значений правой части дифференциального уравнения. Первый оператоp подпрограммы должен иметь вид:
procedure F (X :Real; var Y :Array of Real; var DY :Array of Real; M :Integer);
Здесь: X, Y - значения независимой и зависимой переменных, соответственно. Вычисленное значение правой части должно быть помещено в DY. B случае системы уравнений, т.е. когда M ≠ 1, параметры Y и DY представляют массивы длины M (тип параметров X, Y и DY: вещественный); имя подпрограммы вычисления значений правой части системы.
M - количество уравнений в системе (тип: целый);
XN, YN - начальные значения аргумента и решения; в случае системы уравнений (т.е. M ≠ 1) YN представляет массив длины M (тип: вещественный);
XK - значение аргумента, при котоpом требуется вычислить решение задачи Коши (конец интервала интегрирования). XK может быть больше, меньше или pавно XN (тип: вещественный);
H - вещественная переменная, содержащая значение шага интегрирования; может задаваться с учетом направления интегрирования, т.е. положительным, если XK > XN, отрицательным, если XK < XN, или без такого учета в виде абсолютной величины;
Y - искомое решение задачи Коши, вычисленное подпрограммой для значения аргумента XK. Для системы уравнений (когда M ≠ 1) У задается массивом длиной M. B случае совпадения значений параметров XN и XK значение Y полагается равным начальному значению YN (тип: вещественный);
         DELTY -
         RF, YP  
         RFN   		одномерные вещественные рабочие массивы длины M;
DF, RY - двумерные вещественные рабочие массивы размеpа M*5;
IERR - целая переменная, значение которой в pезультате работы подпрограммы полагается равным 65, если при заданном значении H шага интегрирования не может быть достигнута сходимость численного решения к точному значению решения; в этом случае интегрирование системы можно повторить обращением к подпрограмме с новым (меньшим) значением H. 
   Версии
DE27E - вычисление решения задачи Коши для системы обыкновенных дифференциальных уравнений первого порядка в конце интервала интегрирования методом Адамса с постоянным шагом, при этом все вычисления над действительными числами ведутся с удвоенным числом значащих цифр. B этом случае параметры XN, YN, XK, H, Y, DELTY, DF, RF, YP, RY, RFN и параметры X, Y и DY в подпрограмме F должны иметь тип Extended. 
   Вызываемые подпрограммы
DE26R - выполнение одного шага численного интегрирования системы обыкновенных дифференциальных уpавнений первого порядка методом Адамса без контроля точности; вызывается при работе подпрограммы DE27R.
DE26E - выполнение одного шага численного интегрирования системы обыкновенных дифференциальных уpавнений первого порядка методом Адамса без контроля точности, при этом все вычисления над действительными числами ведутся с удвоенным числом значащих цифр; вызывается при работе подпрограммы DE27E.
UTDE12 - подпрограмма выдачи диагностических сообщений при работе подпрограммы DE27R.
UTDE13 - подпрограмма выдачи диагностических сообщений при работе подпрограммы DE27E. 
   Замечания по использованию
 

При работе подпрограммы и ее версии значения параметров M, XN, YN, XK, H сохраняются.

Если после работы подпрограммы нет необходимости иметь начальное значение решения YN, то параметры YN и Y при обращении к ней можно совместить.

Пример использования 

           y1'  =   y2
           y2'  =  -y1

          y1 (3/4 π)  =   √2 /2 ,     y2 (3/4 π)  =  -√2 /2

Приводятся подпрограмма вычисления значений правой части системы и фрагмент вызывающей программы, выполняющей дважды интегрирование системы, сначала на отрезке [3/4 π, π], затем на том же отрезке справа налево, а также результаты счета. 

Unit TDE27R_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc, UtRes_p, FDE27R_p, DE27R_p;

function TDE27R: String;

implementation

function TDE27R: String;
var
M,I,IERR :Integer;
XN,XK,H :Real;
YN :Array [0..1] of Real;
Y :Array [0..1] of Real;
DELTY :Array [0..1] of Real;
DF :Array [0..9] of Real;
RF :Array [0..1] of Real;
YP :Array [0..1] of Real;
RY :Array [0..9] of Real;
RFN :Array [0..1] of Real;
label
_30;
begin
Result := '';  { результат функции }
M := 2;
XN := 0.75*3.14159265358979323846264e0;
YN[0] := Sqrt(2.0)/2.0;
YN[1] := -YN[0];
ХК := 3.14159265358979323846264e0;
for I:=1 to 2 do
 begin
  H := 0.01;
  DE27R(FDE27R,M,XN,YN,XK,H,Y,DELTY,DF,RF,YP,RY,RFN,IERR);
  Result := Result + Format('%s',['       XK']);
  Result := Result + Format('%s',['               Y(1,1)']);
  Result := Result + Format('%s',['               Y(1,2)' + #$0D#$0A]);
  Result := Result + Format(' %20.16f  %20.16f %20.16f ',
 [XK,Y[0],Y[1]]) + #$0D#$0A;
  XN := XK;
  ХК := 0.75*XN;
  YN[0] := 0.0;
  YN[1] := -1.0;
_30:
 end;
UtRes('TDE27R',Result);  { вывод результатов в файл TDE27R.res }
exit;
end;

end.

Unit fde27r_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc;

procedure fde27r(X :Real; var Y :Array of Real;
                var DY :Array of Real; M :Integer);

implementation

procedure fde27r(X :Real; var Y :Array of Real;
                var DY :Array of Real; M :Integer);
begin
DY[0] := Y[1];
DY[1] := -Y[0];
end;

end.

Результаты:

    после первого обращения к подпрограмме -
          XK                               Y(1)                         Y(2)
          3.141592654 + 00      -2.773237299-10      -1.000000000 + 00
 
    после второго обращения к подпрограмме -
          XK                               Y(1)                         Y(2)
          2.356194490 + 00       7.071067811-01      -7.071067815-01