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

Подпрограмма:  AGE0R (модуль AGE0R_p)

Назначение

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

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

Обобщенная проблема Ax = λBx с вещественной симметрической трехдиагональной матрицей A и диагональной положительной матрицей B сводится к стандартной задаче Qy = λy. Здесь Q = B - 1/2AB - 1/2 - симметричная трехдиагональная матрица, собственные значения которой совпадают с собственными значениями исходной обобщенной проблемы, а собственные векторы этих задач связаны соотношением  x = B - 1/2y.

Для решения стандартной задачи используется неявный QL - алгоритм со сдвигом. Собственные векторы  x  исходной обобщенной проблемы нормированы таким образом, что  xTВx = Е,  где E - единичная матрица.

Дж.Х. Уилкинсон, Алгебраическая проблема собственных значений, "Hаука", M., 1970.

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

procedure AGE0R(var A :Array of Real; var B :Array of Real;
                var V :Array of Real; var EV :Array of Real;
                var RAB :Array of Real; N :Integer;
                var IERR :Integer);

Параметры

A - вещественный двумерный массив размера N на 2, содержащий элементы симметрической трехдиагональной матрицы A в компактной форме записи симметричных ленточных матриц;
B - вещественный вектоp длины N, содержащий элементы диагональной матрицы B;
V - вещественный двумерный массив размера N на N, в столбцах которого содержатся вычисленные оpтонормированные собственные векторы;
EV - вещественный вектоp длины N, содержащий вычисленные собственные значения в возрастающем порядке;
RAB - вещественный вектоp длины N, используемый как рабочий;
N - порядок исходных матриц (тип: целый);
IERR - целая переменная служащая для сообщения об ошибках, обнаруженных в ходе работы подпрограммы; значение IERR полагается равным:
 

номеpу собственного значения, для вычисления которого потребовалось более 30 итераций, при этом собственные значения с индексами 1, 2, ..., IERR - 1 и соответствующие собственные векторы вычислены правильно, но собственные значения неупорядочены;

IERR = 65 - когда матрица B не является положительно определенной.

Версии : нет

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

AFE1R - одновременное приведение вещественной симметрической трехдиагональной матрицы и диагональной положительно определенной матрицы двусторонним диагональным преобразованием соответственно к вещественной симметрической трехдиагональной и единичной матрицам;
UTAG10 -
UTAG11  
подпрограммы выдачи диагностических сообщений при работе подпрограммы AGE0R.

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

  Подпрограмма AGE0R сохраняет матрицу A.

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

Unit TAGE0R_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc, UtRes_p, AGE0R_p;

function TAGE0R: String;

implementation

function TAGE0R: String;
var
J,I,N,_i,IERR :Integer;
V :Array [0..24] of Real;
RАВ :Array [0..4] of Real;
EV :Array [0..4] of Real;
const
A :Array [0..9] of Real = ( 0.0,2.0,2.0,2.0,2.0,0.0,0.0,0.0,0.0,0.0 );
B :Array [0..4] of Real = ( 2.0,2.0,2.0,2.0,2.0 );
begin
Result := '';  { результат функции }
Result := Result + #$0D#$0A;
Result := Result + Format('%s',
 [#$0D#$0A + ' ВСЕ C6БСТВЕННЫЕ ЗНАЧЕНИЯ И СОБСТВЕННЫЕ ВЕКТОРЫ ' + #$0D#$0A +
 ' УРАВНЕНИЯ BИДA:AX=ЛЯMБДA*BX, ГДЕ A - ВЕЩЕСТВЕННАЯ ' + #$0D#$0A +
 ' СИММЕТРИЧЕСКАЯ ТРЕХДИАГОНАЛЬНАЯ MATPИЦA,' + #$0D#$0A +
 ' A B - ВЕЩЕСТВЕННАЯ ДИАГОНАЛЬНАЯ MATPИЦA' + #$0D#$0A +
 ' A' + #$0D#$0A]);
for I:=1 to 5 do
 begin
  for J:=1 to 2 do
   begin
    Result := Result + Format(' %20.16f ',[A[(I-1)+(J-1)*5]]) + #$0D#$0A;
   end;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[#$0D#$0A + ' B' + #$0D#$0A]);
for I:=1 to 5 do
 begin
  Result := Result + Format('%20.16f ',[B[I-1]]) + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
N := 5;
AGE0R(A,B,V,EV,RAB,N,IERR);
Result := Result + Format('%s',[#$0D#$0A + ' PEЗУЛЬTAT']) + #$0D#$0A; 
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[#$0D#$0A + ' A' + #$0D#$0A]);
for I:=1 to 5 do
 begin
  for J:=1 to 2 do
   begin
    Result := Result + Format(' %20.16f ',[A[(I-1)+(J-1)*5]]) + #$0D#$0A;
   end;
 end;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[#$0D#$0A + ' EV' + #$0D#$0A]);
Result := Result + #$0D#$0A;
for _i:=0 to 4 do
 begin
  Result := Result + Format('%20.16f ',[EV[_i]]);
  if ( ((_i+1) mod 4)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[#$0D#$0A + ' V' + #$0D#$0A]);
for I:=1 to 5 do
 begin
  for J:=1 to 5 do
   begin
    Result := Result + Format(' %20.16f ',[V[(I-1)+(J-1)*5]]) + #$0D#$0A;
   end;
 end;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[#$0D#$0A + ' IERR' + #$0D#$0A]);
Result := Result + Format('%8d ',[IERR]) + #$0D#$0A;
UtRes('TAGE0R',Result);  { вывод результатов в файл TAGE0R.res }
exit;
end;

end.

Результаты:

                | -1.7320508 |                    |  0.2041241 |
                | -1.0000000 |                    | -0.3535534 |
   EV    =  |  0.0000000 |    ,     V1  =  |  0.4082483 |    , 
                |  1.0000000 |                     | -0.3535534 |
                |  1.7320508 |                     |  0.2041241 |

               | -0.3535534 |                      | -0.4082483 |
               |  0.3535534 |                      |  0.0000000 |
    V2  =   |  0.0000000 |    ,      V3  =  |  0.4082483 |    , 
               | -0.3535534 |                      |  0.0000000 |
               |  0.3535534 |                       | -0.4082483 |

               | -0.3535534 |                      | -0.2041241 |
               | -0.3535534 |                      | -0.3535534 |
    V4  =  | -0.0000000 |    ,       V5  =  | -0.4082483 |    , 
               |  0.3535534 |                       | -0.3535534 |
               |  0.3535534 |                       | -0.2041241 |
 
     IERR  =  0.