Текст подпрограммы и версий ( Фортран )
qsf9r.zip  qsf9d.zip  		Тексты тестовых примеров ( Фортран )
tqsf9r.zip  tqsf9d.zip 
Текст подпрограммы и версий ( Си )
qsf9r_c.zip  qsf9d_c.zip  		Тексты тестовых примеров ( Си )
tqsf9r_c.zip  tqsf9d_c.zip 
Текст подпрограммы и версий ( Паскаль )
qsf9r_p.zip  qsf9e_p.zip  		Тексты тестовых примеров ( Паскаль )
tqsf9r_p.zip  tqsf9e_p.zip 

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

Назначение

Вычисление определенных интегралов 

                  B
(1)             ∫  F (x) sin( ρ x + φ ) dx ,
                A
                  B
(2)             ∫  F (x) cos( ρ x + φ ) dx
                A

от векторной функции по формулам интерполяционного типа пятой степени точности.

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

Подпрограмма QSF9R вычисляет интегралы (1) и (2), или один из них от векторной функции F (x) = ( f1 (x), ..., fN (x) ) с покомпонентной погрешностью 

                     B
      E ( 1 + |  ∫  fi (x) sin( ρ x + φ) dx | )  или
                   A     
                     B                               
      E ( 1 + |  ∫  fi (x) cos( ρ x + φ) dx | )      
                   A

При ρ h ≥ 0.33 вычисление ведется по формулам интерполяционного типа пятой степени точности, а при ρ h < 0.33 по обобщенной квадратурной формуле Hьютона - Котеса шестой степени точности с шагом h = (B - A)/M, где M - число частичных отрезков разбиения.

Я.М.Жилейкин, М.В.Соколовский. Набор стандартных программ для вычисления интегралов от векторных функций. Сб. "Методы и алгоритмы в численном анализе", Изд-во МГУ, 1981.

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

    SUBROUTINE  QSF9R (RINT, A, B, F, RO, FI, N, E, K, L1, L2,
                                              IND, R, IERR) 
   Параметры
RINT - вещественный вектоp, содержащий вычисленные значения интегралов; его длина pавна 2N, если вычисляются оба интеграла (1) и (2), или N, если вычисляется один из этих интегралов;
A, B - заданные нижний и верхний пределы интегрирования (тип: вещественный);
F - имя подпрограммы, вычисляющей подинтегральные функции;
RO, FI - заданные значения параметров  ρ и  φ (тип: вещественный);
N - размерность вектора - функции F (x) (тип: целый);
E - заданная мера погрешности вычисления интеграла (тип: вещественный);
K - целая переменная, задающая начальное число частичных отрезков разбиения;
L1, L2 - задают режим работы подпрограммы (тип: целый); при этом, если:
 

L1 = 1 и L2 = 0, то вычисляется интеграл (1);

L1 = 0 и L2 = 1, то вычисляется интеграл (2);

L1 = 1 и L2 = 1, то вычисляются интегралы (1) и (2);
IND - целый вектоp длины N, каждая компонента котоpого pавна числу частичных отрезков разбиения, при котоpом достигается заданная точность по этой компоненте;
R - вещественный рабочий вектоp длины 3N, если вычисляются оба интеграла, и 2N, если вычисляется один из них;
IERR - целая переменная, служащая для диагностических сообщений:
IERR = 65 - если заданная точность не может быть достигнута при максимально возможном числе (1048576) частичных отрезков разбиения. 
   Версии
QSF9D - вычисление с удвоенной точностью определенных интегралов (1) и (2) от векторной функции по формулам интерполяционного типа пятой степени точности. 
   Вызываемые подпрограммы
UTQS10 - подпрограмма выдачи диагностических сообщений при работе подпрограммы QSF9R.
UTQS12 - подпрограмма выдачи диагностических сообщений при работе подпрограммы QSF9D. 
   Замечания по использованию
  1. 

Первый оператор подпрограммы F должен иметь вид:
      SUBROUTINE F(X, Y, I)
Здесь:
X - аргумент функции F (тип: вещественный);
Y - вещественный вектоp длины N вычисленных значений функции F;
I - целый вектоp длины N; если его i - ая компонента отлична от нуля, то интеграл от этой компоненты функции F не вычисляется (тип: целый).

  2. 

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

  3. 

Вычисление интеграла от отдельной компоненты прекращается как только на этой компоненте достигается заданная точность.

  4. 

Если для каких - либо компонент векторной функции заданная погрешность интегрирования не может быть достигнута, то соответствующие компоненты вектоpа IND полагаются равными нулю.

  5. 

Погрешность определяется по модулю разности двух просчетов по M и 2M частичным отрезкам разбиения. Если заданная точность не достигнута, то значение M удваивается.

  6.  В подпрограмме QSF9D параметры RINT, A, B, F, RO, FI, E, R имеют тип DOUBLE PRECISION.

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

      SUBROUTINE  F(X, Y, I)
      DIMENSION  Y(1)
      Y = X
      DO 1 J = 1, 4
   1 Y(J + 1) = Y(J)*X
      RETURN
      END

      DIMENSION  RINT(5, 2), R(15), IND(5)
      EXTERNAL  FUN
      PI0 = 1.5707963268
      CALL  QSF9R (RINT, -PI0, PI0, FUN, 1., 0., 5, 1.E-3, 2, 1, 1, IND, R, 
     *                          IERR)

Результаты:

              RINT (I, 1)                           RINT (I, 2)                 IND (I)

      1.9999999618 E+00             -8.661 E-12                            2
     -2.676 E-11                             9.3480217805 E-01              2
      2.8044065499 E+00             -6.428 E-12                            2
     -5.268 E-11                             9.5850994403 E-01              2
      4.7925452374 E+00             -2.098 E-11                            4

      IERR  =  0

      I  =   (1, 2, 3, 4, 5)