Текст подпрограммы и версий
agt2c_p.zip 
Тексты тестовых примеров
tagt2c_p.zip 

Подпрограмма:  AGT2C (модуль AGT2C_p)

Назначение

Вычисление всех собственных значений и собственных вектоpов в обобщенной поблеме AV - 1x = λBV - 1x для комплексных верхней почти треугольной матрицы A, верхней треугольной матрицы B и произвольной невырожденной матрицы V, с помощью LR - алгоритма.

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

Подпрограма AGT2C вычисляет все собственные значения и собственые векторы обобщенной проблемы

     AV - 1 = λBV - 1x 

для комплексных верхней почти треугольной матрицы A, верхней треугольной матрицы B и произвольной невырожденной матрицы V с помощью LR - алгоритма. Матрицы A, B и V имеют размеры N на N.

Информация о вычисленных собственных значениях выдается в вектоpах ALFA и BETA длины N, по которым k - ое собственное значение λk определяется с помощью формулы:

                  λk = ALFA(k) / BETA(k),   при BETA(k) ≠ 0,
 (1)            λk = ∞                            ,   при BETA(k) = 0,          
                  λk = любое число          ,   при ALFA(k) = BETA(k) = 0. 

Вычисленные собственные векторы запоминаются в столбцах матрицы V.

R.S.Mfrtin, J.H.Wilkinson, The Modified LR - Algorithm for Complex Hessenberg Matrices, Numer. Math., 12, 1968.

C.B.Moler, G.W.Stewart, An Algorithm for the Generalized Matrix Eigenvalue Problems, SIAM J. Numer. Anal., 10, 1973.

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

procedure AGT2C(var A :Array of Complex; var B :Array of Complex;
                var V :Array of Complex; var ALFA :Array of Complex;
                var BETA :Array of Complex; N :Integer;
                var IERR :Integer);

Параметры

A, B - комплексные двумерные массивы размера N на N, содержащие соответственно верхнюю почти треугольную и верхнюю треугольную матрицы;
V - комплексный двумерный массив размера N на N, содержащий заданную невырожденную матрицу V; в результате работы подпрограммы в столбцах массива V содержатся вычисленные собственные векторы Xk;
  ALFA -
  BETA  
комплексные векторы длины N, содержащие информацию о собственных значениях λk обобщенной проблемы; при этом собственные значения определяются формулой (1);
N - заданный порядок исходных матриц A, B и V (тип: целый);
IERR - целая переменная, служащая для сообщения об ошибках обнаруженных в ходе работы подпрограммы; значение IERR полагается равным номеpу собственного значения для вычисления которого потребовалось более 30 итераций; при этом собственные значения с индексами IERR+1, ..., N вычислены правильно, а с индексами 1, ..., IERR и соответствующие собственные векторы не вычисляются.

Версии : нет

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

UTAG10 - подпрограмма выдачи диагностических сообщений при работе подпрограммы AGT2C.

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

  1. 

Подпрограмма AGT2C не сохраняет исходную информацию.

  2. 

B массиве V должна задаваться единичная матрица, если решается обобщенная проблема Ax = λBx для комплексных верхней почти треугольной матрицы A и верхней треугольной матрицы B.

  3.  Вычисленные собственные векторы нормализованы так, что наибольшая компонента по модулю pавна 1.

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

Unit TAGT2C_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc, UtRes_p, AFG3C_p, AGT2C_p;

function TAGT2C: String;

implementation

function TAGT2C: String;
var
J,I,N,M,K,IERR :Integer;
V :Array [0..8] of Complex;
ALFA :Array [0..2] of Complex;
ВЕТА :Array [0..2] of Complex;
LAMBDA :Array [0..2] of Complex;
const
A :Array [0..8] of Complex = ( ( re:1.0; im:0.0 ),( re:-10.0; im:0.0 ),( 
re:5.0; im:0.0 ),( re:0.5; im:0.0 ),( re:2.0; 
im:0.0 ),( re:1.0; im:0.0 ),( re:0.0; im:0.0 ),( 
re:0.0; im:0.0 ),( re:0.5; im: 0.0 ) );
B :Array [0..8] of Complex = ( ( re:0.5; im:0.0 ),( re:3.0; im:0.0 ),( re:4.0; 
im:0.0 ),( re:0.0; im:0.0 ),( re:3.0; im:0.0 ),( 
re:0.5; im:0.0 ),( re:0.0; im:0.0 ),( re:0.0; 
im:0.0 ),( re:1.0; im: 0.0 ) );
label
_11,_10;
begin
Result := '';  { результат функции }
Result := Result + Format('%s',
 [' ОПРЕДЕЛЕНИЕ СОБСТВЕННЫХ ЗНАЧЕНИЙ И СОБСТВЕННЫХ BEKTOPOB' + #$0D#$0A +
 ' YРАВНЕНИЯ BИДA:AX=LAMBDA*BX, ГДЕ A - КОМПЛЕКСНАЯ BEPXHЯЯ' + #$0D#$0A +
 ' ПОЧТИ TPEYГОЛЬНАЯ MATPИЦA,' + #$0D#$0A +
 ' B - КОМПЛЕКСНАЯ BEPXHETPEYГОЛЬНАЯ MATPИЦA']) + #$0D#$0A; 
