Текст подпрограммы и версий
qsgce_p.zip 
Тексты тестовых примеров
tqsgce_p.zip  tqsgce1_p.zip 

Подпрограмма:  QSGCE (модуль QSGCE_p)

Назначение

Быстрое вычисление узлов и весов квадратурной формулы Гаусса с расширенной (Extended) точностью.

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

Вычисляются узлы X ( k ) и веса C ( k ), k = 1, 2, ..., N квадратуры Гаусса на отрезке [- 1, 1] при любом N (2 ≤ N < ∞) за O ( N ) операций с расширенной (Extended) точностью.

Л.Г. Васильева, Я.М. Жилейкин, О быстром вычислении узлов и весов квадратуры Гаусса, ЖВМ и МФ, 2003.

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

procedure qsgce(N :Integer; N1 :Integer; var X :Array of Extended;
                var C :Array of Extended);

Параметры

N - число узлов квадратуры Гаусса (тип: целый);
N1 - N1 = N / 2, если N - четное,
N1 = ( N + 1) / 2, если N - нечетное, (тип: целый);
X - вектор длины N1 для вычисленных узлов квадратуры Гаусса (тип: вещественный);
C - вектор длины N1 для вычисленных весов квадратуры Гаусса (тип: вещественный).

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

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

  1. 

Требуется, чтобы N ≥ 2.

  2. 

Узлы квадратур Гаусса располагаются симметрично на отрезке [- 1, 1] относительно начала координат: X ( k ) = - X ( N - k + 1), а веса в симметричных точках равны C ( k ) = C ( N - k + 1), k = 1, 2, ..., N1. Поэтому подпрограмма QSGCD вычисляет N1 узлов и весов Гаусса только для отрезка [0, 1].

  3. 

Узлы и веса X ( k ) и C ( k ) вычисляются с точностью 10- 13.

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

 1.
Unit tqsgce_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc, UtRes_p, QSGCE_p;

function tqsgce: String;

implementation

function tqsgce: String;
var
N,N1,I :Integer;
X :Array [0..2] of Extended;
C :Array [0..2] of Extended;
begin
Result := '';  { результат функции }
{      ЕСЛИ N - ЧЕТНОЕ, ТО N1=N/2; ЕСЛИ N - НЕЧЕТНОЕ, ТО N1=(N+1)/2 }
N := 5;
N1 := 3;
QSGCE(5,3,X,C);
Result := Result + Format('%s',['  EPS=1.D-15' + #$0D#$0A]) + #$0D#$0A; 
Result := Result + Format('%s',['  N=']);
Result := Result + Format('%5d ',[N]) + #$0D#$0A;
Result := Result + Format('%s',
 ['  МАССИВ УЗЛОВ Х(N1)' + #$0D#$0A]) + #$0D#$0A; 
Result := Result + #$0D#$0A;
for I:=1 to N1 do
 begin
  Result := Result + Format('%20.16f ',[X[I-1]]) + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',
 ['  МАССИВ ВЕСОВ С(N1)' + #$0D#$0A]) + #$0D#$0A; 
Result := Result + #$0D#$0A;
for I:=1 to N1 do
 begin
  Result := Result + Format('%20.16f ',[C[I-1]]) + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
UtRes('tqsgce',Result);  { вывод результатов в файл tqsgce.res }
end;

end.

Результаты:

     X(1) = .906179845938664          C(1) = .236926885056189
     X(2) = .538469310105683          C(2) = .478628670499366
     X(3) = .000000000000000          C(3) = .568888888888889

  2. 
Unit tqsgce1_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc, UtRes_p, QSGCE_p;

function tqsgce1: String;

implementation

function tqsgce1: String;
var
N,N1,I :Integer;
X :Array [0..255] of Extended;
C :Array [0..255] of Extended;
begin
Result := '';  { результат функции }
{      ЕСЛИ N - ЧЕТНОЕ, ТО N1=N/2; ЕСЛИ N - НЕЧЕТНОЕ, ТО N1=(N+1)/2 }
N := 512;
N1 := 256;
QSGCE(N,N1,X,C);
Result := Result + Format('%s',['  EPS=1.D-15' + #$0D#$0A]) + #$0D#$0A; 
Result := Result + Format('%s',['  N=']);
Result := Result + Format('%5d ',[N]) + #$0D#$0A;
Result := Result + Format('%s',
 ['  МАССИВ УЗЛОВ Х(N1)' + #$0D#$0A]) + #$0D#$0A; 
Result := Result + #$0D#$0A;
for I:=1 to N1 do
 begin
  Result := Result + Format('%20.16f ',[X[I-1]]) + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',
 ['  МАССИВ ВЕСОВ С(N1)' + #$0D#$0A]) + #$0D#$0A; 
Result := Result + #$0D#$0A;
for I:=1 to N1 do
 begin
  Result := Result + Format('%20.16f ',[C[I-1]]) + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
UtRes('tqsgce1',Result);  { вывод результатов в файл tqsgce1.res }
exit;
end;

end.
 
Результаты:

     X(1) = .999988990984382          C(1) = .000028252637374
     X(2) = .999941994606846          C(2) = .000065765731660
     ......................................................................................................

     X(255) = .009194771386433        C(255) = .006129674838036
     X(256) = .003064962185159        C(256) = .006129905175406