Текст подпрограммы и версий
afp4r_p.zip , afp4c_p.zip , afp4e_p.zip 		Тексты тестовых примеров
tafp4r_p.zip , tafp4c_p.zip , tafp4e_p.zip

Подпрограмма:  AFP4R (модуль AFP4R_p)

Назначение

Приведение прямоугольной вещественной матрицы размера N*M (N ≥ M) к верхнему двухдиагональному виду методом отражений.

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

Для прямоугольной N*М (N ≥ М) матрицы А строятся две последовательности матриц отражения Q1, Q2, ..., QM  и  R1, R2, ... RM - 2  такие, что 

     QM ... Q1AR1 ... RM-2 = D ,

где D - верхняя двухдиагональная матрица.

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

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

procedure AFP4R(var A :Array of Real; N :Integer; M :Integer;
                var D1 :Array of Real; var D2 :Array of Real);


   Параметры
A - двумерный N*М массив, в котором задается исходная матрица (тип: Real); в результате работы подпрограммы в А запоминается информация о матрицах отражения; в наддиагональной части массива А в последовательных строках запоминаются векторы, порождающие матрицы отражения R1, ..., RM - 2 ; в остальной части массива в последовательных столбцах запоминаются векторы, порождающие матрицы отражения Q1, Q2, ..., QM .
N, M - заданные размеры исходной матрицы, причем N ≥ М (тип: целый);
D1 - одномерный массив длины N, используемый подпрограммой как рабочий (тип: Real); в результате работы подпрограммы в первых М компонентах D1 запоминаются элементы главной диагонали матрицы D;
D2 - одномерный массив длины М, используемый подпрограммой как рабочий (тип: Real); в результате работы подпрограммы в первых М - 1 компонентах D2 запоминаются элементы второй диагонали матрицы D. 
   Версии
AFP4C - приведение методом отражений к верхнему двухдиагональному виду прямоугольной матрицы N*М (N ≥ М) комплексной матрицы;
AFP4E - приведение методом отражений к верхнему двухдиагональному виду прямоугольной матрицы N*М (N ≥ М) вещественной матрицы, заданной с расширенной (Extended) точностью. 
   Вызываемые подпрограммы:  нет 
   Замечания по использованию
  1. 

В подпрограмме АFР4С массивы А, D1, D2 имеют тип Complex.

  2.  В подпрограмме АFР4E массивы А, D1, D2 имеют тип Extended.

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

Unit TAFP4R_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc, UtRes_p, AFP4R_p;

function TAFP4R: String;

implementation

function TAFP4R: String;
var
J,I,_i :Integer;
D1 :Array [0..5] of Real;
D2 :Array [0..3] of Real;
const
A :Array [0..23] of Real = ( 1.0,0.0,-3.0,2.0,-1.0,-6.0,4.0,1.0,1.0,-1.0,2.0,
0.0,3.0,1.0,5.0,1.0,-1.0,2.0,2.0,-1.0,6.0,1.0,0.0,
-1.0 );
begin
Result := '';  { результат функции }
AFP4R(A,6,4,D1,D2);
Result := Result + #$0D#$0A;
for I:=1 to 6 do
 begin
  for J:=1 to 4 do
   begin
    Result := Result + Format('%20.16f ',[A[(I-1)+(J-1)*6]]) + #$0D#$0A;
   end;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
for _i:=0 to 5 dО
 begin
  Result := Result + Format('%20.16f ',[D1[_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 ',[D2[_i]]);
  if ( ((_i+1) mod 4)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
UtRes('TAFP4R',Result);  { вывод результатов в файл TAFP4R.res }
exit;
end;

end.


Результат:

                 |  1.068   1.064    0.870    0.332 |
                 |  0.000  -1.045    1.310  -0.534 |
     A  =    | -0.393  -0.904    1.279   0.000 |
                 |  0.262  -0.118   -0.252   1.257 |
                 | -0.131   0.070    0.501  -0.518 |
                 | -0.787  -0.266  -0.221  -0.391 |

      D1  =   (-7.141,  7.608,  -2.625,  -4.142,  -0.518,  -0.391)

      D2  =   (-3.175,  -3.851,  0.412,  -0.534)

Это означает, что

                 | -7.141  -3.175   0.000    0.000 |
                 |  0.000   7.608  -3.851    0.000 |
      D  =    |  0.000   0.000  -2.625    0.412 |
                 |  0.000   0.000   0.000   -4.142 |
                 |  0.000   0.000   0.000     0.000 |
                 |  0.000   0.000   0.000     0.000 |

      Qi  =  I - Ui * UiT ,    i = 1, 2, 3, 4 ,
 где
      U1T  =   (1.068,  0.000,  -0.393,  0.262,  -0.131,  -0.787)

      U2T  =   (0.000,  -1.045, -0.904,  -0.118,  0.070,  -0.266)
   
      U3T  =   (0.000,  0.000,  1.279,  -0.252,  0.501,  -0.221)

      U4T  =   (0.000,  0.000,  0.000,  1.257,  -0.518,  -0.391)

      Ri  =  I - Wi * WiT ,     i = 1, 2 , 
 где
      W1T  =   (0.000,  1.064,  0.870,  0.332)

      W2T  =   (0.000,  0.000,  1.310,  -0.534)