тк+1
^fc+l
wk,rk)
(.Awk,wk)'
тк+1 (wk,rk) 1 \
тк (wk~1,rk-l)akJ
к = 1,2..., ct\ — 1.
Этот метод наиболее широко используется в вычислительной практике при решении задач с симметричной положительно определенной матрицей.
Упражнения
Упражнение 5.1 Напишите программу, реализующую приблиэюеииое решение системы линейных алгебраических уравнений итерационным методом Зейделя. С ее помощью найдите решение системы уравнений с трехдиагональной матрицей
Ах = /,
в которой
вц 2, &г,г+1 1 QJ, {2>г,г—1 — 1 “Ь О!,
а правая часть
fo 1 fi б) i 2,3, ...,П 1, fn — 1 "h qj,
определяет точное решение х^ = 1, г = 1,2, ...,п. Исследуйте зависимость числа итераций от п и параметра а при 0 < а < 1.
В модуле seidel функция seidelO обеспечивает итерационное решение системы линейных уравнений методом Зейделя, результатом выступает приближенное решение и число итераций для достижения необходимой точности (нормы невязки).
:.У::' л::::":;;|;::;: ■
import numpy as np
import math as mt
def seidel(A, f, tol = 1.0e-9):
и и и
Solve the linear system Ax = b by Seidel method
и и и
n = len(f)
x = np.zeros((n), ’float’) r = np.copy(f)
# the maximum number of iterations is 10000 for k in ranged, 10001) : for i in range(n):
x[i] = x[i] + (f [i] - np.dot(A[i,0:n] , x[0:n])) / A[i,i] for i in range(n):
r[i] = f[i] - np.dot(A[i,0:n], x[0:n]) if np.dot(r,r) < tol**2: return x, k print ’Seidel method failed to converge:’
print ’after 10000 iteration residual = ’, mt.sqrt(np.dot(r,r)) return x, к
Решение модельной задачи дается следующей программой.
%хёг -5 1 РУ
import numpy as np from seidel import seidel n = 10 al = 0.5
A = np.zeros((n, n), ’float’) for i in range(0,n):
A [i, i] = 2 if i > 0:
A[i,i-1] = - 1 + al if i < n-1:
A[i,i+1] * - 1 -.al print ’A:\n’, A f = np.zeros((n), ’float’) f CO] = 1 - al f[n-l] = 1 + al print ’f:\n’, f x, iter = seidel(A, f) print ’iter =’, iter print ’x:\n’, x
1
iiilil
|
-1 5
|
0.
|
0.
|
|
: '-M
|
0 ; !
|
0.
|
0 i
|
■Hi
|
piss
|
1 iiiiiii
|
-1 5
|
0
|
|
|
0 ;
|
0 .
|
, 0
|
HiH
|
Iiilil
|
—0.5
|
2
|
-15
|
0 ■
|
|
0
|
йО;;й]
|
0 ■. P
|
■ill!
|
(siiiis
|
issliii
|
-0 5
|
•Mh;
|
-1 5
|
: oSp
|
o. :S
|
ШИ
|
Illlll
|
illlll
|
illlill
|
1ЙК
|
|
-0.5
|
I; 2 pi
|
-1.5
|
-■oil
|
■ill
|
: ip
|
0
|
0
|
0
|
iliis
|
!l:!l
|
-0.5
|
2.
|
lilil
|
Illlill
|
SliSlII
|
0
|
pill
|
illlll
|
mm
|
.■lllll:
|
0 :
|
-0 5
|
’ 2 ■
|
Pil 51
|
siiiis
|
; 0
|
iiilil
|
iiiSlll
|
0
|
■lli
|
0
|
' '0 S »■
|
unis
|
■liii
|
■Ill
|
■III
|
piii
|
iltilllil
|
■Illlll
|
iiilil
|
■ill
|
iiiiiii
|
Illlll
|
Illlll
|
iiisti
|
lill
|
Illlll
|
m§m
|
Illlll
|
Iiilil
|
ills
|
111
|
Illlll
|
111!
|
iilli
|
lill
|
lill
|
illlill
|
IlSlll
|
Illlll
|
Illlll
|
11111
|
Illlll
|
Illlll
|
111
|
lill
|
llllp:
|
iiiiiss
|
|
|
|
|
|
|
|
|
1111!
|
Illlill
|
ISSlll
|
IIIIISS
|
iiilil
|
Illlll
|
Illlll
|
lill
|
|
|
11
Данные по числу итераций для рассматриваемой задачи с диагональным преобладанием суммированы в табл. 5.1
Таблица 5.1 Метод Зейделя
(\
|
n = 25
|
о
Ю
II
e
|
n= 100
|
0
|
1233
|
4482
|
> 1000
|
0.25
|
289
|
410
|
599
|
0.5
|
107
|
159
|
254
|
0.75
|
58
|
93
|
160
|
1
|
25
|
50
|
100
|
При а ф () матрица системы уравнений несимметрична. Наблюдается вполне ожидаемое увеличение числа итераций при увеличении размерности задачи (параметр п) и сильная зависимость от параметра несимметричности а.
Упражнение 5.2 Напишите программу, реализующую приближенное решение системы линейных алгебраических уравнений с симметричной поло- Э1сительио определенной матрицей методом сопряженных градиентов. С ее помощью найдите решение системы
Ах = f
с матрицей Гильберта
Q>ij = ~ j I 7 5 ^ 2,..., 71, j 1,2,п
и правой частью
п
fi = ^ ^ Q>ij ? * 1, 2, ..., 77,
3=1
для которой точное решение есть ж» = 1, i = 1,2, Исследуйте зависимость числа итераций от п.
Расчетные формулы метода сопряженных градиентов возьмем (сравни с (5.24) при В = Е) в виде:
Xм = xk+1 - ak+lsk, 5fc+1 = rk+1 + (3k+lsk, k = 0,1,...
при задании итерационных параметров ak+i,/3k+1 в виде
_ (rk,rk) _ (rfe+i>rfc+i) (rk+1,Ask)
°k+l ~ (sk,Ask)’ Л+1 = (rk,rk) = ~ (sk,Ask) '
В модуле eg функция eg О реализует описанный вариант метода сопряженных градиентов для системы с самосопряженной положительно определенной матрицей А. На выходе мы имеем приближенное решение и число итераций для достижения необходимой точности по невязке.
import mimpy as np
def cg(A, f, tol = 1.0e-9):
к и и
Solve the linear system Ax = b by conjugate gradient method
и и и
n = len(f)
x = np.zeros((n), ’float1) r = np.copy(f) for i in range(n):
r[i] = np.dot(A[i,0:n] , x[0:n]) - f[i] s = np.copy(r)
As = np.zeros((n), ’float’) for к in range(n):
for i in range(n):
As[i] = np.dot(A[i,0:n], s[0:n]) alpha = np.dot(r, r) / np.dot(s, As) x = x - alpha*s for i in range(n):
r[i] = np.dot(A[i,0:n], x[0:n]) - f[i] if np.dot(r, r) < tol**2: break else:
beta = - np.dot(r, As) / np.dot(s, As) s = r + beta*s return x, к
Решение тестовой задачи дается следующей программой.
import numpy as np from eg import eg n = 4
A = np.zeros((n, n), ’float’) for i in range(0,n):
for j in range(0,n):
A[i,j] = l./(i+j+l) print ’A:\n’, A f = np.ones((n), ’float’) for i in range(0,n): f[i] = 0.
for j in range(0,n):
f [i] = f[i] + A [i, j] print ’f:\n’, f x, iter = cg(A, f) print ’iter =’, iter print ’x:\n’, x
| [ 1 : 0.33333333 0 25 ]
[ 0.33333333 0,25 : Q 16666667]
Зависимость числа итераций от размерности матрицы прослеживается по следующим данным: при п = 25 необходимо 13 итераций, при п — 50 — 17, а при п = 100 — 20.
Задачи
Задача 5.1 Напишите программу, реализующую приближенное решение системы линейных алгебраических уравнений итерационным методом Якоби. С ее помощью найдите решение системы уравнений, рассмотренной в упразднении 5.1. Сравните скорости сходимости итерационных методов Якоби и Зейделя при различных параметрах п и а.
Задача 5.2 Напишите программу, реализующую приближенное решение системы линейных алгебраических уравнений методов релаксации (5.9). Исследуйте зависимость скорости сходимости этого итерационного метода от итерационного парсиметра Tk+\ = г при численном решении системы уравнений из упразднения 5.1 при различных п и а.
Задача 5.3 Напишите программу для решения системы линейных алгебраических уравнений явным (В = Е) итерационным методом минимальных поправок (итерационные параметры определяются согласно (5.23)). Сравните скорость сходимости этого метода со скоростью сходимости метода сопряженных градиентов на примере задачи из упражнения 5.2.
Задача 5.4 Пусть
А = А1 + А2 = А*> о, А\ = А, = l-D + L,
В попеременно-треугольном методе переобуславливатель В зададим в виде В = (Е + изА\)(Е +
Напишите программу для решения системы линейных алгебраических уравнений с симметричной положительно определенной матрицей попеременно-треугольным методом с выбором итерационных параметров по (5.23). Исследуйте зависимость скорости сходимости этого итерационного метода от параметра из в задаче из упрсю1спепия 5.2.
Спектральные задачи линейной алгебры
Важной и трудной задачей линейной алгебры является нахождение собственных значений и собственных векторов матриц. Рассматривается проблема устойчивости собственных значений по отношению к малым возмущениям элементов матрицы. Для приближенного нахождения отдельных собственных значений широко используется степенной метод в различных модификациях. Для решения полной проблемы для симметричных матриц применяется итерационный метод Якоби и QR-алгоритм.
Основные обозначения
|
X \Х{ } {^1) 3-2 > • • • > Хп }
|
— n-мерный вектор
|
А={ац)
|
— матрица с вещественными элементами ац
|
Е
|
— единичная матрица
|
D = diag{di,d2,.. .,dn}
|
— диагональная матрица
|
Aj, 2 — 1,2, . • • ,72
|
— собственные значения
|
<р\ г = 1,2,... ,n
n
|
— собственные вектора
|
(x, y) = ^2XiVi
i= 1
|
— скалярное произведение
|
Собственные значения и собственные вектора матриц
Рассматриваются проблемы нахождения собственных значений и собствен' ных векторов квадратной вещественной матрицы А. Собственным числом называется число Л такое, что для некоторого ненулевого вектора (собствен' ного вектора) р имеет место равенство
Собственные вектора определены с точностью до числового множителя. Множество всех собственных значений матрицы А называется спектром матрицы А.
С учетом того, что ищется нетривиальное решение уравнения (6.1), то
det(i4 - АЕ) = 0. (6.2)
Тем самым собственные значения А матрицы А являются корнями характеристического многочлена n-й степени (6.2). Задача отыскания собственных значений и собственных векторов матрицы сводится к построению характеристического многочлена, отысканию его корней и последующему нахождению нетривиальных решений уравнения (6.1) для найденных собственных значений.
В вычислительной практике рассматривается как полная проблема собственных значений, когда необходимо найти все собственные значения матрицы А, так и частичная проблема собственных значений, когда ищутся только некоторые собственные значения, например, максимальные (минимальные) по модулю.
Численные методы решения задач на собственные значения
Начнем с приведением некоторых полезных фактов о свойствах собственных значений и собственных векторов квадратной матрицы. Далее рассматриваются методы решения частичной и полной проблемы собственных значений.
Свойства собственных значений и собственных векторов
Квадратная вещественная матрица порядка п имеет п собственных значений, при этом каждое собственное значение считается столько раз, какова его кратность как корня характеристического многочлена. Для симметричной матрицы А собственные значения вещественны, а собственные вектора, соответствующие различным собственным значениям, ортогональны, т.е. (<рг, ip*) =0, если г ф j.
Отметим также некоторые свойства собственных значений и собственных векторов для сопряженной матрицы А*:
А
(6.3)
*ф = рф.
Для спектральных задач (6.1), (6.3) имеем
Аг Pi) Ъ 1, 2, . . . , 72,
(<р,,ф1)= 0, i ф j.
Две квадратные матрицы А и В одинаковых размеров называются по ными, если существует такая невырожденная матрица Р, что А = Р~1 Подобные матрицы имеют одни и те же собственные значения, так ка] (6.1) непосредственно следует
Вф = \ф, ф = Р^р.
Поэтому вычислительные алгоритмы решения спектральных задач баз] ются на подобном преобразовании матрицы к матрице В, для которой с тральная задача решается проще. В качестве В естественно выбирать ди нальную матрицу, причем в данном случае это будет
A = diag{AbA2,...,An}.
Упорядочим собственные значения симметричной матрицы А по возра нию, т.е. Ai < Л2 < • • • < Лп. Свойства собственных значений и собствен
(Ах х)
функций связаны с отношением Релея ~Отметим, например, что
(ж, х)
любого ненулевого вектора х справедливы неравенства
(х,х)
В
Do'stlaringiz bilan baham: |