Текст подпрограммы и версий
adb1r_p.zip   adb1e_p.zip   adb1c_p.zip
Тексты тестовых примеров
tadb1r_p.zip   tadb1e_p.zip   tadb1c_p.zip

Подпрограмма:  ADB1R (модуль ADB1R_p)

Назначение

Вычисление определителя методом Гаусса с выбором ведущего элемента по столбцу и оценка числа обусловленности вещественной ленточной матрицы, заданной в компактной форме.

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

Для заданной в компактной форме ленточной вещественной матрицы А порядка N выполняется треугольная факторизация L - 1 * А = U, где U - верхняя треугольная ленточная матрица, вычисляется величина, обратная числу обусловленности матрицы А

           RCOND = 1/( || A ||1 * || A-1 ||1 ) ,

 где   || A ||1 = maxj = 1, ..., N { | a1 j | + | a2 j | + ...+ | aN j | } ,

а затем вычисляется определитель матрицы А как произведение диагональных элементов матрицы U, умноженное на (-1)I, где I - число выполненных в процессе факторизации перестановок.

Определитель записывается в виде:

     det A = D1*10D2 ,   где 1.0 ≤ D1 < 10.0 

Дж. Форсайт, М. Малькольм, К. Моулер, Машинные методы математических вычислений. М.: Мир, 1980.

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

procedure ADB1R(var A :Array of Real; var MA :Integer; var N :Integer;
                var ML :Integer; var MU :Integer;
                var NLEAD :Array of Integer; var DET1 :Real;
                var DET2 :Real; var RCOND :Real; var Z :Array of Real;
                var IERR :Integer); 

Параметры

A - вещественный двумерный массив размера МА на N, в первых МL+МU+1 столбцах которого задается в компактном виде исходная ленточная матрица порядка N; на выходе в первых МL столбцах массива находятся нижние кодиагонали ленточной матрицы L1* ... *LN - 1 ( Li ,  i = 1, ..., N - 1 суть элементарные матрицы исключения метода Гаусса ), в следующих МL+МU+1 столбцах содержится в компактном виде матрица U (см. замечания по использованию);
MA - первая размерность массива А в вызывающей программе (тип: целый);
N - порядок матрицы А (тип: целый);
ML - число нижних кодиагоналей матрицы А (тип: целый);
MU - число верхних кодиагоналей матрицы А (тип: целый);
NLEAD - целый вектор длины N, содержащий на выходе информацию о выполненных в ходе факторизации перестановках (см. замечания по использованию);
DET1 - вещественная переменная, содержащая на выходе мантиссу определителя;
DET2 - вещественная переменная, содержащая на выходе десятичный порядок определителя;
RCOND - вещественная переменная, содержащая на выходе вычисленное значение величины, обратной числу обусловленности матрицы А;
Z - вещественный рабочий вектор длины N;
IERR - целая переменная, содержащая на выходе информацию о прохождении счета, при этом:
IЕRR=65 - если МА ≤ 0 или N ≤ 0;
IЕRR=66 - если в процессе работы произошло переполнение (это говорит о том, что либо ||А||1, либо некоторые элементы матрицы U превосходят по абсолютной величине максимальное представимое на данной машине число);
IЕRR=-К - если в результате факторизации в К-ой строке матрицы U диагональный элемент равен нулю (свидетельствует о вырожденности матрицы А). Если таких строк у матрицы U несколько, то значение K полагается равным номеру последней из них.

Версии

ADB1E - вычисление определителя методом Гаусса с выбором ведущего элемента по столбцу и оценка числа обусловленности вещественной ленточной матрицы, заданной с расширенной (Extended) точностью в компактном виде;
ADB1C - вычисление определителя методом Гаусса с выбором ведущего элемента по столбцу и оценка числа обусловленности комплексной ленточной матрицы, заданной в компактном виде.

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

AFB2R - подпрограмма треугольной факторизации и оценки числа обусловленности ленточной матрицы.
UTAFSI - подпрограмма выдачи диагностических сообщений.

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

  1. 

В подпрограмме АDВ1E массивы А, Z и переменные RСОND, DЕТ1 и DЕТ2 имеют тип Extended, для треугольного разложения и оценки числа обусловленности ленточной матрицы А вызывается подпрограмма АFВ2E.

  2. 

