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

117
Из первой строки уравнений находим сначала
11 11
a
u
=
, затем
11 1
1
u
a
u
j
j
=
при
n
j
,...,
2
=
. Из второй –
2 12 22 22
u
a
u
−
=
, затем
22 1
12 2
2
u
u
u
a
u
j
j
j
−
=
при
n
j
,...,
3
=
и т.д. Завершается процесс вычислением
∑
−
=
−
=
1 1
2
n
k
kn
nn
nn
u
a
u
Таким образом, матрица U может быть определена совокупностью формул
∑
−
=
−
=
1 1
2
i
k
ki
ii
ii
u
a
u
при
n
i
,...,
2
,
1
=
;
ii
i
k
kj
ki
ij
ij
u
u
u
a
u
∑
−
=
−
=
1 1
при
n
i
,...,
2
=
,
i
j
>
(
0
=
ij
u
при
i
j
<
).
(3.24)
Осуществимости вещественного
U
U
T
– разложения вещественной симметрической матрицы A могут помешать два обстоятельства: обращение в нуль элемента
ii
u при каком-либо
{
}
n
i
,...,
2
,
1
∈
и отрицательность подкоренного выражения. Известно, что для важного в приложениях класса
симметричных положительно определенных матриц разложение по
формулам (3.24) выполнимо.
При наличии
U
U
T
- разложения решение системы с симметрической матрицей
b
x
A
r r
=
сводится к последовательному решению двух систем с треугольными матрицами
b
y
U
T
r r
=
и
y
x
U
r r
=
Первая из них имеет вид








=
+
+
+
=
+
=
,
....
..........
..........
,
,
2 2
1 1
2 2
22 1
12 1
1 11
n
n
nn
n
n
b
y
u
y
u
y
u
b
y
u
y
u
b
y
u
откуда получаем вспомогательные неизвестные
n
y
y
y
,...,
,
2 1
по формулам
,
11 1
1
u
b
y
=
ii
i
k
k
ki
i
i
u
y
u
b
y
∑
−
=
−
=
1 1
(
( )
1
>
i
).
(3.25)

118
Из второй системы








=
=
+
=
+
+
+
n
n
nn
n
n
n
n
y
x
u
y
x
u
x
u
y
x
u
x
u
x
u
..........
..........
,
,
2 2
2 22 1
1 2
12 1
11
находим искомые значения
i
x в обратном порядке, т.е. при
1
,...,
1
,
−
=
n
n
i
, по формулам
,
nn
n
n
u
y
x
=
ii
n
i
k
k
ik
i
i
u
x
u
y
x
∑
+
=
−
=
1
(
)
n
i
<
(3.26)
Решение симметричных СЛАУ по формулам (3.24)-(3.26) называют
методом квадратных корней или схемой Холецкого. В случае систем с положительно определенными матрицами можно ожидать хороших результатов применения такого метода (особенно если в процессе решения делать проверку на немалость
ii
u , чтобы избежать большого роста погрешностей). В противном случае нет, например, гарантий, что в процессе разложения не появятся чисто мнимые числа, что, кстати, может не отразиться на результатах, если в алгоритме реализации метода квадратных корней предусмотреть возможность появления мнимых чисел.

3.7. Метод вращений

Как и в методе Гаусса, цель прямого хода преобразований в этом методе – приведение системы к треугольному виду последовательным обнулением поддиагональных элементов сначала первого столбца, затем второго и т.д. Делается это следующим образом.
Пусть
12
c и
12
s – некоторые отличные от нуля числа. Умножим первое уравнение системы







=
+
+
+
=
+
+
+
=
+
+
+
n
n
nn
n
n
n
n
n
n
b
x
a
x
a
x
a
b
x
a
x
a
x
a
b
x
a
x
a
x
a
......
..........
..........
,
,
2 2
1 1
2 2
2 22 1
21 1
1 2
12 1
11
(3.27)

