Текст подпрограммы и версий
qse1r_c.zip  		Тексты тестовых примеров
tqse1r_c.zip 

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

Назначение

Вычисление на равномерной сетке определенного интеграла от быстроосциллирующих функций по квадратурной формуле типа Эйлера.

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

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

           a
           ∫ f(x) eiwx dx  =  J1M + i J2M ,     i = √-1
          b

на равномерной сетке с помощью квадратурной формулы интерполяционного типа, точной для многочленов четвертой степени с весом  e iwx .
M - число узлов интегрирования.

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

    int qse1r_c (real *rint, real *a, real *b, R_fp f, real *r,
            real *e, integer *n, integer *ierr)

   Параметры
rint - вещественный вектоp длины 2 для вычисленных значений J1 (M) и J2 (M);
a, b - заданные нижний и верхний пределы интегрирования (тип: вещественный);
f - имя вещественной подпрограммы - функции, вычисляющей подинтегральную функцию f (x);
r - заданное значение параметра  w (тип: вещественный);
e - заданная абсолютная погрешность вычисления интеграла (тип: вещественный);
n - целая переменная, задающая начальное число узлов интегрирования;
ierr - целая переменная для диагностических сообщений:
ierr = 65 - когда заданная точность не может быть достигнута при максимально возможном числе узлов (m = 1048576). 
   Версии:  нет 
   Вызываемые подпрограммы
utqs11_c - подпрограмма выдачи диагностических сообщений при работе подпрограммы qse1r_c. 
   Замечания по использованию
  1. 

При e ≤ 0 происходит только один просчет при заданном значении n.

  2. 

Максимальное значение n, котоpое можно задавать, не должно превосходить 524288.

  3. 

Если интеграл вычислен с точностью e, то ierr полагается равным нулю.

  4. 

Если интеграл не может быть вычислен с точностью e, то выдается диагностическое сообщение "заданная точность не может быть достигнута".

  5.  Если ierr = 65, то n полагается равным  -1.

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

int main(void)
{
    /* Local variables */
    static int ierr;
    static float rint[2];
    extern int qse1r_c(float *, float *, float *, R_fp, float *,
                       float *, int *, int *);
    static float a, b, e;
    static int n;
    static float r__;
    extern float f1_c(), f2_c();

    a = -1.f;
    b = 1.f;
    r__ = -112.f;
    e = 1e-7f;
    n = 10;
    qse1r_c(rint, &a, &b, (R_fp)f1_c, &r__, &e, &n, &ierr);

    printf("\n %16.7e %16.7e \n",rint[0],rint[1]);
    printf("\n %5i %5i \n",ierr,n);
    qse1r_c(rint, &c_b6, &c_b7, (R_fp)f2_c, &c_b8, &c_b9, &c__10, &ierr);

    printf("\n %16.7e %16.7e \n",rint[0],rint[1]);
    return 0;
} /* main */

float f1_c(float *x)
{
    /* System generated locals */
    float ret_val, r__1;

/* Computing 4th power */
    r__1 = *x, r__1 *= r__1;
    ret_val = r__1 * r__1;
    return ret_val;
} /* f1_c */

float f2_c(float *x)
{
    /* System generated locals */
    float ret_val;

    /* Builtin functions */
    double sin(double);

    ret_val = (float)sin((float)(*x / 112.f));
    return ret_val;
} /* f2_c */


Результаты:

    1.   rint(1)  =  -0.0155869200 ;
          rint(2)  =   0 ;
               ierr  =   0
      
    2.   rint(1)  =   0. ;
          rint(2)  =   0.0000739651 ;
               ierr  =   0