В подпрограмме АDВ1С массивы А и Z и переменная DЕТ1 имеют тип Сomplex, для треугольного разложения и оценки числа обусловленности ленточной матрицы А вызывается подпрограмма АFВ2С.

  3. 

На выходе К - ый элемент вектора NLЕАD равен номеру строки, переставленной на К - ом шаге факторизации с К - ой строкой матрицы А. Поскольку факторизация Гаусса требует N - 1 шагов, то NLЕАD [N] := N

  4. 

Так как в результате выполненных в ходе факторизации перестановок число верхних кодиагоналей матрицы U равно МL+МU, а также в силу некоторых конструктивных особенностей подпрограммы, для правильной ее работы необходимо выполнение условия МА ≥ N > МU+2*МL+1.

   

Если же МU+2*МL+1 ≥ N, то имеет смысл, задав матрицу А не в компактной, а в полной форме, обратиться к подпрограмме АDG2R.

  5. 

Если вырабатывается значение переменной IЕRR, отличное от нуля, то выдается соответствующее диагностическое сообщение, полагается RСОND := 0.0, DЕТ1 := 0.0, DЕТ2 := 0.0 и происходит выход из подпрограммы.

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

Unit TADB1R_p;
interface
uses
SysUtils, Math, { Delphi }
LStruct, Lfunc, UtRes_p, ADB1R_p;

function TADB1R: String; 

implementation

function TADB1R: String;
var
MA,N,I,J,ML,MU,J0,J1,K,IERR :Integer;
E,DET1,DET2,RCONE :Real;
A :Array [0..24] of Real;
Z :Array [0..4] of Real;
NLEAE :Array [0..4] of Integer;
label
_2,_1,_3,_4;
begin
Result := '';
МА := 5;
N := 5;
for I:=1 to МА do
 begin
  for J:=1 to N do
   begin
    A[(I-1)+(J-1)*5] := 0.0;
_2:
   end;
_1:
 end;
ML := 1;
MU := 1;
for I:=1 to МА do
 begin
  J0 := Max0(1,I-ML);
  J1 := Min0(N,I+MU);
  for J:=J0 to J1 do
   begin
    K := J-I+ML+1;
    A[(I-1)+(K-1)*5] := (I*10+J);
_3:
   end;
_4:
 end;
Result := Result + Format('%s',[' A=']);
Result := Result + #$0D#$0A;
for J:=1 to N do
 begin
  for I:=1 to МА do
   begin
    Result := Result + Format('%20.16f',[A[(I-1)+(J-1)*5]]) + #$0D#$0A;
   end;
 end;
Result := Result + #$0D#$0A;
ADB1R(A,MA,N,ML,MU,NLEAE,DET1,DET2,RCONE,Z,IERR);
Result := Result + Format('%s',['   NLEAE=']);
Result := Result + #$0D#$0A;
for J:=1 to N do
 begin
  Result := Result + Format('%3d',[NLEAD[J-1]]) + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[' A=']);
Result := Result + #$0D#$0A;
for J:=1 to N do
 begin
  for I:=1 to МА do
   begin
    Result := Result + Format('%20.16f',[A[(I-1)+(J-1)*5]]) + #$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.16f',[RCOND]) + #$0D#$0A;
Result := Result + Format('%s',['   DET1=']);
Result := Result + Format('%20.16f',[DET1]) + #$0D#$0A;
Result := Result + Format('%s',['   DET2=']);
Result := Result + Format('%20.16f',[DET2]) + #$0D#$0A;
UtRes('TADB1R',Result);  { вывод результатов в файл TADB1R.res }
exit;
end;

end.


Результаты:

           |    0         21.0     22.0     23.0     0 |
           | -0.524    32.0     33.0     34.0     0 |
    A = | -0.015    43.0     44.0     45.0     0 |
           |  0.292    54.0     55.0        0       0 |
           | -0.228    0.569      0         0       0 |

       NLEAD = (2, 3, 4, 5, 5)
       RCOND = 1.47362E-03
       DET1 = 8.88360
       DET2 = 5.0