Текст подпрограммы и версий
aeh3c_p.zip , aeh3z_p.zip 		Тексты тестовых примеров
taeh3c_p.zip , taeh3z_p.zip

Подпрограмма:  AEH3C (модуль AEH3C_p)

Назначение

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

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

Комплексная эрмитова матрица  А последовательностью унитарных подобных преобразований приводится к симметричной трехдиагональной форме F = V*Q*AQV (где  Q - произведение унитарных матриц,  V - унитарная диагональная матрица), используемой для вычисления собственных значений, заданных в интервале, и их номеров методом бисекций. Собственные векторы, соответствующие вычисленным собственным значениям вычисляются методом обратной итерации, при этом связь между векторами матриц  А и  F определяется следующим соотношением: x = QVy.

Уилкинсон, Pайнш. "Справочник алгоритмов на языке АЛГОЛ. Линейная алгебра", М., "Машиностроение", 1976.

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

procedure AEH3C(N :Integer; var AR :Array of Real;
                var AI :Array of Real; var VR :Array of Real;
                var VI :Array of Real; var EV :Array of Real;
                var RLB :Real; var RUB :Real; MM :Integer;
                var M :Integer; var IRAB1 :Array of Integer;
                var RAB1 :Array of Real; var RAB2 :Array of Real;
                var IERR :Integer);


   Параметры
N - порядок исходной матрицы (тип: целый);
AR, AI - вещественные двумерные массивы размера N на N, содержащие соответственно действительную и мнимую части исходной матрицы; в подпрограмме используется только полный нижний треугольник матрицы АI и строгий нижний треугольник матрицы АR;
VR, VI - вещественные двумерные массивы размера N на MM, содержащие соответственно действительную и мнимую части вычисленных ортонормированных собственных векторов;
EV - вещественный вектор длины ММ, содержащий вычисленные в возрастающем порядке собственные значения, принадлежащие заданному интервалу;
      RLB -
      RUB   		заданные нижняя и верхняя границы интервала собственных значений (тип: вещественный); если RLВ > RUВ, то собственные значения не вычисляются;
MM - оценка сверху числа собственных значений исходной матрицы, принадлежащих заданному интервалу, ММ ≤ N (тип: целый); если вычисленное число  М собственных значений, принадлежащих заданному интервалу, больше, чем ММ, то собственные значения и соответствующие собственные векторы не вычисляются;
M - целая переменная, в которой запоминается вычисленное число собственных значений, принадлежащих заданному интервалу;
IRAB - целый вектор длины ММ, содержащий индексы расположенных в возрастающем порядке  М собственных значений;
RAB1 - вещественный двумерный массив размера 2 на N, используемый как рабочий;
RAB2 - вещественный вектор длины 8 на N, используемый как рабочий;
IERR - целая переменная, служащая для сообщения об ошибках, обнаруженных в ходе работы подпрограммы; при этом:
IERR = 3*N + 1 - если значение ММ меньше числа  М вычисленных собственных значений на интервале; при этом собственные значения и соответствующие собственные векторы не вычисляются;
IERR = - К - если для вычисления собственного вектора с индексом  К потребовалось более 5 итераций; при этом компоненты этого вектора полагаются равными нулю. Если таких собственных векторов несколько, то значение IЕRR полагается равным индексу последнего из них. 
   Версии
AEH3Z - вычисление собственных значений, принадлежащих заданному интервалу, их номеров и соответствующих собственных векторов комплексной эрмитовой матрицы, заданной с расширенной (Extended) точностью. Массивы АR, АI, VR, VI, ЕV, RАВ1, RАВ2 и переменные RLВ, RUВ имеют тип Extended. 
   Вызываемые подпрограммы
UTAE10 - подпрограмма выдачи диагностических сообщений при работе подпрограмм AEH3C и AEH3Z. 
   Замечания по использованию
  Подпрограмма AEH3C сохраняет полный верхний треугольник массива АR и строгий верхний треугольник массива АI, остальные элементы массивов АR и АI используются как рабочие.

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

Unit TAEH3C_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc, UtRes_p, AEH3C_p;

function TAEH3C: String;

implementation

