Подпрограмма:  de28r_c

Назначение

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

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

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

          Y ' = F (X, Y) ,

          Y = ( y1, ..., yM ) ,
          F = ( f1 (X, y1,..., yM), ... , fM (X, y1,..., yM) ) 

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

По заданному значению YX решения в узле X_n вычисляется значение этого решения в узле X_n + H. Все компоненты решения вычисляются с контролем точности по меpе погрешности, который заключается в следующем. Если некоторая компонента приближенного решения по абсолютной величине не меньше некоторой заданной константы P (называемой границей перехода), то контроль точности для этой компоненты ведется по относительной погрешности, иначе - по абсолютной. Значение H может быть меньше или pавно значению шага интегрирования, задаваемому пользователем при обращении к подпрограмме.

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

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

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

    int de28r_c (S_fp f, integer *m, integer *jstart, real *hmin,
                 real *hmax, real *eps, real *p, real *yx, real *x, real *h,
                 logical *bul, real *delty, real *df, real *rf, real *yp, real *ry,
                 real *rfn, integer *ierr)

   Параметры 

   Версии 

   Вызываемые подпрограммы 

   Замечания по использованию  

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

           y1'  =   y2
           y2'  =  -y1

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

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

int main(void)
{
    /* Builtin functions */
    double sqrt(double), sin(double), cos(double);

    /* Local variables */
    extern int de28r_c(U_fp, int *, int *, float *, float *, float *,
                       float *, float *, float *, float *, logical *,
                       float *, float *, float *, float *, float *, float *,
                       int *);
    static float hmin, hmax;
    static int ierr;
    extern int f_c();
    static float h__;
    static int m;
    static float p, x, y[2], delty[2], y1, y2, df[10];
    static int ih;
    static float rf[2], yp[2], ry[10];
    static int jstart;
    static logical bul;
    static float rfn[2], eps;

    m = 2;
    x = 2.3561944901925003f;
    h__ = .01f;
    y[0] = (float)sqrt(2.f) / 2.f;
    y[1] = -y[0];
    hmin = 1e-9f;
    jstart = 0;
    eps = 1e-5f;
    p = 100.f;
    bul = FALSE_;
    ih = 0;
    hmax = 1.f;
l6:
    de28r_c((U_fp)f_c, &m, &jstart, &hmin, &hmax, &eps, &p, y, &x, &h__, &bul,
            delty, df, rf, yp, ry, rfn, &ierr);
    ++ih;
    if (ih == 1) { goto l6; }
    y1 = (float)sin(x);
    y2 = (float)cos(x);

    printf("\n %16.7e %16.7e %16.7e \n", x, y[0], y[1]);
    printf("\n %16.7e %16.7e %16.7e \n\n", h__, y1, y2);
    if (ih == 4) {
        goto l20;
    }
    jstart = -1;
    if (ih == 3) {
        goto l7;
    }
    h__ = .005f;
    goto l6;
l7:
    h__ = .02f;
    goto l6;
l20:
    return 0;
} /* main */

int f_c(float *x, float *y, float *dy, int *m)
{
    /* Parameter adjustments */
    --dy;
    --y;

    /* Function Body */
    dy[1] = y[2];
    dy[2] = -y[1];
    return 0;
} /* f_c */


Результаты:

    после второго обращения к подпрограмме -
          x                                y(1)                         y(2)
          2.376194490+00        6.928241717-01     -7.211065574-01
          h__                              y1                            y2
          1.000000000-02         6.928241717-01     -7.211065574-01
 
    после третьего обращения к подпрограмме -
          x                                y(1)                         y(2)
          2.371194490+00        6.964210292-01     -7.176334371-01
          h__                              y1                            y2
          5.000000000-03         6.964210292-01     -7.176334371-01

    после четвертого обращения к подпрограмме -
          x                                y(1)                         y(2)
          2.386194490+00        6.855785854-01     -7.279986286-01
          h__                              y1                            y2
          2.000000000-02         6.855785854-01     -7.279986286-01