Текст подпрограммы и версий aig3r_p.zip aig3e_p.zip aig3c_p.zip |
Тексты тестовых примеров taig3r_p.zip taig3e_p.zip taig3c_p.zip |
Обращение вещественной матрицы общего вида методом Гаусса с выбором ведущего элемента по столбцу и оценка числа обусловленности матрицы.
Для заданной вещественной квадратной матрицы А порядка N выполняется треугольная факторизация L- 1А = U, где L - нижняя треугольная матрица, U - верхняя треугольная матрица, вычисляется величина, обратная числу обусловленности матрицы А
RCOND = 1/(|| A ||1*|| A-1 ||1) , где || A ||1 = maxj = 1, ..., N{ | a1 j | + | a2 j | + ...+ | aN j | }.
Далее вычисляется матрица U-1, и затем вычисляется А-1 = U-1L-1.
Дж. Форсайт, М. Малькольм, К. Моулер. Машинные методы математических вычислений. М.: Мир, 1980.
procedure AIG3R(var A :Array of Real; var M :Integer; var N :Integer; var NLEAD :Array of Integer; var RCOND :Real; var Z :Array of Real; var IERR :Integer);
Параметры
A - | вещественный двумерный массив размера М на N, в котором задается исходная матрица; на выходе на месте исходной матрицы находится вычисленная обращенная матрица; |
M - | первая размерность массива А в вызывающей программе (тип: целый); |
N - | порядок обращаемой матрицы (тип: целый); |
NLEAD - | целый вектор длины N, содержащий на выходе информацию о выполненных в ходе факторизации перестановках (см. замечания по использованию); |
RCOND - | вещественная переменная, содержащая на выходе вычисленное значение величины, обратной числу обусловленности матрицы; |
Z - | вещественный рабочий вектор длины N; |
IERR - | целая переменная, содержащая на выходе информацию о прохождении счета; при этом: |
IЕRR=65 - | если М ≤ 0 или N ≤ 0; |
IЕRR=66 - | если в процессе счета возникло переполнение (это говорит о том, что либо ||А||1, либо некоторые элементы матрицы U, либо некоторые элементы обращенной матрицы превосходят по абсолютной величине максимальное представимое на данной машине число); |
IЕRR=-К - | если в результате факторизации в К-й строке матрицы U диагональный элемент равен нулю (это свидетельствует о вырожденности матрицы А). Если таких строк у матрицы U несколько, то значение K полагается равным номеру последней из них. |
Версии
AIG3E - | обращение методом Гаусса с выбором ведущего элемента по столбцу и оценка числа обусловленности вещественной матрицы А, заданной с расширенной (Extended) точностью. |
AIG3C - | обращение методом Гаусса с выбором ведущего элемента по столбцу и оценка числа обусловленности комплексной матрицы А. |
Вызываемые подпрограммы
AFG4R - | подпрограмма треугольной факторизации и оценки числа обусловленности матрицы А; |
UTAFSI - | подпрограмма выдачи диагностических сообщений. |
Замечания по использованию
1. |
В подпрограмме АIG3E массивы А и Z и переменная RСОNE имеют тип External, для треугольного разложения и оценки числа обусловленности матрицы А вызывается подпрограмма АFG4E. | |
2. |
В подпрограмме АIG3С массивы А и Z имеют тип Сomplex, для треугольного разложения матрицы А вызывается подпрограмма АFG4С. | |
3. |
На выходе К - й элемент вектора NLЕАE равен номеру строки, переставленной на К - м шаге факторизации с К - й строкой матрицы А. Так как факторизация Гаусса требует N - 1 шагов, то NLЕАD (N) = N. | |
4. |
Если вырабатывается значение переменной IЕRR, отличное от нуля, то выдается соответствующее диагностическое сообщение и происходит выход из подпрограммы. |
Unit TAIG3R_p; interface uses SysUtils, Math, { Delphi } LStruct, Lfunc, UtRes_p, AIG3R_p; function TAIG3R: String; implementation function TAIG3R: String; var J,I,JI,JJ,II,IERR,M,N,K :Integer; RCONE :Real; B :Array [0..15] of Real; C :Array [0..15] of Real; Z :Array [0..3] of Real; NLEAE :Array [0..3] of Integer; const A :АRRАУ [0..15] ОF RЕАL = ( 7.9,8.5,4.3,3.2,5.6,-4.8,4.2,-1.4,5.7,0.8,-3.2, -8.9,-7.2,3.5,9.3,3.3 ); labeL _1,_2,_3,_4,_5,_6,_7,_8; begin Result := ''; M := 4; N := M; for I:=1 to M do begin for J:=1 to N do begin B[(I-1)+(J-1)*4] := A[(I-1)+(J-1)*4]; _1: end; _2: end; Result := Result + Format('%s',[' A=']); Result := Result + #$0D#$0A; for J:=1 to N do begin for I:=1 to M do begin Result := Result + Format('%20.12f ',[A[(I-1)+(J-1)*4]]) + #$0D#$0A; end; end; Result := Result + #$0D#$0A; AIG3R(A,M,N,NLEAE,RCONE,Z,IERR); Result := Result + Format('%s',[' NLEAE=']); Result := Result + #$0D#$0A; for JI:=1 to N do begin Result := Result + Format('%3d ',[NLEAD[JI-1]]) + #$0D#$0A; end; Result := Result + #$0D#$0A; Result := Result + Format('%s',[' A=']); Result := Result + #$0D#$0A; for JJ:=1 to N do begin for II:=1 to M do begin Result := Result + Format('%20.12f ',[A[(II-1)+(JJ-1)*4]]) + #$0D#$0A; end; end; Result := Result + #$0D#$0A; Result := Result + Format('%s',[' IERR=']); Result := Result + Format('%3d',[IERR]) + #$0D#$0A; Result := Result + Format('%s',[' RCONE=']); Result := Result + Format('%20.12f ',[RCOND]) + #$0D#$0A; for I:=1 to M do begin for J:=1 to N do begin C[(I-1)+(J-1)*4] := 0.0; for K:=1 to N do begin C[(I-1)+(J-1)*4] := C[(I-1)+(J-1)*4]+B[(I-1)+(K-1)*4]*A[(K-1)+(J-1)*4]; _3: end; _4: end; _5: end; Result := Result + Format('%s',[' A*INVERSE(A)=']); Result := Result + #$0D#$0A; for J:=1 to N do begin for I:=1 to M do begin Result := Result + Format('%20.12f ',[C[(I-1)+(J-1)*4]]) + #$0D#$0A; end; end; Result := Result + #$0D#$0A; for I:=1 to M do begin for J:=1 to N do begin C[(I-1)+(J-1)*4] := 0.0; for K:=1 to N do begin C[(I-1)+(J-1)*4] := C[(I-1)+(J-1)*4]+A[(I-1)+(K-1)*4]*B[(K-1)+(J-1)*4]; _6: end; _7: end; _8: end; Result := Result + Format('%s',[' INVERSE(A)*A=']); Result := Result + #$0D#$0A; for J:=1 to N do begin for I:=1 to M do begin Result := Result + Format('%20.12f ',[C[(I-1)+(J-1)*4]]) + #$0D#$0A; end; end; Result := Result + #$0D#$0A; UtRes('TAIG3R',Result); {вывод результатов в файл TAIG3R.res} exit; end; end. Результаты: | 0.05056 0.05429 0.00629 0.03500 | A = | 0.05189 -0.08460 0.07212 -0.00030 | | -0.00841 0.04319 0.02021 -0.12113 | | -0.04971 0.02797 0.07900 -0.05773 | NLEAD = (2, 2, 4, 4) RCOND = 0.41764