Методическое пособие для практических и лабораторных работ по курсу «Численные методы»

160 0
lim lim
=
=
−
∞
→
∞
→
A
k
k
A
k
k
z
x
x
r r
r
(8.16)
Доказательство. Доказательство проведем, следуя [12]. Перейдем от
(8.14) к задаче для вектора погрешности
0 1
=
+
−
+
k
k
k
z
A
z
z
B
r r
r
τ
(8.17)
Запишем тождество
τ
τ
k
k
k
k
k
z
z
z
z
z
r r
r r
r
−
−
+
=
+
+
1 1
2
)
(
2 1
, тогда (8.17) перепишется в виде
0
)
(
2 1
)
2
(
1 1
=
+
+
−
−
+
+
k
k
k
k
z
z
A
z
z
A
B
r r
r r
τ
τ
(8.18)
Так как
T
A
A
=
, то
2 2
1 1
1
))
(
),
(
(
A
k
A
k
k
k
k
k
z
z
z
z
z
z
A
r r
r r
r r
−
=
−
+
+
+
+

Умножив (8.18) скалярно на
)
(
2 1
k
k
z
z
r r
−
+
и используя последнее равенство, получаем
2 2
1 1
1
,
)
2
(
2
A
k
A
k
k
k
k
k
z
z
z
z
z
z
A
B
r r
r r
r r
=
+






−
−
−
+
+
+
τ
τ
τ
τ
(8.19)
В силу (8.15) и
0
>
τ
из (8.19) вытекает, что
2 0
2 2
1
A
A
k
A
k
z
z
z
r r
r
≤
≤
≤
+
Поэтому последовательность
2
A
k
zr не возрастает, ограничена снизу нулем и, следовательно, является сходящейся. Это дает основание перейти к пределу в
(8.19) и заключить, что
0
)
(
),
)(
2
(
lim
1 1
=






−
−
−
+
+
∞
→
k
k
k
k
k
z
z
z
z
A
B
r r
r r
τ

Последнее выражение и условие (8.15) влекут за собой предельный переход

161 0
lim
)
,
(
lim
2 2
1 1
1
=
−
=
−
−
+
∞
→
+
+
∞
→
k
k
k
k
k
k
k
k
z
z
z
z
z
z
r r
r r
r r

Обратимся теперь к (8.17) и перепишем его так
τ
k
k
k
z
z
B
A
z
A
r r
r
−
−
=
+
−
1 2
1 2
1

Тогда
2 2
2 1
2 1
2
τ
k
k
A
k
z
z
B
A
z
r r
r
−
≤
+
−

Переход к пределу в последнем неравенстве и завершает доказательство теоремы.
Эта теорема позволила получить доказательство сходимости ряда итерационных методов, в том числе и релаксационных методов, в которых на каждом шаге происходит «подавление» одной, наиболее сильно меняющейся компоненты.
Среди специалистов по прикладной математике широко известен метод последовательной верхней релаксации – SOR-метод:
i
i
j
n
i
j
j
k
ij
j
k
ij
i
k
ii
i
k
ii
b
x
a
x
a
x
a
x
a
ω
ω
ω
ω
+
−
−
−
=
∑
∑
−
=
+
=
+
+
1 1
1
)
(
)
(
1
)
(
)
(
1
)
1
(
(8.20)
Здесь
ω
−
это параметр, такой, что
2 0
<
< ω
. При
2 1
<
< ω
говорят о верхней релаксации, при
1 0
<
< ω
– о нижней. Значение
1
=
ω
соответствует полной релаксации или методу Зейделя. В последнее время при любых
2 0
<
<ω
его часто называют методом последовательной верхней релаксации.
Приведем (8.20) к каноническому виду (8.14). Вспомним разложение матрицы
2 1
A
D
A
A
+
+
=
. Тогда (8.20) перепишется в виде
b
x
A
x
x
A
D
k
k
k
r r
r r
=
+
−
+
+
ω
ω
1 1
)
(
(8.21)
Теорема 5. Если в (8.21)
0
>
=
T
A
A
,
2 0
<
< ω
,
(8.22)
то SOR-метод (8.21) сходится:

162 0
lim
=
∞
→
A
k
k
zv

Доказательство. Для доказательства этой теоремы нужно показать, что
A
A
D
ω
ω
2 1
1
>
+
и воспользоваться Теоремой 4. В силу (8.22)
2 1
A
A
T
=
,
0
>
D
и поэтому
(
)
(
)
(
)
(
)
0
,
2 1
1
,
2 1
,
2 1
1
,
2 1
2 1
2 1
,
2 1
2 2
1 2
1 1
>





 −
=
−





 −
=
=










−
+





 −
=












+
+
−
+
x
x
D
x
x
A
x
x
D
x
x
A
A
D
x
x
A
D
A
A
D
r r
r r
r r
r r
r r
ω
ω
ω
ω
ω
ω
ω
ω

Теорема доказана.
В последнее время получили широкое применение попеременно- треугольные методы. Теория этих методов связана с работами
А.А.Самарского. Основная идея методов заключается в следующем.
Представим матрицу A в виде суммы двух треугольных матриц
1
R и
2
R :
2 1
2 1
2 1
2 1
R
R
A
D
D
A
A
+
=






+
+





 +
=

По матрицам
1
R и
2
R строятся треугольные матрицы
1
B и
2
B , например так:
1 1
B
E
R
ω
= +
,
2 2
B
E
R
ω
= +
,
0
ω >
конструируют двухслойный итерационный метод
1 1
2 1
k
k
k
k
x
x
B B
Ax
b
τ
+
+
−
+
=
r r
r r
(8.23)
Переход от k -ой итерации к (
1)
k
+
-ой в (8.23) реализуется путем последовательного обращения треугольных матриц
1
B и
2
B

163 1
k
B y
b
Ax
= −
r r
r
,
1 2
1
k
k
k
x
x
B
y
τ
+
+
−
=
r r
r

Среди трехслойных итерационных методов выделяют метод сопряженных градиентов. Этот метод в случае симметричной и положительно определенной матрицы сходится за конечное число итераций.
Приведем этапы и формулы, реализующие метод сопряженных градиентов
(иногда его называют методом Ланцоша):
1.
Задаем
0
xr и
ε
2.
Вычисляем
0 0
x
A
b
r r
r
−
=
ξ
(невязка начального приближения).
3.
Полагаем
0 0
ξ
r r
=
p
,
0
=
k
(
k – номер итерации).
4.
Вычисляем скаляр
)
,
/(
)
,
(
k
k
k
k
k
p
A
p
r r
r r
ξ
ξ
α =
5.
Вычисляем вектор
k
k
k
k
p
x
x
r r
r
α
+
=
+
1
(очередное приближение).
6.
Вычисляем
k
k
k
k
p
Ar r
r
α
ξ
ξ
−
=
+
1
(невязка
)
1
(
+
k
-го приближения.
Действительно:
k
k
x
A
b
r r
r
−
=
ξ
,
k
k
k
k
k
k
k
k
p
A
p
A
x
A
b
x
A
b
r r
r r
r r
r r
α
ξ
α
ξ
−
=
−
−
=
−
=
+
+
1 1
. Такое выражение невязки
1
+
k
ξ
r позволяет обходиться без вычисления вектора
1
+
k
x
Ar
).
7.
Проверяем выполнение неравенства
ε
ξ
≤
+
1
k
r
. Если «да», то останавливаем работу алгоритма и выводим на печать результаты. В противном случае переходим к следующему пункту.
8.
Вычисляем скаляр
)
,
/(
)
,
(
1
k
k
k
k
k
p
A
p
p
A
r r
r r
+
= ξ
β
9.
Вычисляем вектор
k
k
k
k
p
p
r r
r
β
ξ −
=
+
+
1 1
(новое направление минимизации).
10.
Полагаем
1
+
=
k
k
и возвращаемся к вычислению скаляра
1
k
α
+
В настоящее время наиболее развиты итерационные методы для решения систем с симметричными и положительно определенными матрицами. Иногда этого достаточно, т.к., домножая слева Ax b
=
r r
на
T
A , мы приходим к эквивалентной системе
T
Ax
A b
=
r
)r с симметричной и положительно определенной матрицей
T
A
A A
=
)
. Однако, указанным образом поступают сравнительно редко. Дело в том, что при переходе от A к A
)
может быть утеряно свойство разряженности. Кроме того, вообще говоря, существенно ухудшается обусловленность системы. Имеется реальная опасность, что преобразованная система окажется очень плохо обусловленной.

