<< Предыдущая

стр. 113
(из 131 стр.)

ОГЛАВЛЕНИЕ

Следующая >>

12 1,n?1 1,n
A = A0 = ? · · ?.
· ··· ·
an,1 an,2 · · · an,n?1 an,n
0 0 0 0


Здесь и ниже верхний индекс нумерует последовательность матриц, возникающих
в процессе вычислений.
Алгоритм приведения включает следующие действия:
1) Начинаем с предпоследнего элемента нижней строки a0n,n?1 . Если этот эле-
мент не равен нулю — делим на него все элементы столбца, а из остальных
столбцов вычитаем (n ? 1)-й столбец, умноженный на такие величины, чтобы в
нижней строке получались нули. В результате получаем матрицу вида
? ?
b11 b12 · · · b1,n?1 b1,n
B=? · · ?.
· ··· ·
0 ···
0 1 0

Эта матрица не подобна матрице A0 . Чтобы восстановить подобие, нужно еще
выполнить некоторое, также элементарное, преобразование, эквивалентное умно-
жению слева на матрицу
? ?
···
1 0 0 0
?0 0?
···
1 0
? ?
M =? · · ?.
· ··· ·
?0 ?
? an,1 a0 a0 ?
· · · a0 n,n
n,2 n,n?1
···
0 0 0 1
Это преобразование изменяет лишь (n ? 1)-ю строку матрицы B. Возникшая ма-
трица
?1 ?
a11 a1 · · · a1 a1
12 1,n?1 1,n
A1 = M B = ? · ·?
· ··· ·
0 ···
0 1 0
500 В.И. Фущич, В.В. Корняк

уже будет подобна матрице A0 . После этого переходим к (n ? 1)-й строке и повто-
ряем операцию с элементом, расположенным в позиции (n ? 1, n ? 2).
2) Если очередной элемент в позиции (n?k, n?k ?1) оказался равным нулю —
ищем слева от него в (n ? k)-й строке ненулевой элемент. Допустим, что он распо-
ложен в l-м столбце, — переставляем (n ? k ? 1)-й и l-й столбец и одновременно,
для сохранения подобия, (n ? k ? 1)-ю и l-ю строки.
3) Если ненулевой элемент в (n ? k)-й строке не найден, то матрица Ak прио-
бретает блочно-треугольный вид
C D
Ak = ,
0 F
где
? ?
···
ak ak ak
n?k,n?k n?k,n?1 n?k,n
? ?
···
1 0 0
F =? ?,
? ?
· ··· · ·
···
0 1 0
т.е. F имеет фробениусову форму.
Поскольку

det(A ? ?E) = det(Ak ? ?E) = det(C ? ?E) det(F ? ?E),

сразу можно выделить det(F ? ?E) как сомножитель характеpистического поли-
нома и далее применить алгоритм к матрице C.
В соответствии с этим алгоритмом была написана программа NLP. Эта про-
грамма вычисляет и корни полиномов, если их степень оказывается меньше пяти.
Для более экономного представления разреженных матриц в памяти ЭВМ мы при-
меняем представление строк матрицы в виде сумм, слагаемыми которых является
произведения элементов матрицы на координаты “пробного” вектора. Ниже будет
приведен поясняющий пример.
После определения собственных чисел матрицы A и их алгебраических кратно-
стей необходимо найти жорданову форму J и матрицу преобразования W , такую,
что

W AW ?1 = J.

Для этих целей была написана программа NLJW. Алгоритм работает следующим
образом:
Выполняется цикл последовательного перебора корней. Выбрав определенный
корень, мы должны установить структуру жорданова блока, соответствующего
этому корню. Жордановым блоком мы здесь будем называть совокупность всех
жордановых клеток с данным корнем ?. Жордановой клеткой порядка p называют
матрицу размера p ? p вида
? ?
? 1 0 ··· 0 0
? 0 ? 1 ··· 0 0 ?
? ?
? · · · ··· · · ?.
? ?
? 0 0 0 ··· ? 1 ?
0 0 0 ··· 0 ?
Реализация на ЭВМ алгоритма вычисления нелокальных симметрий 501

В частности, диагональному случаю соответствуют жордановы клетки порядка 1.
Существует общая формула, позволяющая найти число np жордановых клеток
порядка p в блоке с собственным числом ?, именно
np = rank (A ? ?E)p?1 ? 2 rank (A ? ?E)p + rank (A ? ?E)p+1 .
По определению полагается rank (A ? ?E)0 = N , где N — размер матрицы A. Не-
посредственное использование этой формулы может привести к излишним вычи-
слениям. Во многих случаях (на практике — в большинстве) для определения
числа np нет необходимости вычислять rank (A ? ?E)p+1 . Вычисление высоких
степеней матрицы A ? ?E является самым узким местом алгоритма. Например,
если A ? ?E — разреженная матрица, то уже (A ? ?E)2 как правило не будет
разреженной. Поэтому важно минимизировать количество вычислений степеней
матриц. Мы получим некоторые соотношения, позволяющие достичь этого.
Поскольку все рассматриваемые ниже величины не зависят от выбора системы
координат, рассмотрим матрицу A ? ?E в жордановом базисе. Она будет иметь
вид
?
? ?
?
?
0
? ?
.. n
? ?
. ?1
? ?
?
? ?
0
? ?
?
? ?
?
01
? ?
?
?
? ?
?
?
00
? ?
?
?
A ? ?E = ? ?,
..
? ?
. n2
? ?
?
? ?
01? ?
? ?
?
?
? ?
?
?
? ?
00
? ?
? ?
..
? ?
.
B
здесь B — матрица полного ранга rank B = N ? kalg , kalg = n1 + 2n2 + · · · +
mnm — алгебраическая кратность корня ?. Последовательно возведение в степень
нильпотентной матрицы вида
? ?
0 1 0 ··· 0 0
? 0 0 1 ··· 0 0 ?
? ?
? · · · ··· · · ?.
? ?
? 0 0 0 ··· 0 1 ?
0 0 0 ... 0 0
приводит к сдвигам на 1 позицию вправо диагонали из единиц и уменьшению ран-
га на 1. Исхода из этого, рассматривая последовательно ранги степеней матрицы
A ? ?E можно получить необходимые нам величины. Допустим, что уже найдены
числа n1 , n2 , . . . , np?1 и вычислен rank (A ? ?E)p , тогда нам известны величины
a = pnp + (p + 1)np+1 + · · · + mnm = kalg ? [n1 + 2n2 + · · · + (p ? 1)np?1 ],
b = np+1 + 2np+2 + · · · + (m ? p)nm = kalg ? N + rank (A ? ?E)p .
Поскольку a, b, nl — неотрицательные числа, то в случаях когда b = 0 или
b = 1 величины np , np+1 , . . . определяются однозначно. Именно, если b = 0, то
502 В.И. Фущич, В.В. Корняк

