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

Подпрограмма:  AGJ1R (модуль AGJ1R_p)

Назначение

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

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

Обобщенная проблема Ax = λBx с якобиевой матрицей А и диагональной положительной матрицей B последовательностью преобразований сводится к стандартной задаче Qz = λz. Здесь Q - симметричная трехдиагональная матрица, собственные значения которой совпадают с собственными значениями исходной обобщенной проблемы. Для решения этой задачи применяется неявный QL - алгоритм со сдвигом.

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

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

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


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

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

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

IERR = 66 - когда в исходной несимметрической трехдиагональной матрице A не все попарные произведения соответствующих элементов побочных диагоналей неотрицательны. 
   Версии :   нет 
   Вызываемые подпрограммы
AFJ0R - приведение матрицы Якоби к симметрической трехдиагональной матрице преобразованием подобия.
AFE1R - одновременное приведение вещественной симметрической трехдиагональной матрицы и диагональной положительно определенной матрицы двусторонним диагональным преобразованием соответственно к вещественной симметрической трехдиагональной и единичной матрицам.
UTAG10 -
UTAG11   		подпрограммы выдачи диагностических сообщений при работе подпрограммы AGJ1R. 
   Замечания по использованию
  Подпрограмма AGJ1R сохраняет матрицу A.

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

Unit TAGJ1R_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc, UtRes_p, AGJ1R_p;

function TAGJ1R: String;

implementation

function TAGJ1R: String;
var
J,I,N,_i,IERR :Integer;
EV :Array [0..4] of Real;
RАВ :Array [0..9] of Real;
const
A :Array [0..14] of Real = ( 0.0,1.0,1.0,1.0,1.0,0.0,0.0,0.0,0.0,0.0,4.0,4.0,
4.0,4.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 + ' ВСЕ СОБСТВЕННЫЕ ЗНАЧЕНИЯ УРАВНЕНИЯ BИДA:AX=LAMBDA*BX,' + #$0D#$0A +
 ' ГДЕ A - ВЕЩЕСТВЕННАЯ МАТРИЦА ЯКОБИ,' + #$0D#$0A +
 ' B - ВЕЩЕСТВЕННАЯ ДИАГОНАЛЬНАЯ МАТРИЦА ' + #$0D#$0A +
 ' A' + #$0D#$0A]);
for I:=1 to 5 do
 begin
  for J:=1 to 3 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;
AGJ1R(A,B,EV,RAB,N,IERR);
Result := Result + Format('%s',
 [#$0D#$0A + ' PEЗУЛЬTAT' + #$0D#$0A]) + #$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 3 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;
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 3)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[#$0D#$0A + ' IERR' + #$0D#$0A]);
Result := Result + Format('%8d ',[IERR]) + #$0D#$0A;
UtRes('TAGJ1R',Result);  { вывод результатов в файл TAGJ1R.res }
exit;
end;

end.

Результаты:

                    | -1.7320508 |
                    | -1.0000000 |
       EV    =  |  0.0000000 |    , 
                    |  1.0000000 |
                    |  1.7320508 |

       IERR  =  0