164
Лекция 9. Полная и частичная проблема собственных значений. Прямые и итерационные методы. Степенной метод вычисления максимального по модулю собственного числа. Метод скалярных произведений. Методы исчерпывания.

9.1.
Прямые методы. Обусловленность задач на собственные значения

Численные методы решения проблемы собственных значений
x
x
A
r r
λ
=
,
(9.1) использовавшиеся до конца 1940-х годов, сводились в конечном счете к решению характеристического или иначе «векового» уравнения
1 2
1 2
1 0
n
n
n
n
n
p
p
p
p
λ
λ
λ
λ
−
−
−
−
−
− −
−
=
(9.2)

При реализации такого подхода основные усилия были направлены на разработку эффективных методов быстрого вычисления коэффициентов характеристического уравнения. Методы такого класса получили название
прямых, к ним относятся пользовавшиеся популярностью методы Крылова,
Данилевского, Леверье и др.
В методе Крылова, исходя из некоторого вектора
0
cr , строится последовательность векторов
1
i
i
c
Ac
−
=
r r
,
1,...,
i
n
=
и формируется система линейных алгебраических уравнений
1
n
i n i
n
i
p c
c
−
=
=
∑
r r
. Решая эту систему уравнений, получаем в зависимости от ее ранга либо коэффициенты характеристического уравнения, либо коэффициенты минимального
многочлена вектора
0
cr . Метод опирается на следующие теоретические факты: а) корнями минимального многочлена ( )
ψ λ
матрицы
A являются все различные корни характеристического многочлена, т.е. все различные собственные числа матрицы A; б) минимальный многочлен
0
( )
c
ψ λ
r любого вектора
0
cr является делителем минимального многочлена ( )
ψ λ
матрицы
A .
В методе Данилевского матрица
A
подобными преобразованиями за
(
1)
n
−
шаг приводится к канонической форме Фробениуса

