Текст подпрограммы и версий
afe1c_p.zip , afe1z_p.zip 		Тексты тестовых примеров
tafe1c_p.zip , tafe1z_p.zip

Подпрограмма:  AFE1C (модуль AFE1C_p)

Назначение

Одновременное приведение эрмитовой трехдиагональной и диагональной положительно определенной матриц двусторонним диагональным преобразованием соответственно к вещественной симметрической трехдиагональной и единичной матрицам.

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

Пусть Т - заданная эрмитова трехдиагональная матрица, D - диагональная положительно определенная матрица. Подпрограмма AFE1C получает вещественную симметрическую трехдиагональную матрицу Н и диагональную матрицу D1, такие, что  Н = D2ТD1 и  D2DD1 = Е , где D2 - диагональная матрица.

Эта подпрограмма может быть использована для сведения обобщенной проблемы собственных значений вида Тx = λDx , где D > 0, к стандартной проблеме Нz = λz с вещественной симметрической трехдиагональной матрицей Н; при этом собственные значения этих задач совпадают, а соответствующие собственные векторы связаны следующим образом:  x = D1z .

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

procedure AFE1C(var A :Array of Real; var BR :Array of Real;
                var BI :Array of Real; var D :Array of Real;
                var FR :Array of Real; var FI :Array of Real;
                N :Integer; var IERR :Integer);


   Параметры
A - вещественный вектор длины N, содержащий диагональные элементы исходной трехдиагональной матрицы Т; на выходе из подпрограммы на месте вектора А запоминаются вычисленные диагональные элементы матрицы Н;
BR, BI - вещественные векторы длины N, в последних N - 1 компонентах которых содержатся соответственно вещественные и мнимые части наддиагональных элементов исходной трехдиагональной матрицы Т; на выходе из подпрограммы в последних N - 1 компонентах вектора ВR запоминаются вычисленные наддиагональные элементы трехдиагональной матрицы Н;
D - вещественный вектор длины N, содержащий диагональные элементы исходной диагональной матрицы D;
FR, FI - вещественные векторы длины N, содержащие соответственно вещественные и мнимые части диагональных элементов диагональной матрицы перехода D1;
N - заданный порядок исходных матриц Т и D (тип: целый);
IERR - целая переменная, служащая для сообщения об ошибках, обнаруженных в ходе работы подпрограммы; при этом:
IЕRR=65 - если исходная диагональная матрица D не является положительно определенной. 
   Версии
AFE1Z - одновременное приведение эрмитовой трехдиагональной и диагональной положительно определенной матриц, заданных с расширенной (Extended) точностью, двусторонним диагональным преобразованием соответственно к вещественной симметрической трехдиагональной и единичной матрицам. 
   Вызываемые подпрограммы
UTAF12 - подпрограмма выдачи диагностических сообщений при работе подпрограмм АFЕ1С и АFЕ1Z. 
   Замечания по использованию
  В подпрограмме АFЕ1Z параметры А, ВR, ВI, D, FR, FI должны иметь тип Extended.

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

Unit tafe1c_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc, UtRes_p, AFE1C_p;

function tafe1c: String;

implementation

function tafe1c: String;
var
_i,IERR :Integer;
FR :Array [0..3] of Real;
FI :Array [0..3] of Real;
const
A :Array [0..3] of Real = ( 4.0,4.0,4.0,4.0 );
BR :Array [0..3] of Real = ( 0.0,0.0,-4.0,2.0 );
BI :Array [0..3] of Real = ( 0.0,2.0,0.0,2.0 );
D :Array [0..3] of Real = ( 1.0,4.0,4.0,2.0 );
begin
Result := '';  { результат функции }
AFE1C(A,BR,BI,D,FR,FI,4,IERR);
Result := Result + Format('%s',[' A,BR,FR,FI,IERR:']) + #$0D#$0A; 
Result := Result + #$0D#$0A;
for _i:=0 to 3 do
 begin
  Result := Result + Format('%20.16f ',[A[_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 ',[BR[_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 ',[FR[_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 ',[FI[_i]]);
  if ( ((_i+1) mod 4)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + Format(' %5d ',[IERR]) + #$0D#$0A;
UtRes('tafe1c',Result);  { вывод результатов в файл tafe1c.res }
exit;
end;

end.


Результаты:

      A   =   (4., 1., 1., 2.) ,     BR  =   (0., 1., 1., 1.) ,

      FR  =   (1., 0., 0., 0.5) ,    FI  =   (0., -0.5, 0.5, 0.5) ,

      IERR  =  0 .