Текст подпрограммы и версий
idb9r_c.zip , idb9d_c.zip
Тексты тестовых примеров
tidb9r_c.zip , tidb9d_c.zip

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

Назначение

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

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

Пусть заданы узлы неравномерной прямоугольной сетки:

   X1A = { x1(1), x2(1), ..., xM(1) } ;   X2A  =  { x1(2), x2(2), ..., xN(2) }

и значения функции  f ( x(1), x(2) ) в узлах этой сетки:

  YA(I, J)  =  f ( X1A(I), X2A(J) ) ,   I = 1, 2, ..., M ;   J = 1, 2, ..., N.

Предполагается, что  x1(1) < x2(1) <...< xM(1)  и  x1(2) < x2(2) <...< xN(2). 

Подпрограмма idb9r_c вычисляет вторые частные производные

2 f (x(1), x(2)) / ∂x(2) 2 

в узлах заданной сетки методом одномерных натуральных кубических сплайнов по направлению x (2).

Н.С.Бахвалов. Численные методы. Изд - во "Наука", 1973.

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

    int idb9r_c (real *x1a, real *x2a, real *ya, integer *m,
            integer *n, real *y2a, real *yt, real *y2t, real *rab)

Параметры

x1a -
x2a  
вещественные векторы длины m и n, содержащие узлы  { x1 (1), x2 (1), ...,xm (1)}  и  { x1 (2), x2 (2), ...,xn (2)} соответственно;
ya - вещественный двумерный массив размеров m на n, содержащий значения функции  f (x (1), x (2)) двух переменных в узлах заданной сетки;
m, n - длины векторов x1a и x2a соответственно (тип: целый);
y2a - вещественный двумерный массив размеров m на n, содержащий вычисленные значения второй частной производной функции  f по второй переменной в узлах заданной сетки;
yt, y2t -
rab  
вещественные векторы длины n, используемые в подпрограмме в качестве рабочих.

Версии

idb9d_c - вычисление второй частной производной функции двух переменных, заданной на прямоугольной неравномерной сетке, по второй переменной методом одномерных натуральных кубических сплайнов в режиме удвоенной точности; все формальные параметры, за исключением m и n, должны иметь тип double.

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

ids8r_c -
ids8d_c  
вычисление вторых производных таблично - заданной функции в узлах неравномерной сетки методом кубических или натуральных кубических сплайнов при заданных первых производных функции в концевых узлах сетки в режимах одинарной и удвоенной точности соответственно; используются в подпрограммах idb9r_c и idb9d_c.

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

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

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

    /* Local variables */
    extern int idb9r_c(float *, float *, float *, int *, int *, float *,
                       float *, float *, float *);
    static int i__, j, m, n;
    static float r__, ya[20] /* was [5][4] */, yt[4], x1a[5], x2a[4],
                     y2a[20] /* was [5][4] */, y2t[4], rab[4],
                    y2at[20] /* was [5][4] */;
    int i__1, i__2;

#define ya_ref(a_1,a_2) ya[(a_2)*5 + a_1 - 6]
#define y2at_ref(a_1,a_2) y2at[(a_2)*5 + a_1 - 6]

    m = 5;
    n = 4;
    r__ = 0.f;
    i__1 = m;
    for (i__ = 1; i__ <= i__1; ++i__) {
        x1a[i__ - 1] = r__;
/* l1: */
        r__ += .1f;
    }
    r__ = 0.f;
    i__1 = n;
    for (i__ = 1; i__ <= i__1; ++i__) {
        x2a[i__ - 1] = r__;
/* l2: */
        r__ += .1f;
    }
    i__1 = m;
    for (i__ = 1; i__ <= i__1; ++i__) {
        i__2 = n;
        for (j = 1; j <= i__2; ++j) {
            ya_ref(i__, j) = (float)sin(x1a[i__ - 1]) * (float)sin(x2a[j - 1]);
/* l3: */
            y2at_ref(i__, j) = -ya_ref(i__, j);
        }
    }
    idb9r_c(x1a, x2a, ya, &m, &n, y2a, yt, y2t, rab);

    for (i__ = 0; i__ <= 15; i__+= 5) {
        printf("\n %15.6e%15.6e%15.6e%15.6e%15.6e \n",
                    y2a[i__], y2a[i__+1], y2a[i__+2], y2a[i__+3], y2a[i__+4]);
    }
    for (i__ = 0; i__ <= 15; i__+= 5) {
        printf("\n %15.6e%15.6e%15.6e%15.6e%15.6e \n",
               y2at[i__], y2at[i__+1], y2at[i__+2], y2at[i__+3], y2at[i__+4]);
    }
    return 0;
} /* main */


Результаты:      

                   |  0.0      0.0                    0.0                    0.0
                   |  0.0     -0.80066e-02   -0.27724e-01    0.0
       y2a =   |   0.0     -0.15933e-01   -0.55172e-01    0.0
                   |  0.0     -0.23700e-01   -0.82068e-01    0.0
                   |  0.0     -0.31230e-01   -0.10814            0.0