Документ взят из кэша поисковой машины. Адрес
оригинального документа
: http://num-anal.srcc.msu.ru/lib_na/cat/q_htm_p/qsgce_p.htm
Дата изменения: Tue Nov 17 12:59:15 2015 Дата индексирования: Sun Apr 10 02:06:05 2016 Кодировка: Windows-1251 |
Текст подпрограммы и версий qsgce_p.zip |
Тексты тестовых примеров tqsgce_p.zip tqsgce1_p.zip |
Быстрое вычисление узлов и весов квадратурной формулы Гаусса с расширенной (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