Текст подпрограммы и версий
qtsxr_c.zip  qtsxd_c.zip  		Тексты тестовых примеров
tqtsxr_c.zip  tqtsxd_c.zip 

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

Назначение

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

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

Вычисляются значения неопрeдeлeнного двойного интеграла по прямоугольникам x1 ≤ x ≤ xi,  y1 ≤ y ≤ yj,   от табличной функции f (x, y), заданной на неравномерной сетке (xi, yj),   i = 1, 2, ..., N,   j = 1, 2, ..., M, по квадратурной формуле, точной для многочленов второй степени.

Н.С.Бахвалов. Численные методы, "Hаука", 1975.

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

    int qtsxr_c (real *rint, real *x1, real *x2, real *f,
            integer *n, integer *m)

   Параметры
rint - вещественный двумерный массив размера n на m, содержащий вычисленные значения неопределенного интеграла;
x1 - вещественный вектоp длины n, содержащий значения узлов неравномерной сетки по x;
x2 - вещественный вектоp длины m, содержащий значения узлов неравномерной сетки по y;
f - вещественный двумерный массив размера n на m, содержащий значения подинтегральной функции f (x, y);
n - заданное число узлов сетки по x (тип: целый);
m - заданное число узлов сетки по y (тип: целый). 
   Версии
qtsxd_c - вычисление с удвоенной точностью неопределенного двукратного интеграла от табличной функции, заданной на неравномерной сетке, по квадратурной формуле, точной для многочленов второй степени. 
   Вызываемые подпрограммы:  нет 
   Замечания по использованию
 

Требуется, чтобы  n ≥ 4  и  m ≥ 4.

В подпрограмме qtsxd_c параметры rint, x1, x2, f имеют тип double.

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

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

    /* Local variables */
    static float rint[2000] /* was [50][40] */,
                    f[2000] /* was [50][40] */;
    static int i__, j, m, n;
    static float x, h1, h2, x1[50], x2[40];
    extern int qtsxr_c(float *, float *, float *, float *, int *, int *);
    int i__1, i__2;

#define rint_ref(a_1,a_2) rint[(a_2)*50 + a_1 - 51]
#define f_ref(a_1,a_2) f[(a_2)*50 + a_1 - 51]

    h1 = .064114081632653058f;
    h2 = .08055358974358974f;
    n = 50;
    m = 40;
    x1[0] = 0.f;
    x2[0] = 0.f;
    x1[1] = h1 / 4.f;
    x2[1] = h2 / 4.f;
    i__1 = n;
    for (i__ = 3; i__ <= i__1; ++i__) {
/* l1: */
        x1[i__ - 1] = x1[i__ - 3] + h1 * 2.f;
    }
    i__1 = m;
    for (j = 3; j <= i__1; ++j) {
/* l2: */
        x2[j - 1] = x2[j - 3] + h2 * 2.f;
    }
    i__1 = n;
    for (i__ = 1; i__ <= i__1; ++i__) {
        i__2 = m;
        for (j = 1; j <= i__2; ++j) {
            x = x1[i__ - 1] + x2[j - 1];
/* l3: */
            f_ref(i__, j) = (float)sin(x);
        }
    }
    qtsxr_c(rint, x1, x2, f, &n, &m);

        printf("\n  %16.7f %16.7f %16.7f %16.7f \n",
            rint_ref(1, 1), rint_ref(2, 2), rint_ref(n-1, m), rint_ref(n, m));
    return 0;
} /* main */


Результаты:

       rint(1, 1)        =  0.00000
       rint(2, 2)        =  0.00001
       rint(n-1, m)   =  0.24858
       rint(n, m)      =  0.21648