function TAEH3C: String;
var
J,I,N,MM,M,IERR :Integer;
RLB,RUB :Real;
VR :Array [0..11] of Real;
VI :Array [0..11] of Real;
EV :Array [0..2] of Real;
IRAB1 :Array [0..2] of Integer;
RAB1 :Array [0..7] of Real;
RAB2 :Array [0..31] of Real;
const
AR :Array [0..15] of Real = ( 1.0,0.42,0.54,0.66,0.42,1.0,0.32,0.44,0.54,0.32,
1.0,0.22,0.66,0.44,0.22,1.0 );
AI :Array [0..15] of Real = ( 0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
0.0,0.0,0.0,0.0 );
begin
Result := '';  { результат функции }
Result := Result + Format('%s',
 [' ВЫЧИСЛЕНИЕ НЕКОТОРЫХ СОБСТВЕННЫХ ЗНАЧЕНИЙ И COOTBETCTBYЮЩИХ COБCTBEHHЫX' + #$0D#$0A +
 ' ВЕКТОРОВ КОМПЛЕКСНОЙ ЭРМИТОВОЙ MATPИЦЫ']) + #$0D#$0A; 
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[' AR' + #$0D#$0A]);
for I:=1 to 4 do
 begin
  for J:= 1 to 4 do
   begin
    Result := Result + Format('%20.16f ',[AR[(I-1)+(J-1)*4]]) + #$0D#$0A;
   end;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[' AI' + #$0D#$0A]);
for I:=1 to 4 do
 begin
  for J:=1 to 4 do
   begin
    Result := Result + Format('%20.16f ',[AI[(I-1)+(J-1)*4]]) + #$0D#$0A;
   end;
 end;
Result := Result + #$0D#$0A;
N := 4;
ММ := 3;
RLB := 0.0;
RUB := 1.0;
AEH3C(N,AR,AI,VR,VI,EV,RLB,RUB,MM,M,IRAB1,RAB1,RAB2,IERR);
Result := Result + Format('%s',[' PEЗYЛЬTAT']) + #$0D#$0A;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[' AR' + #$0D#$0A]);
for I:=1 to 4 do
 begin
  for J:=1 to 4 do
   begin
    Result := Result + Format('%20.16f ',[AR[(I-1)+(J-1)*4]]) + #$0D#$0A;
   end;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[' AI' + #$0D#$0A]);
for I:= 1 to 4 do
 begin
  for J:=1 to 4 do
   begin
    Result := Result + Format('%20.16f ',[AI[(I-1)+(J-1)*4]]) + #$0D#$0A;
   end;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[' EV' + #$0D#$0A]);
for I:=1 to 3 do
 begin
  Result := Result + Format('%20.16f ',[EV[I-1]]) + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[' VR' + #$0D#$0A]);
for I:=1 to 4 do
 begin
  for J:=1 to 3 do
   begin
    Result := Result + Format('%20.16f ',[VR[(I-1)+(J-1)*4]]) + #$0D#$0A;
   end;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[' VI' + #$0D#$0A]);
for I:=1 to 4 do
 begin
  for J:=1 to 3 do
   begin
    Result := Result + Format('%20.16f ',[VI[(I-1)+(J-1)*4]]) + #$0D#$0A;
   end;
 end;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[' IERR' + #$0D#$0A]);
Result := Result + Format('%8d ',[IERR]) + #$0D#$0A;
UtRes('TAEH3C',Result);  { вывод результатов в файл TAEH3C.res }
exit;
end;

end.

Результаты:

Собственные значения в интервале (0., 1.):

                  |  0.242260708  |       
       EV  =  |  0.638283803  |   ,
                  |  0.796706689  |

Собственные векторы, соответствующие вычисленным
в интервале (0., 1.) собственным значениям:
 
                  |  - 0.718845953  |                             |  - 0.380449882   |
       V1  =  |  - 0.095698981  |   ,               V2  =  |    0.850275473   |   ,
                  |    0.387435463   |                             |    0.035889606   |
                  |    0.569206432   |                             |  - 0.361941215  |

                  |    0.050328449  |
       V3  =  |  - 0.237226458  |   ,
                  |    0.812846171  |
                  |  - 0.529595844  |

       IERR  =  0