np = a/p, np+1 = 0, . . .; если b = 1, то np = (a ? p ? 1)/p, np+1 = 0, . . . (Можно
объединить формулу вычисления np для обоих случаев b = 0 и b = 1: np = [a?b(p+
1)]/p). При p = 1, например, случай b = 0 соответствует равенству геометрической
кратности алгебраической, что приводит к диагонализации блока, а равенство b =
1 означает, что геометрическая кратность на единицу меньше алгебраической, что
указывает на наличие в жордановом блоке диагональной части длины kalg ? 2 и
одной жордановой клетки вида
?1
.
0?
Таким образом, в подобных случаях для определения np достаточно знать только
rank (A ? ?E)p . Отдельно необходимо выделить случай алгебраической кратности
равной единице — этот случай вообще не требует вычисления рангов.
В процессе определения жордановой структуры матрицы можно по частям ре-
шать матричное уравнение W A ? JW = 0 для матрицы преобразования W . Имен-
но, можно показать, что каждой жордановой метке соответствует автономная си-
стема линейных уравнений, содержащая только те элементы матрицы W , которые
входят в строки W , расположенные против соответствующих строк жордановой
клетки. Более того, если в блоке имеется несколько жордановых клеток одного
размера, то им соответствуют системы уравнений одинаковой структуры. Поэтому
достаточно решать систему уравнений один раз и, заменив переменные, получить
для всех клеток данного размера решения. Проиллюстрируем это примером для
матриц размера 2 ? 2. Пусть матрица J имеет вид
? 0
J= ,
0 ?
тогда уравнение W A ? JW = 0 будет иметь вид
w11 w12 a11 a12 ? 0 w11 w12
? = 0.
w21 w22 a21 a22 0 ? w21 w22
Для верхней строки матрицы W получаем систему уравнений:
w11 (a11 ? ?) + w12 a21 = 0,
w11 a12 + w12 (a22 ? ?) = 0,
для нижней — систему
w21 (a11 ? ?) + w22 a21 = 0,
w21 a12 + w22 (a22 ? ?) = 0.

Мы видим, что эти системы одинаковы с точностью до замены неизвестных. Ана-
логичная ситуация и в случае жордановых клеток общего вида.
Таким образом алгоритм производит следующие действия:
1) Выполняется цикл по различным корням.
2) После выбора корня выполняется цикл по размерам жордановых клеток. Для
каждого размера определяется количество клеток, и, если оно не равно нулю —
методом Гаусса решается уравнение для определения соответствующих элементов
матрицы свободными параметрами, а остальные выражаются через них.
Реализация на ЭВМ алгоритма вычисления нелокальных симметрий 503

3) Последовательно выводятся на печать по клеткам строки матриц J и W .
Производится подсчет числа произвольных параметров.
Такая последовательная организация вычислений может оказаться полезной
в случае, если известны не все корни. В таком случае мы получим частичное
приведение матрицы к жорданову виду, что может дать некоторую информацию о
симметриях системы уравнений.
В конце работы программа печатает число параметров матрицы преобразова-
ния. Параметры должны удовлетворять условиям невырожденности матрицы. Не-
вырожденность матрицы соответствует случаю общего положения, т.е. из любой
вырожденной матрицы малым изменением параметров можно получить невыро-
жденную.
Мы не вычисляем здесь обратную матрицу W ?1 по следующим причинам:
1) Элементы обратной матрицы будут слишком громоздкими, т.к. они представ-
ляют собой сложные нелинейные функции параметров.
2) Даже в численной математике непосредственное вычисление любого объе-
кта, получаемого с помощью обратных матриц требует в 3 раза меньше времени
и в 2 раза меньше памяти, чем обращение матрицы. Поэтому с вычислительной
точки зрения обратные матрицы бесполезны.
Для получения выражений, содержащих обратные матрицы целесообразней на-
писать программу, вычисляющую эти выражения непосредственно, например, ме-
тодом Гаусса.

2. Примеры применения программ
Приведем примеры, иллюстрирующие работу программ NLP и NLJW. Все
вычисления проводились на ЭВМ ЕС–1022 с доступной оперативной памятью 390
килобайтов.
Рассмотрим вначале систему уравнений с минимальным нетривиальным разме-

<< Предыдущая

стр. 113
(из 131 стр.)

ОГЛАВЛЕНИЕ

Следующая >>