Текст подпрограммы и версий
asg1r_p.zip , asg1e_p.zip , asg1c_p.zip 		Тексты тестовых примеров
tasg1r_p.zip , tasg1e_p.zip , tasg1c_p.zip

Подпрограмма:  ASG1R (модуль ASG1R_p)

Назначение

Решение вещественной системы линейных алгебраических уравнений А*x = b с матрицей A и вектором правой части b методом Гаусса с выбором ведущего элемента по столбцу. При повторном решении системы с той же матрицей и другой правой частью эффективно используется выполненная ранее факторизация матрицы.

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

Матрица системы А N - го порядка приводится к верхнетреугольному виду U последовательностью элементарных преобразований Гаусса L1, L2, ..., LN - 1 и матрицей перестановок Q так, что 

     LN-1*LN-2...L1*Q*A = U ,

и затем решается система с треугольной матрицей U и правой частью
LN - 1*LN - 2...L1*Q*b . Матрицы Li ,  i = 1, 2, ..., N - 1 имеют единичные диагональные элементы, матрица Q осуществляет перестановку строк матрицы А и обеспечивает стратегию выбора ведущего элемента по столбцам.

В.В.Воеводин, Р.В.Петрина, Комплекс алгоритмов, основанных на преобразованиях типа Гаусса, в пакете линейной алгебры, Сб. "Численный анализ на ФОРТРАНе", вып,3, Изд-во МГУ, 1973.

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

procedure ASG1R(var A :Array of Real; var B :Array of Real;
                var X :Array of Real; var S :Array of Integer;
                N :Integer; L :Integer);


   Параметры
A - двумерный N на N массив, в котором задается исходная матрица (тип: Real); по окончании работы подпрограммы в массиве А на сответствующих местах запоминаются элементы матрицы U и векторы-столбцы, порождающие матрицы Li ;
B - вектор длины N, в котором задается вектор правой части системы (тип: вещественный);
X - вектор длины N, в котором запоминается вектор вычисленного решения (тип: вещественный);
S - вектор длины N, в котором запоминается вектор, порождающий матрицу перестановок; при этом в S (К) запоминается номер строки переставленной на К - ом шаге с К - ой строкой (тип: целый);
N - заданный порядок матрицы системы (тип: целый);
L - признак решаемой задачи (тип: целый), причем
L = 1 - если система с данной матрицей решается впервые;
L ≠ 1 - если система с данной матрицей решается повторно. 
   Версии
ASG1E - решение системы линейных алгебраических уравнений А*x = b методом Гаусса с выбором ведущего элемента по столбцу для вещественных А и b, заданных с расширенной (Extended) точностью.
ASG1C - решение системы линейных алгебраических уравнений А*x = b методом Гаусса с выбором ведущего элемента по столбцу для комплексных А и b. 
   Вызываемые подпрограммы:  нет 
   Замечания по использованию
  1. 

В подпрограмме АSG1E массивы А, В, Х имеют тип Extended.

  2. 

В подпрограмме АSG1С массивы А, В, Х имеют тип Complex.

  3. 

Поскольку факторизация Гаусса требует N - 1 шагов, то S (N) = 0.

  4.  Подпрограммы АSG1R, АSG1E, АSG1С позволяют помещать вычисленное решение Х на место правой части В.

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

Unit TASG1R_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc, UtRes_p, ASG0R_p, ASG1R_p;

function TASG1R: String;

implementation

function TASG1R: String;
var
S,I,J,_i :Integer;
B :Array [0..3] of Real;
S2 :Array [0..3] of Integer;
РВ :Array [0..15] of Real;
B1 :Array [0..3] of Real;
B2 :Array [0..3] of Real;
X1 :Array [0..3] of Real;
X2 :Array [0..3] of Real;
const
РА :Array [0..15] of Real = ( 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
_14,_3,_11;
begin
Result := '';  { результат функции }
for I:=1 to 4 do
 begin
  for J:=1 to 4 do
   begin
    PB[(I-1)+(J-1)*4] := PA[(I-1)+(J-1)*4];
_14:
   end;
 end;
for I:=1 to 4 do
 begin
_3:
  B[I-1] := 0.0;
 end;
B[0] := 1.0;
for I:=1 to 4 do
 begin
  B1[I-1] := B[I-1];
  B2[I-1] := B[I-1];
_11:
 end;
ASG0R(PA,B1,X1,4,1);
ASG1R(PB,B2,X2,S2,4,1);
Result := Result + #$0D#$0A;
for _i:=0 to 15 do
 begin
  Result := Result + Format('%20.16f ',[PA[_i]]);
  if ( ((_i+1) mod 4)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
for _i:=0 to 3 do
 begin
  Result := Result + Format('%20.16f ',[B1[_i]]);
  if ( ((_i+1) mod 4)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
for _i:=0 to 3 do
 begin
  Result := Result + Format('%20.16f ',[X1[_i]]);
  if ( ((_i+1) mod 4)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
for _i:=0 to 15 do
 begin
  Result := Result + Format('%20.16f ',[PB[_i]]);
  if ( ((_i+1) mod 4)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
for _i:=0 to 3 do
 begin
  Result := Result + Format('%20.16f ',[B2[_i]]);
  if ( ((_i+1) mod 4)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
for _i:=0 to 3 do
 begin
  Result := Result + Format('%20.16f ',[X2[_i]]);
  if ( ((_i+1) mod 4)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
UtRes('TASG1R',Result);  { вывод результатов в файл TASG1R.res }
end;

end.


Результат:

      X  =  (0.051,  0.052,  -0.008,  0.050)