165 1
2 3
1 1
0 0
0 0
0 1
0 0
0
...
...
...
0 0
0 1
0
n
n
p
p
p
p
p
−








Φ =









Элементами первой строки матрицы
Φ
являются коэффициенты характеристического многочлена матрицы A .
Метод Леверье заключается в следующем. Пусть
1
,...,
n
λ
λ
– корни характеристического многочлена, среди которых могут быть и равные.
Обозначим
1
n
k
i
k
i
s
λ
=
=
∑
. Тогда справедливы соотношения, известные под названием формул Ньютона
1 1
1 1
k
k
k
k
kp
s
p s
p s
−
−
= −
− −
,
1,...,
k
n
=
. Так как с другой стороны известно, что
k
k
s
SpA
=
, то приведенные формулы Ньютона позволяют найти все коэффициенты
i
p . (Напомним, что SpA – это след
матрицы A .)
Однако подход, основанный на решении характеристического уравнения становится неудовлетворительным, если речь идет о вычислении собственных значений матриц, имеющих порядок n в несколько десятков и тем более сотен, т.е. матриц довольно скромных по нынешним понятиям размеров.
Одна из основных причин неудовлетворительности указанного подхода состоит в том, что задачи (9.1) и (9.2), хотя формально эквивалентны, имеют разную обусловленность. Корни многочлена
)
(x
P
n
высокой степени чрезвычайно чувствительны к погрешностям в коэффициентах. Поэтому уже на этапе вычисления коэффициентов характеристического уравнения может быть значительно потеряна информация о собственных значениях матрицы.
С появлением ЭВМ широкое распространение получили итерационные методы решения проблемы собственных значений, не использующие вычисление характеристического многочлена.
Отметим, что задача вычисления собственных значений симметричных матриц очень хорошо обусловлена, что вытекает из следующего факта.
Пусть A
– матрица с приближенно заданными элементами
ij
ij
a
a
≈
,
j
λ
- собственные числа матрицы A
. Справедливо следующее утверждение:
Если A и A
– симметричные матрицы, а
j
λ
и
j
λ
– их собственные числа, упорядоченные по возрастанию, то справедливы оценки погрешности

166 2
1
|
|
max
A
A
j
j
n
j
−
≤
−
≤
≤
λ
λ

В то же время нет аналогичного утверждения для несимметричных матриц, что показывает следующий пример, принадлежащий Дж. Х.
Уилкинсону.
Собственными числами верхней треугольной матрицы 20-го порядка






















=
1 0
0 0
0 0
20 2
0 0
0 0
.
.
.
0 0
17 0
0 0
0 0
20 18 0
0 0
0 0
20 19 0
0 0
0 0
20 20
A
являются числа
1 1
=
λ
,
2 2
=
λ
, …,
20 20
=
λ
. Добавим малое
10 10
−
=
ε
к элементу
0 1
,
20
=
a
. Собственные числа возмущенной матрицы таковы:
996 0
1
≈
λ
,
11 2
2
≈
λ
,
57 2
3
≈
λ
,
i
09 1
97 3
5
,
4
±
≈
λ
,
i
95 1
89 5
7
,
6
±
≈
λ
,
i
53 2
12 8
9
,
8
±
≈
λ
,
i
73 2
5 10 11
,
10
±
≈
λ
,
i
53 2
9 12 13
,
12
±
≈
λ
,
i
95 1
1 15 15
,
14
±
≈
λ
,
i
09 1
0 17 17
,
16
±
≈
λ
,
4 18 18
≈
λ
,
9 18 19
≈
λ
,
0 20 20
≈
λ

Большинство собственных значений оказалось полностью искаженным. Матрица очень плохо обусловлена по отношению к проблеме собственных значений.
Отметим, что число обусловленности
)
( A
cond
матрицы A не характеризует обусловленность матрицы по отношению к проблеме собственных значений.
Задачи вычисления собственных чисел и собственных векторов матрицы можно разделить на две группы: полная проблема собственных значений (поиск всех собственных чисел и соответствующих им собственных векторов) и частичная проблема собственных значений (вычисление одного или нескольких собственных чисел и собственных векторов).
9.2. Степенной метод вычисления максимального по модулю
собственного числа

167
Рассмотрим степенной (итерационный) метод решения частичной проблемы собственных значений.
Изложим его в простейшем случае. Пусть матрица A действительная с положительными коэффициентами.
Обозначим через
1 2
,
, ...,
n
λ λ
λ
собственные числа матрицы, а через
1 2
,
, ...,
n
x x
x
r r r
– соответствующие им собственные векторы.
Пусть справедливы неравенства
1 2
n
λ
λ
λ
>
≥ ≥
(9.3)
Требуется найти максимальное по модулю собственное число матрицы
A , т.е.
1
λ
и соответствующий ему собственный вектор
1
xr .
Выбираем произвольно вектор
0
yr и строим последовательность
1
k
k
y
Ay
−
=
r r
,
1, 2,...
k
=
На каждой итерации вычисляем отношение
( )
( )
( )
1 1
i
k
k
i
k
y
y
λ
−
=
, i – любое.
Пусть разложение вектора
0
yr по собственным векторам
1 2
,
, ...,
n
x x
x
r r r
следующее
0 1
n
j
j
j
y
c x
=
=
∑
r r
, тогда
0 1
1
n
n
j
j
j
j
j
j
j
Ay
c Ax
c
x
λ
=
=
=
=
∑
∑
r r
r
. Обозначим через
{
}
1
, ...,
n
e
e
r r
базис из единичных векторов в
n -мерном векторном пространстве.
Пусть
1
n
j
ji i
i
x
x e
=
=
∑
r r
,
1, ...,
j
n
=
– разложение собственных векторов по этому базису. Тогда
i
n
j
n
j
n
i
n
i
n
j
k
j
ji
j
i
ji
k
j
j
j
k
j
j
k
k
e
x
c
e
x
c
x
c
y
A
y
r r
r r
r
∑
∑
∑
∑ ∑
=
=
=
=
=






=
=
=
=
1 1
1 1
1 0
λ
λ
λ
Поэтому
∑
=
=
n
j
k
j
ji
j
i
k
x
c
y
1
)
(
λ
, а
∑
=
+
+
=
n
j
k
j
ji
j
i
k
x
c
y
1 1
)
(
1
λ
. Рассмотрим отношение этих компонент
k
n
i
ni
n
k
i
i
k
n
i
ni
n
k
i
i
k
n
ni
n
k
i
k
i
k
n
ni
n
k
i
k
i
i
k
i
k
x
c
x
c
x
c
x
c
x
c
x
c
x
c
x
c
x
c
x
c
x
c
x
c
x
c
x
c
y
y




+
+




+




+
+




+
=
+
+
+
+
+
+
=
+
+
+
+
+
+
1 1
1 1
2 1
1 2
2 1
1 1
1 1
1 2
1 1
2 2
1 2
2 2
1 1
1 1
1 2
2 2
1 1
1 1
)
(
)
(
1 1
1
...
λ
λ
λ
λ
λ
λ
λ
λ
λ
λ
λ
λ
λ
λ
λ

Отсюда в силу (9.3) получаем, что

168
)
(
)
(
1 1
lim
i
k
i
k
k
y
y
+
→∞
=
λ
(9.4)
Итак, для достаточно больших k








+
=
−
k
i
k
i
k
O
y
y
1 2
1
)
(
1
)
(
λ
λ
λ
(9.5)
Чтобы найти собственный вектор, рассмотрим
k
yr .
,
2 1
1 1
1 1
1 2
1 1
1 0












+
=
=
+
=
=
=
∑
∑
∑
=
=
=
n
j
j
k
j
j
k
j
n
j
n
j
k
j
j
k
j
k
j
j
k
k
x
c
c
x
c
x
c
x
c
x
c
y
A
y
r r
r r
r r
r
λ
λ
λ
λ
λ
λ
поэтому, для достаточно больших k
1 1
1
x
c
y
k
k
r r
λ
≈

Так как собственный вектор находится с точностью до множителя, то в качестве
1
xr можно взять
k
yr .
Замечание. Иногда предела нет. Получаются прыгающие значения отношения
)
(
)
(
1
i
k
i
k
y
y
+
. В этом случае расчеты можно прекратить с этим
0
yr .
Чтобы не появлялось больших чисел, рекомендуется на каждой итерации или через несколько итераций производить нормировку векторов
k
yr . Метод работает, когда
1
λ
– комплексное, а также в случае, когда
1
λ
– кратное.
Изложенный метод и есть степенной метод вычисления максимального по модулю собственного числа (иногда его называют PM-алгоритм – от англ.
Power method).
Чтобы избавиться от некоторых недостатков степенного метода
(возможное появление больших чисел, малых чисел в знаменателе, большое число итераций) на практике чаще применяют его модификацию: метод
скалярных произведений (SP-метод, от англ. Scalar product).
Строим две последовательности векторов
1
−
=
k
k
y
A
y
r r
и
1
−
′
=
′
k
T
k
y
A
y
r r
,
,...
2
,
1
=
k
,
0 0
y
y
′
=
r r
(9.6)
Нетрудно показать, что

169
(
)
(
)
(
)
(
)
0 0
1 0
0 1
1
,
,
,
,
y
A
y
A
y
A
y
A
y
y
y
y
k
T
k
k
T
k
k
k
k
k
r r
r r
r r
r r
−
−
=
′
′
≈
λ
,
(9.7) причем
(
)
(
)








+
=
′
′
−
k
k
k
k
k
O
y
y
y
y
2 1
2 1
1
,
,
λ
λ
λ
r r
r r
,
(9.8) т.е. этот метод сходится в два раза быстрее, чем степенной.
Для симметрических матриц
(
)
(
)
0 0
1 0
0 1
,
,
y
A
y
A
y
A
y
A
k
k
k
k
r r
r r
−
≈
λ
Если
1 2
3
n
λ
λ
λ
λ
>
>
≥ ≥
, то приближенное значение второго по модулю собственного числа
2
λ
можно найти по формуле
( )
( )
1 1
2
( )
( )
1 1
2
i
i
k
k
i
i
k
k
y
y
y
y
λ
λ
λ
−
−
−
−
≈
−
Степенной метод позволяет также вычислять наименьшее по модулю собственное значение матрицы. Для этого после того, как было определено максимальное по модулю собственное значение, необходимо решить аналогичную задачу для матрицы
E
A
B
max
λ
−
=
. Обозначим максимальное собственное значение этой матрицы через
B
λ
. Тогда минимальное по модулю собственное значение матрицы A будет равно max min
λ
λ
λ
+
=
B
. Вектора, соответствующие значениям
B
λ
матрицы B и min
λ
матрицы A , будут при этом равны.
Рассмотрим еще один способ построения наименьшего по модулю собственного числа матрицы A – метод обратных итераций

Скачать 1.29 Mb.

Поделитесь с Вашими друзьями: