Та¬кого же ко¬л¬и¬чес¬т¬¬ва дейс¬твий по-тре¬бу¬ет дальнейшее умно¬жение на мат¬¬рицу V*(?kl) в со¬от¬вет¬ствии с фор¬му¬ла¬ми (1.44), ког¬да про¬ис¬хо¬дит из¬ме¬не¬ние k-й и пер¬вой строк.
Общее количество арифметических операций для плу¬¬че¬ния каждой последующей матрицы A(n) со¬с¬та¬вит при этом ~ 6n.
Рассмотренный алгоритм (с некоторой мо¬ди¬фи¬ка¬ци¬ей, предложенной в работе [Greenstadt, 1961]) р嬬а¬ли¬зо¬ван в прце¬дуре EIGJAC, используя которую можно вы¬чис¬лить все собственные зна¬че¬ния и соб¬ствен¬ные век¬то-ры за¬дан¬¬ной сим¬мет¬рич¬ной матрицы A[1:n, 1:n]. На вы¬хде из про¬ц嬬ду¬ры l-й столбец мат¬рицы S (на¬чаль¬ное зна¬че¬ние мас¬си¬ва S не¬су¬щест¬венно) со¬дер¬жит первый из n соᬬ¬ст¬¬вен¬ных веꬬ¬то¬ров дан¬ной мат¬рицы, ди¬а¬го¬наль¬ный эле¬мент a[l, l] мас¬¬сива A сот¬вет¬ст¬вует первому соᬬст¬вен¬ному знବ¬че¬нию. Ит嬬ра¬ци¬онный процесс зବка¬н¬чи¬ва¬ет-ся, ког¬да для каж¬дого н嬬¬ди¬а¬го¬наль¬нго элемента a[i, j] вы¬пол¬ня¬¬ется не¬ра¬вен¬ст¬во
,
где параметр ?? связанный с оценкой по¬греш¬нос-ти rn в фор¬му¬ле (1.46) и введенный вместо rn для оценки ка¬чес¬т¬ва ре¬ш嬬ния в работе [Gre¬en-stadt, 1961], и есть за¬да¬ва¬емая точ¬ность.
Формальные параметры про¬цедуры. Входные: n (тип in¬teger) - порядок матрицы A; eps (тип real) - за¬дବваемая точ¬¬¬ность; A[1:n,1:n] (тип real) - ис¬ход¬ная си쬬метричная мат¬¬¬рица. Вы-ход¬ные: A[1:n, 1:n] (тип re¬al) - матрица, ди¬а¬го-наль¬¬ные элементы ко¬трой A[i, i] яв¬ляются со¬от-вет¬ст¬вен¬но ??-ми соб¬ст¬вен¬ны¬ми зна¬ч嬬н謬я¬ми; S[1:n, 1:n] (тип re¬al) - мат¬¬ри¬ца, k-й столбец кто-рой содержит k-й из n соб¬ст¬вен¬ных векторов ис-хо䬬ной матрицы A.
PROCEDURE EIGJAC (N:INTEGER;EPS:REAL;
VAR A, S:MAS1);
VAR NORM1,NORM2,THR,MU,OMEGA,SINT,COST,
INT1,V1, V2,V3,T : REAL;
I,J,P,Q,IND : INTEGER;
LABEL MAIN, MAIN1;
BEGIN
{ *** ФОРМИРОВАНИЕ В МАССИВЕ S
ЕДИНИЧНОЙ МАТРИЦЫ *** }
FOR I:=1 TO N DO
FOR J:=1 TO I DO
IF I=J THEN
S[I,J]:=1
ELSE
BEGIN
S[I,J]:=0;
S[J,I]:=0;
END;
{ *** ВЫЧИСЛЕНИЕ НАЧАЛЬНОЙ НОРМЫ NORM1,
КОНЕЧНОЙ НОРМЫ NORM2 И ПОРОГА THR *** };
INT1:=0;
FOR I:=2 TO N DO
FOR J:=1 TO I-1 DO
INT1:=INT1+2*A[I,J]*A[I,J];
IF INT1=0 THEN
EXIT;
NORM1:=SQRT(INT1);
THR:=SQRT(INT1);
NORM2:=(EPS/N)*NORM1;
IND:=0;
MAIN:
THR:=THR/N;
{**НАЧАЛО ОБРАБОТКИ НЕДИАГОНАЛЬНЫХ ЭЛЕМЕНТОВ**};
MAIN1:
FOR Q:=2 TO N DO
FOR P:=1 TO Q-1 DO
IF ABS(A[P,Q]) >= THR THEN
BEGIN
IND:=1;
V1:=A[P,P];
V2:=A[P,Q];
V3:=A[Q,Q];
MU:=0