Result := Result + Format('%s',[' A' + #$0D#$0A]);
for I:=1 to 3 do
 begin
  for J:=1 to 3 do
   begin
    Result := Result + Format('%20.16f %20.16f ',
 [A[(I-1)+(J-1)*3].re,A[(I-1)+(J-1)*3].im]) + #$0D#$0A;
   end;
 end;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[' B' + #$0D#$0A]);
for I:=1 to 3 do
 begin
  for J:=1 to 3 do
   begin
    Result := Result + Format('%20.16f %20.16f ',
 [B[(I-1)+(J-1)*3].re,B[(I-1)+(J-1)*3].im]) + #$0D#$0A;
   end;
 end;
Result := Result + #$0D#$0A;
N := 3;
M := 0;
AFG3C(A,B,V,N,M);
AGT2C(A,B,V,ALFA,BETA,N,IERR);
Result := Result + Format('%s',[' PEЗYЛЬTAT' + #$0D#$0A]) + #$0D#$0A; 
Result := Result + Format('%s',[' A' + #$0D#$0A]);
for I:=1 to 3 do
 begin
  for J:=1 to 3 do
   begin
    Result := Result + Format('%20.16f %20.16f ',
 [A[(I-1)+(J-1)*3].re,A[(I-1)+(J-1)*3].im]) + #$0D#$0A;
   end;
 end;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[' B' + #$0D#$0A]);
for I:=1 to 3 do
 begin
  for J:=1 to 3 do
   begin
    Result := Result + Format('%20.16f %20.16f ',
 [B[(I-1)+(J-1)*3].re,B[(I-1)+(J-1)*3].im]) + #$0D#$0A;
   end;
 end;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[' ALFA' + #$0D#$0A]);
for I:=1 to 3 do
 begin
  Result := Result + Format('%20.16f %20.16f ',
 [ALFA[I-1].re,ALFA[I-1].im]) + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[' BETA' + #$0D#$0A]);
for I:=1 to 3 do
 begin
  Result := Result + Format('%20.16f %20.16f ',
 [BETA[I-1].re,BETA[I-1].im]) + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[' V' + #$0D#$0A]);
for I:=1 to 3 do
 begin
  for J:=1 to 3 do
   begin
    Result := Result + Format('%20.16f %20.16f ',
 [V[(I-1)+(J-1)*3].re,V[(I-1)+(J-1)*3].im]) + #$0D#$0A;
   end;
 end;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[' IERR' + #$0D#$0A]);
Result := Result + Format('%3d ',[IERR]) + #$0D#$0A;
for K:=1 to 3 do
 begin
  if ( CReal(BETA[K-1]) = 0.0 ) 
   then goto _11;
  LAMBDA[K-1] := DivC(ALFA[K-1],BETA[K-1]);
  goto _10;
_11:
  LAMBDA[K-1] := Cmplx(3.4E38,3.4E38);
_10:
 end;
Result := Result + Format('%s',[' LAMBDA=ALFA/BETA' + #$0D#$0A]);
for I:=1 to 3 do
 begin
  Result := Result + Format('%20.16f %20.16f ',
 [LAMBDA[I-1].re,LAMBDA[I-1].im]) + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
UtRes('TAGT2C',Result);  { вывод результатов в файл TAGT2C.res }
exit;
end;

end.

Результаты:

      IERR  =  0

                      |  2.79963 - 3.65018i  |
      ALFA  =  |  0.42142 + 0.54945i |
                      | -1.09877 + 0.0i         |
 
                      | -1.05899 - 1.89748i  |
      BETA  =  | -0.15941 + 0.27809i | 
                      | -2.19753 + 0.0i         |

               | -0.25205 + 0.19169i |
               | -0.25205 - 0.19169i  |
               |  0.0 + 0.0i                 |
               | -0.08799 - 0.72598i |
      V  =  | -0.08799 - 0.72598i |
               |  0.0 + 0.0i                 |
               |  1.00000 + 0.0i         |
               |  1.00000 + 0.0i         |
               |  1.0 + 0.0i                 |

Собственные значения   λk = ALFA(k) / BETA(k)

      λ1  =  0.83333 + 1.9930i
      λ2  =  0.83331 - 1.9931i
      λ3  =  0.50000 + 0.0i

Собственные векторы

                 | -0.25205 + 0.19169i |
      X1  =  | -0.08799 - 0.72598i  |
                 |  1.00000 + 0.0i         |
 
                 | -0.25205 - 0.19169i |
      X2  =  | -0.08799 - 0.72598i |
                 |  1.00000 + 0.0i        |

                 |  0.0 + 0.0i |
      X3  =  |  0.0 + 0.0i |
                 |  1.0 + 0.0i |