119
на
12
c , второе – на
12
s и сложим их. Полученным уравнением заменяем первое уравнение системы. Затем первое уравнение исходной системы умножаем на
12
s
−
, второе – на
12
c и результатом их сложения заменяем второе уравнение. Таким образом, первые два уравнения системы (3.27) заменяются уравнениями
,
)
(
)
(
)
(
2 12 1
12 2
12 1
12 2
22 12 12 12 1
21 12 11 12
b
s
b
c
x
a
s
a
c
x
a
s
a
c
x
a
s
a
c
n
n
n
+
=
+
+
+
+
+
+
.
)
(
)
(
)
(
2 12 1
12 2
12 1
12 2
22 12 12 12 1
21 12 11 12
b
c
b
s
x
a
c
a
s
x
a
c
a
s
x
a
c
a
s
n
n
n
+
−
=
+
−
+
+
+
−
+
+
−

На введенные два параметра
12
c и
12
s наложим два условия:
0 21 12 11 12
=
+
−
a
c
a
s
– условие обнуления (т.е. исключения
1
x из второго уравнения);
1 2
12 2
12
=
+
s
c
– условие нормировки.
Легко проверить, что за
12
c и
12
s , удовлетворяющие этим условиям, можно принять, соответственно,
,
2 21 2
11 11 12
a
a
a
c
+
=
.
2 21 2
11 21 12
a
a
a
s
+
=
(3.28)
Эти числа можно интерпретировать как синус и косинус некоторого угла
1
α
(отсюда название «Метод вращений», так как один промежуточный шаг прямого хода такого метода может рассматриваться как преобразование вращения на угол
1
α
расширенной матрицы системы в плоскости, определяемой индексами обнуляемого элемента).
После проведенных преобразований система (3.27) принимает вид








=
+
+
+
=
+
+
+
=
+
+
=
+
+
+
,
...
...
...
...
...
,
,
,
2 2
1 1
3 3
2 32 1
31
)
1
(
2
)
1
(
2 2
)
1
(
22
)
1
(
1
)
1
(
1 2
)
1
(
12 1
)
1
(
11
n
n
nn
n
n
n
n
n
n
n
n
b
x
a
x
a
x
a
b
x
a
x
a
x
a
b
x
a
x
a
b
x
a
x
a
x
a
(3.29) где
j
j
j
a
s
a
c
a
2 12 1
12
)
1
(
1
+
=
,
n
j
,...,
1
=
,
2 12 1
12
)
1
(
1
b
s
b
c
b
+
=
,
j
j
j
a
c
a
s
a
2 12 1
12
)
1
(
2
+
−
=
,
n
j
,...,
3
,
2
=
,
2 12 1
12
)
1
(
2
b
c
b
s
b
+
−
=

120
Далее первое уравнение системы (3.29) заменяем новым, полученным сложением результатов умножения первого и третьего уравнений (3.29), соответственно, на
2 31 2
)
1
(
11
)
1
(
11 13
)
(
a
a
a
c
+
=
и
2 31 2
)
1
(
11 31 13
)
(
a
a
a
s
+
=
, а третье заменяем уравнением, полученным сложением результатов умножения тех же уравнений, соответственно, на
13
s
−
и
13
c . Получаем систему










=
+
+
+
=
+
+
+
=
+
+
=
+
+
=
+
+
+
,
...
...
...
...
...
,
,
,
,
2 2
1 1
4 4
2 42 1
41
)
1
(
3
)
1
(
3 2
)
1
(
32
)
1
(
2
)
1
(
2 2
)
1
(
22
)
2
(
1
)
2
(
1 2
)
2
(
12 1
)
2
(
11
n
n
nn
n
n
n
n
n
n
n
n
n
n
b
x
a
x
a
x
a
b
x
a
x
a
x
a
b
x
a
x
a
b
x
a
x
a
b
x
a
x
a
x
a
где
j
j
j
a
s
a
c
a
3 13
)
1
(
1 13
)
2
(
1
+
=
,
n
j
,...,
1
=
,
3 13
)
1
(
1 13
)
2
(
1
b
s
b
c
b
+
=
,
j
j
j
a
c
a
s
a
3 13
)
1
(
1 13
)
1
(
3
+
−
=
,
n
j
,...,
3
,
2
=
,
3 13
)
1
(
1 13
)
1
(
3
b
c
b
s
b
+
−
=

Проделав такие преобразования
)
1
(
−
n
раз, приходим к системе







=
+
+
=
+
+
=
+
+
+
−
−
−
−
)
1
(
)
1
(
2
)
1
(
2
)
1
(
2
)
1
(
2 2
)
1
(
22
)
1
(
1
)
1
(
1 2
)
1
(
12 1
)
1
(
11
...
...
...
...
,
,
n
n
nn
n
n
n
n
n
n
n
n
n
b
x
a
x
a
b
x
a
x
a
b
x
a
x
a
x
a
(3.30)
такого же вида, какой приняла бы система (3.27) после первого этапа преобразований прямого хода метода Гаусса:








=
+
+
+
=
+
+
+
=
+
+
+
=
+
+
+
+
...
...
...
...
,
,
,
)
1
(
)
1
(
)
1
(
3 2
)
1
(
2
)
1
(
3
)
1
(
3
)
1
(
33 2
)
1
(
32
)
1
(
2
)
1
(
2
)
1
(
23 2
)
1
(
22 1
1 3
13 2
12 1
n
n
nn
n
n
n
n
n
n
n
n
b
x
a
a
x
a
b
x
a
a
x
a
b
x
a
a
x
a
y
x
c
x
c
x
c
x
(3.31)

121
Однако в отличие от (3.31) система (3.30) обладает замечательным свойством: длина любого вектора-столбца расширенной матрицы системы
(3.30) остается такой же, как у соответствующего столбца исходной системы.
Чтобы убедиться в этом, достаточно посмотреть на результаты следующих преобразований, учитывающих условия нормировки:
.
)
(
)
(
2 2
)
(
)
(
2 2
2 1
2 2
2 12 2
12 2
1 2
12 2
12 2
2 2
12 2
1 12 12 2
1 2
12 2
2 2
12 2
1 12 12 2
1 2
12 2
)
1
(
2 2
)
1
(
1
j
j
j
j
j
j
j
j
j
j
j
j
j
j
a
a
a
s
c
a
s
c
a
c
a
a
s
c
a
s
a
s
a
a
s
c
a
c
a
a
+
=
+
+
+
=
=
+
−
+
+
+
=
+

Это равенство показывает, что сохраняется величина суммы квадратов, изменяемых на данном промежуточном шаге преобразований, пары элементов любого столбца. Так как остальные элементы столбцов при этом остаются неизменными, то, значит, на любом этапе преобразований длина столбца будет одной и той же, т.е. не будет роста элементов матрицы и, тем самым, не будет нарастать вычислительная погрешность, несмотря на то, что коэффициенты пересчитываются
)
1
(
−
n
раз.
Дальше точно так же за
)
2
(
−
n
промежуточных шага преобразуем подсистему





=
+
+
=
+
+
,
...
...
...
...
,
)
1
(
)
1
(
2
)
1
(
2
)
1
(
2
)
1
(
2 2
)
1
(
22
n
n
nn
n
n
n
b
x
a
x
a
b
x
a
x
a
системы (3.30), создавая нули под элементом
)
1
(
22
a и т.д.
В результате
)
1
(
−
n
таких этапов прямого хода исходная система (3.27) будет приведена к треугольному виду







=
=
+
+
=
+
+
+
−
−
−
−
−
−
−
−
−
...
...
,
,
)
1
(
)
1
(
)
1
(
2
)
1
(
2 2
)
1
(
22
)
1
(
1
)
1
(
1 2
)
1
(
12 1
)
1
(
11
n
n
n
n
nn
n
n
n
n
n
n
n
n
n
n
n
b
x
a
b
x
a
x
a
b
x
a
x
a
x
a

Нахождение отсюда неизвестных
1 1
,...,
,
x
x
x
n
n
−
не отличается от рассмотренного ранее обратного хода метода Гаусса.
Метод вращений применим к любой системе
b
x
A
r r
=
с невырожденной матрицей A и обладает хорошей вычислительной устойчивостью.
Метод вращений требует в четыре раза больше операций умножения, чем схема единственного деления.

122
Замечание. Для больших n схема единственного деления требует






3 3
n
O
умножений и делений, метод квадратных корней
−






6 3
n
O
, метод вращений
−






3 3
4
n
O
, метод прогонки (он будет рассмотрен ниже)
−
)
5
( n
O

123
Лекция 4. Векторные и матричные нормы. Согласованность норм. Обусловленность СЛАУ. Число обусловленности матрицы.
Вычисление определителей. Обращение матриц.
4.1. Векторные и матричные нормы. Согласованность норм.
Пусть Н – линейное пространство n -мерных векторов.
Напомним, что в Н задана норма, если каждому вектору xr из Н сопоставлено вещественное число xr , удовлетворяющее аксиомам:
1.
0
x
≥
r и
0
x
=
r тогда и только тогда, когда
0
x
=
r r
2.
x
x
α
α
=
⋅
r r
для любого числа
α
и любого x
H
∈
r
3.
x
y
x
y
+ ≤
+
r r r
r
Аналогично определяется норма матриц. Напомним, что норма вектора может быть вычислена, в частности, по одной из следующих формул:
,
1 1
∑
=
=
n
i
i
x
xr
,
1 2
2
∑
=
=
n
i
i
x
xr
,
max
1
n
i
i
x
x
≤
≤
∞
=
r а норма матрицы по формулам:
,
max
1 1
1
∑
=
≤
≤
=
n
i
ij
n
j
a
A
,
|
)
(
|
max
1 2
A
A
A
T
j
n
j
⋅
=
≤
≤
λ
.
max
1 1
∑
=
≤
≤
∞
=
n
j
ij
n
i
a
A

Определение 1. Данная норма матриц называется согласованной с данной нормой векторов, если для любой матрицы A и для любого вектора
xr справедливо неравенство
Ax
A
x
≤
⋅
r r

Нетрудно проверить, что приведенные нормы согласованы. В дальнейшем, не оговаривая особо, будем считать, что нормы матриц и векторов согласованы.
Определение
2.
Вектор
в
в
x
x
x
∆ = −
r r r называется
вектором
погрешности. Здесь через xr обозначили точное решение системы Ax b
=
r r
, а через
в
xr – результат вычислений (говорят «
x вычисленное»).
Определение 3. Вектор
в
r
Ax
b
=
−
r r
r называется вектором невязки.
Вектор невязки вычисляют для проверки правильности решения.
Близость нормы вектора невязки к нулю, очевидно, необходима для того, чтобы система линейных уравнений была решена верно. С другой стороны близость нормы вектора невязки к нулю не гарантирует правильность решения.

124
Итак, малость вектора невязки необходима для того, чтобы вычисленные значения неизвестных были близки к точным, но не достаточна. Докажем последнее утверждение. Действительно, так как
в
в
в
x
A
x
x
A
b
x
A
r
r r
r r
r r
∆
=
−
=
−
=
)
(
, то
r
A
x
в
r r
1
−
=
∆
и
r
A
x
в
r r
⋅
≤
∆
−
1

Из последнего неравенства вытекает, что из малости rr не обязательно следует малость
в
xr
∆
, т.е. близость вычисленного решения к точному решению системы
b
x
A
r r
=
Видно, что
в
x
∆
r может быть большой, если велика норма обратной матрицы
1
A
−
Последнее наблюдается у так называемых плохо обусловленных матриц.
4.2.
Обусловленность систем линейных уравнений

Под обусловленностью вычислительной задачи понимают чувствительность ее решения к малым погрешностям входных данных.
Задачу называют хорошо обусловленной, если малым погрешностям входных данных отвечают малые погрешности решения, и плохо обусловленной, если возможны большие изменения решения.
Часто оказывается возможным ввести количественную меру степени обусловленности вычислительной задачи – число обусловленности. Эту величину можно интерпретировать как коэффициент возможного возрастания погрешностей в решении по отношению к вызвавшим их погрешностям входных данных.
Пусть между погрешностями входных данных x
и решениями y
установлены неравенства
)
(
)
(
x
v
y
∆
≤
∆
∆
,
)
(
)
(
x
v
y
δ
δ
δ
≤
Величины
∆
v и
δ
v называют, соответственно, абсолютным и относительным числом обусловленности. Чаще всего под числом обусловленности
v
задачи понимают относительное число обусловленности.
Соотношение
1
>>
v
свидетельствует о реальной возможности существенного роста ошибок. Грубо говоря, если
N
v
10
≈
δ
, то порядок
N

125
показывает число верных цифр, которое может быть утеряно в результате, по сравнению с числом верных цифр входных данных.
Каково же значение v , при котором следует признать задачу плохо обусловленной? Ответ на этот вопрос существенно зависит, с одной стороны, от предъявляемых требований к точности решения, с другой
−
от уровня обеспеченности точности исходных данных. Например, если требуется найти решение с точностью 0.1%, а входная информация задается с точностью
0,02%, то уже значение
10
=
v
сигнализирует о плохой обусловленности.
Однако (при тех же требованиях к точности результата), гарантия того, что исходные данные задаются с точностью не ниже 0.0001%, означает, что и при
3 10
=
v
задача хорошо обусловлена.
При решении систем линейных алгебраических уравнений, как и при решении других задач, возникает вопрос о влиянии погрешностей во входных данных: в элементах
ij
a матрицы A и компонентах
i
b правой части
b
r
, – на результат.
В реальных условиях вычисления на ЭВМ практически всегда сопровождаются погрешностями, вызванными погрешностями перевода числовых исходных данных из десятичной системы счисления в двоичную, и погрешностями округления уже при записи информации в память ЭВМ.
Предположим, что в системе
b
x
A
r r
=
возмущены (допущены ошибки при измерениях, произведены округления при вычислениях) как правая часть, так и коэффициенты. Рассмотрим возмущённую систему
b
x
A
=
и обозначим
A
A
A
−
=
∆
,
x
x
x
r r
−
=
∆
,
b
b
b
r r
−
=
∆
. Если
A
имеет обратную матрицу
1
−
A
и выполнено условие
1 1
−
<
∆
A
A
, то матрица
A
A
A
∆
+
=
имеет обратную и справедливы оценки относительной погрешности








∆
+
∆
−
≤
∆
−
b
b
A
A
A
A
A
cond
A
cond
x
x
r r
r r
)
/
)(
(
1
)
(
1
,
)
(
)
(
1
A
A
A
x
δ
δ
−
≤
,
)
(
)
(
1
b
A
A
x
δ
δ
−
≤
, где через cond(A) обозначили произведение
1
−
⋅
A
A
Определение 4. Число
1
)
(
−
⋅
=
A
A
A
cond
называют числом обусловленности матрицы, а вернее, стандартным числом обусловленности матрицы.
Очевидно, что число обусловленности зависит от выбора нормы.

126
Чем больше число обусловленности матрицы, тем больше может быть погрешность решения системы.
Определение 5. Систему линейных алгебраических уравнений
(матрицу
A
) называют плохо обусловленной, если число обусловленности велико, в противном случае – хорошо обусловленной.
Свойства
)
( A
cond
:
1.
1
)
(
=
E
cond
2.
1
)
(
≥
A
cond
3.
)
(
)
(
)
(
B
cond
A
cond
AB
cond
≤
4. Число обусловленности не меняется при умножении матрицы
A
на ненулевое число.
5. Для симметрической матрицы
A
:
|
)
(
|
min
|
)
(
|
max
)
(
1 1
A
A
A
cond
j
n
j
j
n
j
λ
λ
≤
≤
≤
≤
=
Не существует численного метода, с помощью которого можно было бы устранить чувствительность плохо обусловленной системы к возмущениям элементов матрицы и правой части.
Решение систем линейных уравнений с плохо обусловленными матрицами может оказаться некорректной задачей.
В случае систем второго порядка понятие обусловленности матрицы допускает наглядную геометрическую интерпретацию. В этом случае каждое из уравнений графически изображается прямой линией, и если система плохо обусловлена, то прямые почти параллельны. При этом небольшое изменение угла наклона или сдвиг одной из прямых очень изменяет положение точки пересечения прямых, т.е. значительно меняет решение соответствующей линейной системы (рис. 1).

Рис.1
Реальное вычисление числа обусловленности требует нахождения
1
−
A и ее нормы. Это влечет за собой около
2 3
2n
n
+
дополнительных операций и приблизительно в четыре раза увеличивает затраты на решение системы линейных уравнений
b
x
A
r r
=

127
Поэтому можно грубо оценить
1
−
A
. Заметим, что если
y
A
z
r r
1
−
=
, то
y
A
z
r r
⋅
≤
−
1
. Следовательно,
y
z
A
r r
≥
−
1
. Можно выбрать
k
векторов
i
yr ,
k
i
,...,
1
=
, решить системы
i
i
y
z
A
r r
=
и положить
i
i
i
y
z
A
r r
max
1
≈
−
Если
k
мало, то это вычисление потребует около
3
kn
операций, что немного по сравнению с общими трудозатратами на решение
b
x
A
r r
=
Эвристически установлено, что если yr выбирать случайным образом, то ожидаемой оценкой
y
z
r r
будет
1 2
1
−
A
. Представляется надежным использовать небольшие значения k .
Традиционным примером очень плохо обусловленной матрицы является матрица Гильберта A с элементами
1 1
−
+
=
j
i
a
ij
. Из табл. 1 видно, что для матрицы A даже сравнительно невысокого порядка число обусловленности оказывается чрезвычайно большим.

Скачать 1.29 Mb.

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