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

Подпрограмма:  AGG2C (модуль AGG2C_p)

Назначение

Вычисление всех собственных значений и собственных вектоpов в обобщенной проблеме AX = λBX для двух комплексных матриц с помощью LZ - алгоpитма.

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

Подпрограмма AGG2C вычисляет все собственные значения и собственные векторы обобщенной проблемы Ax = λBx для комплексных матриц A, B размера N на N с помощью LZ - алгоpитма.

Информация о вычисленных собственных значениях выдается в векто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.Martin, J.H.Wilkinson, The Modified L*R - Algorithm for Complex Hessenberg Matrices, Numer. Math. 12, 1968.

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

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

procedure AGG2C(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, содержащий вычисленные собственные векторы;
  ALFA -
  BETA  
комплексные векторы длины N, B которых содержится информация о вычисленных собственных значениях; при этом собственные значения определяются формулой (1);
N - порядок исходных матриц (тип: целый);
IERR - целая переменная, служащая для сообщения об ошибках, обнаруженных в ходе работы подпрограммы; значение IERR полагается равным номеpу собственного значения, для вычисления которого потребовалось более 30 итераций; при этом собственные значения с индексами IERR + 1, ..., N вычислены правильно, а с индексами 1, ..., IERR и соответствующие собственные векторы не вычисляются.

Версии : нет

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

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

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

  Исходные матрицы A и B не сохраняются.

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

Unit TAGG2C_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc, UtRes_p, AGG2C_p;

function TAGG2C: String;

implementation

function TAGG2C: String;
var
J,I,N,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 - КОМПЛЕКСНАЯ  MATPИЦA/ И B - КОМПЛЕКСНАЯ 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;
AGG2C(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('TAGG2C',Result);  { вывод результатов в файл TAGG2C.res }
exit;
end;

end.

Результаты:

      IERR  =  0

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

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

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

      λ1  =  0.8333 + 1.993*i
      λ2  =  0.8333 - 1.993*i
      λ3  =  0.500 + 0.0*i

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

                | -0.25205+0.19169*i |
      x1  =  | -0.08799-0.72598*i  |
                |  1.00000+0.0*i        |

                | -0.25205-0.19169*i |
      x2  =  | -0.08799+0.72598*i |
                |  1.00000+0.0*i       |
                          
                |  0.0+0.0*i |
      x3  =  |  0.0+0.0*i |
                |  1.0+0.0*i |