Если у вас возникли сложности с курсовой, контрольной, дипломной, рефератом, отчетом по практике, научно-исследовательской и любой другой работой - мы готовы помочь.
Предоплата всего
Подписываем
Лабораторная работа №2
Система m линейных алгебраических ура в не ний с n не известными в общем виде может быть за писана сле ду ющим образом:
(1.22)
или в матричном виде AX = B, где A - прямоу голь ная мат рица размерности mn, X - век тор n-го по ряд ка, B - век тор m-го порядка. Ре ше нием сис темы (1.22) на зы ва ет ся такая упо ря до ченная со во куп ность чисел которая об ра ща ет все уравнения сис те мы в вер ные ра вен ства. Две системы называются эк ви ва лентными (рав но силь ны ми), если множества их ре ше ний сов падают.
Система линейных уравнений называется сов мест ной, ес ли она имеет хотя бы одно ре ше ние, и не сов мест ной - в про тивном случае. Сов мест ная сис те ма на зы ва ет ся определенной, если она имеет един ст вен ное ре ше ние, и не оп ре де лен ной - в про тив ном слу чае. Система яв ля ется опре де ленной, ес ли rang A = rang B, где матри ца B, по лученная из матрицы A добавлением стол бца сво бод ных чле нов, на зы вается рас ши рен ной.
Если матрица A - квадратная и det A ¹ 0, то она на зы вается неособенной (невырожденной), при этом сис те ма уравнений, имеющая не о со бен ную мат ри цу A, сов мес тна и имеет единственное ре шение.
МЕТОД ИСКЛЮЧЕНИЯ ГАУССА ДЛЯ СИСТЕМ ЛИНЕЙНЫХ УРАВНЕНИЙ
Из курса линейной алгебры [Крылов и др., 1972; Ку рош, 1962; Фадеев, Фадеева, 1963 и др.] из вестно, что ре ше ние системы линейных урав не ний можно прос то най ти по правилу Крамера - че рез отношение оп ре де ли те лей. Но этот способ не очень удобен для ре шения сис тем урав не ний с чис лом неизвестных > 5, т.е. когда най ти опре де ли тель сложно, а при чис ле не известных > 10 на хож де ние оп ре де ли теля с до ста точ но высокой сте пенью точ нос ти ста но вится са мо сто ятельной вы чис ли тель ной за дачей. В этих слу чаях применяют иные методы ре шения, среди ко то рых самым рас про стра нен ным яв ля ет ся метод Га ус са.
Запишем систему линейных уравнений (1.22) в ви де
(1.24)
Если матрица системы верхняя треугольная, т.е. ее эле мен ты ниже главной диагонали рав ны ну лю, то все хj можно найти по сле до ва тель но, на чиная с хn, по фор му ле
. (1.25)
При j > k и аjj ¹ 0 этот метод дает возможность най ти решение системы.
Метод Гаусса для произвольной системы (1.22) ос но ван на приведении матрицы А сис те мы к верх ней или ниж ней треугольной. Для это го выч и та ют из второго уравнения системы первое, умно женное на такое число, при котором а21 = 0, т.е. ко эф фи ци ент при х1 во второй строке должен быть ра вен ну лю. Затем та ким же образом ис клю чают по сле до ва тельно а31 , а41 , ..., аm1 . После за вер ше ния вы чис ле ний все эле мен ты первого столбца, кроме а11, бу дут равны ну лю.
Продолжая этот процесс, исключают из мат ри цы А все коэффициенты аij, лежащие ниже главной ди а го на ли.
Построим общие формулы этого процесса. Пусть ис клю чены коэффициенты из k - 1 столб ца. Тог да ниже глав ной диагонали долж ны остаться урав не ния с нену ле выми эле ментами:
Умножим k-ю строку на число
и вычтем ее из m-й строки. Первый не ну левой элемент этой строки обратится в нуль, а ос таль ные изменятся по формулам
(1.26)
Произведя вычисления по формулам (1.26) при всех ука занных индексах, исключим эле мен ты k-го столб ца. Та кое исключение назовем цик лом, а вы полнение всех цик лов назовем пря мым ходом ис ключения.
После выполнения всех циклов получим верх нюю тре угольную матрицу А системы (1.24), ко то рая легко ре шается обратным ходом по фор му лам (1.25). Если при прямом ходе коэффициент аjj ока зал ся ра вен нулю, то пе ре ста нов кой строк пе ремеща ем на глав ную ди а го наль не ну ле вой эле мент и про дол жа ем расчет.
Если предположить, что аjj ¹ 0, то тогда можно при ме нить подпрограмму gaus, ре ша ю щую сис те му из n ли ней ных уравнений.
Формальные параметры процедуры. Входные: n (тип in te ger) - порядок системы; а (тип real) - мас сив ко эф фи ци ентов системы размером nn; b (тип real) - массив-стро ка свободных чле нов. Выход ные: х (тип real) - мас сив-строка, в который по ме ща ется ре шение системы.
procedure gaus (n:integer; a : mas; b : mas1;
var x : mas1);
type mst = array [1..n] of real;
mss = array [1..n] of mst;
var a1 : mss; b1 : mst;
i, j, m, k : integer; h : real;
begin
for i := 1 to n do
begin
b1[i] := b[i];
for j := 1 to n do
a1 [i,j] := a[i,j];
end;
for i := 1 to n-1 do
for j := i+1 to n do
begin
a1[j,i] :=- a1[j,i] / a1[i,i];
for k := i+1 to n do
a1[j,k] := a1 [j,k] + a1[j,i]*a1[i,k];
b1[j] := b1[j] + b1[i]*a1[j,i];
end;
x[n] := b1[n] / a1[n,n];
for i := n-1 downto 1 do
begin
h := b1[i];
for j := i+1 to n do
h := h - x[j]*a1[i,j];
x[i] := h / a1[i,i];
end;
end.
Если контроль за невырожденностью мат ри цы А не об ходим, то можно предложить воспользоваться другой про це ду рой gaus1.
В отличие от первой процедуры здесь в вы ход ных па раметрах есть переменная tol, воз вра ща ю щая 0 при нор мальном завершении работы про це ду ры, или 1, ес ли на главной диагонали один из эле ментов ра вен 0, или 2, если матрица А раз мер ностью боль ше, чем 5050.
procedure gaus1 (n:integer; a : mas; b : mas1;
var x : mas1; var tol : integer);
type mst = array [1..50] of real;
mss = array [1..50] of mst;
var a1 : mss; b1 : mst;
i, j, m, k : integer; h : real;
begin
for i := 1 to n do
begin
b1[i] := b[i];
for j := 1 to n do a1 [i,j] := a[i,j];
end;
tol := 0;
if n > 50 then
begin
tol := 2;
exit;
end;
for i := 1 to n-1 do
for j := i+1 to n do
begin
if abs(a1[i,i])< 1.0e-10 then
begin
tol := 1;
exit;
end;
a1[j,i] :=- a1[j,i] / a1[i,i];
for k := i+1 to n do
a1[j,k] := a1 [j,k] + a1[j,i]*a1[i,k];
b1[j] := b1[j] + b1[i]*a1[j,i]
end;
if abs(a1[n,n]) < 1.0e-10 then
begin
tol := 1;
exit;
end;
x[n] := b1[n] / a1[n,n];
for i := n-1 downto 1 do
begin
if abs(a1[i,i]) < 1.0e-10 then
begin
tol := 1;
exit;
end;
h := b1[i];
for j := i+1 to n do h := h - x[j]*a1[i,j];
x[i] := h / a1[i,i];
end;
end.
По поводу приведенных ал го рит мов можно ска зать, что они не самые эф фек тив ные и при год ны для ре ше ния систем, имеющих невырожденные матрицы А и В с рангом не бо лее 50, составленные из при бли зи тельно оди на ко вых ко эф фи ци ентов. Если меж ду на и большим и на и мень шим из ко эф фи ци ен тов рас стоя ние более чем 104, то эти ал го ритмы при ме нять не рекомендуется, так как по греш ность вы чис ле ний бу дет та кова, что может ли бо привести к вы рож ден ной мат ри це А, ли бо дать в ре зультате век тор Х с боль шой вы чис ли тель ной по греш нос тью.
Для проверки работы процедур ре шим сис тему ли ней ных уравнений методом Гаусса с точ нос тью до 0.0001:
Коэффициенты, промежуточные результаты и око н чательное решение приводятся в таб л. 1.14. Подобные таб лицы удобно ис поль зо вать для руч ного счета и про вер ки хода ре ше ни я.
Таблица 1.14
|
Коэффициенты при неизвестных |
Cвободные |
Cтрочные |
|||
|
x1 |
x2 |
x3 |
x4 |
члены |
суммы |
|
0.68 0.21 0.11 -0.08 |
0.55 -0.13 0.84 0.15 |
-0.11 0.27 0.28 0.50 |
0.88 0.80 -0.06 0.12 |
2.15 0.44 0.83 1.16 |
2.85 -0.01 1.44 0.61 |
|
1 |
0.0735 |
-0.1618 |
0.1176 |
3.1618 |
4.1912 |
|
-0.1454 -0.18319 0.15590 |
0.30398 0.26220 -0.51290 |
-0.8247 0.0729 -0.1106 |
-0.22398 -0.48220 1.41290 |
-0.88015 -0.97847 0.94530 |
|
|
1 |
-2.09060 |
5.6719 |
1.54040 |
6.1217 |
|
|
-1.47643 -0.18697 |
4.79139 -0.99480 |
0.79920 1.17230 |
4.1136 -0.0095 |
||
|
1 |
-3.24410 |
-0.54110 |
-2.7851 |
||
|
-1.60130 |
1.07110 |
-0.5302 |
|||
|
1 |
-0.66890 |
0.3311 |
|||
|
2.8264 |
-0.3337 |
-2.7110 |
-0.6689 |
: РЕШЕНИЕ |
На практике применяют множество ва ри ан тов вы чис ли тельных схем метода Гаусса. Например, при при ве дении мат рицы к верхней тре у голь ной фор ме вы би ра ют на и боль ший элемент (в строке или стол бце), умень шая вы чис ли тель ную по греш ность за счет де ле ния на не са мый ма лень кий эле мент. Та кая вы чис ли тель ная схема на зы ва ет ся ме то дом Га ус са с выбо ром ве ду щего элемента. Ес ли же вы би рать при при ведении мат рицы са мый боль шой (по мо ду лю) элемент из всех ос тав ших ся, то такая схе ма бу дет на зы ваться ме то дом Га усса с вы бором глав но го эле мента. По след няя схема от но сит ся к на и бо лее по пулярным. Глав ное ее от личие от метода Гаусса, рас смот рен ного в п. 3.1, сос то ит в том, что при при ведении мат рицы А к верх ней (или ниж ней) тре угольной фор ме ее стро ки и стол б цы пе ре став ля ют так, чтобы на и боль ший из всех ос тав ших ся эле мен тов матрицы стал ве дущим, и на него вы пол ня ет ся де ление. Ес ли мат ри ца хоро шо обу с лов лена, то в ме тоде Га усса с вы бо ром глав но го элемента по греш нос ти ок руг ле ния не ве ли ки.
Формальные параметры процедуры. Входные: n (тип in te ger) - целое положительное число, рав ное по рядку ис ходной системы; aa (тип real) - массив из nn дейс т ви тель ных чисел, содержащий мат ри цу ко эф фи циентов сис те мы (aa[1] = а11; aa[2] = а21; aa[3] = a31; ...; aa[n]= an1; aa[n +1]= =a12; aa[n + 2] = a22; ...; aa[2*n] = аn2; ...; aa[n*n] = аnn); b (тип real) - мас сив из n действительных чисел, со дер жа щий стол бец сво бодных членов исходной системы (b[1] = =b1; b[2] = b2; ...; b[n] = bn). Выходные: b (тип real) - мас сив из n действительных чисел (он же был вход ным) при выходе из подпрограммы, со дер жа щий ре ше ния ис ход ной системы (b[1]= x1; b[2] = x2; ...; b[n] = xn); ks (тип in te ger) - признак пра виль нос ти ре ше ния или код ошиб ки: если ks = =0, то в ма с си ве b со дер жится ре шение ис ход ной сис те мы; если ks = 1, то исход ная сис тема не име ет ре ше ния (глав ный опре де ли тель ее равен нулю).
procedure sinq (var a: mas11; var b : mas1;
var ks : integer; n : integer);
var tol, biga, aa, save : real;
jj,jy,ij,it,j,i, i1, imax, k, iqs, ixj, ixjx : integer;
jjx, ny, ia, ib, io, i2, ix, jx, ky : integer;
label Cont10;
begin tol := 0.0;
ks := 0;
jj := -n;
for j := 1 to n do
begin jy := j + 1;
jj := jj + n + 1;
biga := 0.0;
it := jj - j;
for i := j to N do
begin
ij := it + i;
aa := abs (a[ij]);
if (abs(biga)-aa) < 0.0 then
begin
biga := a[ij];
imax := i; end;
end;
if (abs(biga)-tol) <= 0.0 then
begin
ks := 1;
exit;
end;
i1 := j + n*(j-2);
it := imax - j;
for k := j to n do
begin
inc (i1,n);
I2:= i1 + it;
save := a[i1];
a[i1] := a[i2];
a[i2] := save;
a[i1] := a[i1] / biga;
end;
save := b[imax];
b[imax] := b[j];
b[j] := save / biga;
if j=n then goto Cont10;
iqs := n*(j-1);
for ix := jy to n do
begin
ixj := iqs + ix;
it := j - ix;
for jx := jy to n do
begin
ixjx := n*(jx-1) + ix;
jjx := ixjx + it;
a[ixjx] := a[ixjx] - a[ixj]*a[jjx];
end;
b[ix] := b[ix] - b[j]*a[ixj];
end;
end;
Cont10: ny := n - 1;
it := n*n;
i := 1;
j :=1;
for ky := 1 to ny do
begin ia := it - ky;
ib := n - ky;
io := n;
for k := 1 to ky do
begin
b[ib] := b[ib] - a[ia]*b[io];
Dec (ia,n);
Dec (io);
end;
end;
end.
Для проверки процедуры и сравнения с ра бо той про цедур в п. 3.1 решалась система урав не ний, при ве ден ная в качестве примера в п. 3.1. Вы числения про во ди лись с точностью до 10-5. Ре зуль таты работы про грам мы приведены далее.
|
Исходная матрица А |
Матрица В |
|
0.680000 0.050000 -0.110000 0.080000 0.210000 -0.130000 0.270000 -0.800000 -0.110000 -0.840000 0.280000 0.060000 -0.080000 0.150000 -0.500000 -0.120000 |
2.150000 0.440000 -0.830000 1.160000 |
|
Преобразованн ая матрица А |
Матрица В |
|
1.000000 0.210000 -0.110000 -0.080000 0.000000 1.000000 -0.145441 0.155882 0.000000 0.000000 1.000000 0.258130 0.000000 0.000000 0.000000 1.000000 |
3.161765 0.579636 -2.851572 -0.669070 |
РЕШЕНИЕ СИСТЕМЫ
2.826351 -0.333733 -2.711759 -0.669070 .
Метод квадратного корня основан на пред став лении матрицы А, составленной из ко эф фи ци ентов сис темы в форме произведения тре у голь ных мат риц, что позволят свести решение за дан ной сис те мы к по сле довательному решению двух сис тем с тре у голь ны ми матрицами.
Метод квадратного корня применяется для ре ше ния сис темы линейных уравнений, ко эф фи ци ен ты ко торой об ра зуют эрмитову сим мет ри чес кую мат ри цу (эр ми то ва мат рица совпадает с ком плек с но-со пря женной тран с по ни ро ванной A* = A).
Представим матрицу А в виде А = S*D S, где S - верх няя треугольная матрица с по ло жи тель ны ми эле мен тами на главной диагонали; S*- транс по ни ро ван ная к матрице S; D- ди а го наль ная мат ри ца, на ди агонали ко торой на хо дят ся чис ла (+1) или (— 1).
Если матрицы S и D найдены, то заданная сис те ма AХ= = F может быть решена следующим пу тем:
AX = S* D S X = (S* D) S X = B Y = F, (1.28)
где S*D = В есть нижняя треугольная матрица и Y = S X - вспомогательный вектор.
Таким об ра зом, ре шение системы АX = F рав но силь но ре ше нию двух треугольных систем ВY = F и SX =Y.
Пусть S = (sik) при i > k и sik ¹ 0, sii > 0; S* = (); D = =(dik), d =1; i ¹ k. Тогда из срав не ния мат риц А и S*DS по лучим
.
Ограничение в сумме получается из учета то го фак та, что в матрице S ниже главной ди а го на ли эле менты об ра ща ются в нуль. Последнее ра венс т во можно пе ре пи сать не сколько иначе, при няв для оп ределенности k = min(k, l):
;
,
откуда окончательно получим формулы для вы чис ления элементов матриц S и D:
(1.29)
Единственным условием возможности опре де ле ния sik является skk ¹ 0. Для по стро е ния мат риц по лагаем k = 1 и последовательно вы чис ля ем все эле менты первой строки s по формуле (1.29); за т ем по лагаем k = 2 - и определяем эле мен ты вто рой стро ки и т.д. Когда k = n, тогда най де ны все эле мен ты матриц S и D, а следовательно и S*. За тем последовательно вы пол ня ем вы чис ления (1.28):
S* Z = F; D Y = Z; S X = Y
обыч ным ходом по формулам
(1.30)
при i = 2, 3, ..., n.
Попутно заметим, что определитель матрицы А можно вычислить из выражения
. (1.31)
Этот метод экономичен, требует n3/3 ариф ме тических действий и при больших n ( n > 50) вдвое бы стрее метода Га усса. Если за основу про цедуры при нять алгоритм П5.6 Дья конова [1987], то тогда метод квад ратного кор ня мо жет быть ре а ли зован c помощью сле дующей про це ду ры:
procedure KVK (M, N : integer;
aa: array of real; bb : array of real;
var C: array of real);
type mas1=array [0..4] of real;
mas=array [0..4] of mas1;
var a : mas; j, k, i : integer; s, c1 : real;
begin
i := 0;
for j := 1 to m do
for k := 1 to n do
begin Inc (i); a[j,k] := aa [i]; end;
for j := 1 to N do
begin
for k := j to N do
begin
s := 0.0;
for i := 1 to M do s := s + a[i,j] * a[i,k];
c[k] := s;
end;
c1 := 0.0;
for i := 1 to M do с1 := c1 + a[i,j]*Bb[i];
for i := j to N do a[i,j] := c[i];
c[J] := c1;
end;
a[1,1] := sqrt (a[1,1]);
for j := 2 to N do a[1,j] := a[j,1] / a[1,1];
for i := 2 to n do
begin
s := 0.0;
for k := 1 to i-1 dо s := s + a[k,i]*a[k,i];
a[i,i] := sqrt (a[i,i] - S);
for j := i+1 to N do
begin
s := 0.0;
for k := 1 to i-1 do
S := S + A[K,I]*A[K,J];
A[I,J] := (A[J,I] - S) / A[I,I];
end;
end;
c[1] := c[1] / a[1,1];
for i := 2 to N do
begin s := 0.0;
for k := 1 to i-1do
s := s + a[k,i] * c[k];
c[i] := (c[i] - s) / a[i,i];
end;
c[n] := c[n] / a[n,n];
for i := n-1 downto 1 do
begin s := 0.0;
for k := i+1 to N do s := s + a[i,k] * c[k];
c[i] := (c[i] - s) / a[i,i];
end;
еnd.
Формальные параметры процедуры.Входные: m, N (тип integer) - m уравнений и n не из вест ных оп ре де ляют раз мер матрицы А. Для этой про це ду ры обя за тельно вы пол не ние m > n, т.е. ко ли чество урав не ний долж но быть боль ше ко личества не из вест ных (сис те ма пе ре о пре делена), пре дель ный раз ре шен ный слу чай m = n; а - мат рица, со став лен ная из ко эф фи ци ен тов при не из вест ных; В - мас сив, со став ленный из стол бца сво бод ных чле нов. Вы ход ные: С - массив, в ко то ром со дер жит ся решение сис темы.
Для проверки правильности работы про це ду ры ре ша лась система 44 линейных уравнений ме то дом квад рат ных корней с точностью 0.0001:
0.68X1 + 0.05X2 + 0.11X3 + 0.08X4 = 2.15 ;
0.05X1 + 0.13X2 + 0.27X3 + 0.80X4 = 0.44 ;
0.11X1 + 0.27X2 + 0.28X3 + 0.06X4 = 0.83 ;
0.08X1 + 0.80X2 + 0.06X3 + 0.12X4 = 1.16 ,
решение которой, полученное методом квад рат ных кор ней, будет следующим
Х 1 = 2.97; Х2 = 1.11; Х3 = 0.74; Х4 = -0.07.
Итерационные методы решения систем ли ней ных урав нений дают возможность вычислить ре ше ние сис те мы как предел бесконечной по сле до ва тель нос ти про ме жуточных решений. Причем каж дое по сле ду ющее ре ше ние в слу чае схо ди мости ите рационного про цесса счи та ет ся более точ ным. В этих методах, в от личие от точ ных (см. п. 3.1 - 3.4), ошиб ки в на чаль ном при бли жении и по сле ду ю щих вы чис ле ниях ком пен сируются, т.е. ите ра ци он ные методы (в слу чае схо ди мости) по зволяют по лу чить решение бо лее точ ное, чем пря мые. По э то му ите ра ци онные ме то ды от носят к са мо ис прав ля ю щим ся.
Условия и скорость сходимости процесса в боль шей сте пени зависят от свойств уравнений, т.е. от свойств мат ри цы системы и от выбора на чаль ного при бли же ния.
Пусть дана система линейных ал геб ра и чес ких урав не ний (1.22) с неособенной матрицей. В ме то де простой ите ра ции если аii ¹ 0, то ис ход ная сис тема мо жет быть пре об ра зо вана к виду хi = bi + aij хj , i ¹ j, т.е. из каждого уравнения по сле до ва тельно вы ра жа ют хi.
Здесь bi = Fi / аii; aij = - аij / аii. Таким образом, в мат рич ном виде имеем Х = В + AХ. Полученную сис тему бу дем решать методом по сле до ва тель ных при ближений. За ну левое приближение Х(0) мож но при нять матрицу В: Х(0)= = B, и далее, под ста вив най денные значения в исходную систему, по лу чим Х (1) = В + A Х(0) .
При бесконечном повторении этой вы чис ли тель ной схемы имеем
,
где и будет искомое решение системы.
Условия сходимости итерационного процесса опре де ля ются теоремами, которые приводятся на ми без до ка за тельств.
Теорема 1. Для того, чтобы по сле до ва тель ность при бли жений Х(n) сходилась, до ста точ но, что бы все соб ст вен ные значения матрицы A были по мо ду лю меньше еди ни цы: | li | < 1, i = 1, 2, ..., n.
Теорема 2. Если требовать, чтобы по сле до ва тель ность Х(n) сходилась к при любом на чаль ном при бли же нии Х(0) , то условие те о ре мы 1 яв ля ет ся не об хо ди мым.
Применение теорем 1 и 2 требует знания всех соб ст вен ных значений матрицы A, нахождение ко то рых яв ля ет ся очень не простой за да чей. По э то му на прак ти ке огра ни чи ва ются бо лее прос той те о ре мой, дающей до статочные ус ло вия схо димости.
Теорема 3. Если для системы Х = В + AХ в ы пол няется хотя бы одно из условий :
;
,
то итерационный процесс сходится к един ствен но му ре ше нию этой системы независимо от вы бо ра на чаль но го при ближения.
Для многих приложений важно знать, какой яв ля ет ся ско рость сходимости про цес са, и оценить по греш ность по лу ченного ре ше ния.
Теорема 4. Если какая-либо норма мат ри цы A, со гла со ванная с рассматриваемой нормой век то ра Х, меньше еди ницы, то верна следующая оцен ка по греш нос ти при бли жения в методе прос той итерации:
.
В библиотеках стандартного ма те ма ти чес ко го обес пе че ния ЭВМ всегда можно найти не сколь ко ва риантов про грам мы, выполняющей ре шение сис те мы линейных урав не ний методом прос той ите ра ции.
Так, можно предложить к использованию сле ду ю щую удобную и достаточно эффективную про це дуру, взятую из БСП БЭСМ для транслятора ALFA. Перевод на язык PASCAL выполнен ав то ра ми.
procedure ITER (n, ik :integer; eps : real;
a : mas1; b : mas; var x : mas);
var x1 : mas; s : real; i, j, k : integer;
begin x1 := b;
x := x1; k := 0;
repeat s := 0.0;
Inc(k);
for i := 1 to n do
begin
for j := 1 to n do x[i] := a[i,j]*x1[j] + b[j];
s := s + abs (x[i]-x1[i]);
end;
s := s / n; x1 := x;
until (s<eps) and (k>ik);
end.
Формальные параметры процедуры. Входные: А (тип real) - матрица, составленная из коэф фи ци ен тов при Х преобразованного уравнения; В (тип real) - мат ри ца, сос тав лен ная из свободных членов; N (тип integer) - раз мер ность мат риц А (N N) и В(N); IK (тип integer) - пре дель но воз мож ное количество итераций (введено для то го, что бы в слу чае рас хож дения процесса выйти из под про грам мы. Обы ч но решение достигается за 3-6 ите ра ций); ЕРS (тип real) - за дан ная погрешность ре шения. Вы ходные: Х (тип real) - матрица, в ко торой на хо дится ре ше ние сис темы.
Для примера методом простых итераций ре шена система 44 линейных уравнений с точ ностью до 0.0001:
Х1 = 0.08Х1 + 0.05Х2 + 0.11Х3 + 0.08Х4 + 2.15
Х2 = 0.05Х1 + 0.13Х2 + 0.27Х3 + 0.28Х4 + 0.44
Х3 = 0.11Х1 + 0.27Х2 + 0.28Х3 + 0.06Х4 + 0.83
Х4 = 0.08Х1 + 0.18Х2 + 0.06Х3 + 0.12Х4 + 1.16
Результаты вычислений представлены в табл. 1.16.
Таблица 1.16
|
X[1] |
X[2] |
X[3] |
X[4] |
№ итерации |
|
2.150000 |
0.440000 |
0.830000 |
1.160000 |
iter = 0 |
|
1.252800 |
1.484800 |
1.229600 |
1.299200 |
iter = 1 |
|
1.263936 |
1.523776 |
1.237952 |
1.315904 |
iter = 2 |
|
1.265272 |
1.528453 |
1.238954 |
1.317908 |
iter = 3 |
|
1.265433 |
1.529014 |
1.239075 |
1.318149 |
iter = 4 |
|
1.265452 |
1.529082 |
1.239089 |
1.318178 |
iter = 5 |
|
1.265454 |
1.529090 |
1.239091 |
1.318181 |
iter = 6 |
|
1.265455 |
1.529091 |
1.239091 |
1.318182 |
iter = 7 |
|
1.265455 1.529091 1.239091 1.318182 :РЕШЕНИЕ |
Этот метод относится к итерационным, имеет бо лее бы струю сходимость по сравнению с ме то дом простых ите раций (см. п. 3.5), но часто при во дит к громоздким вы чис лениям.
Метод Зейделя применяют в двух расчетных схе мах. Рассмотрим каноническую форму (для схе мы ите раций). Пусть система приведена к ви ду Х = Вх + b. В схеме про с той итерации сле ду ю щее при бли же ние Х(k +1) находит ся пу тем под ста нов ки пре ды ду ще го приближения Х(k) в правую часть выражения X(k +1) = B X (k) + b.
Для удобства рассуждений полагаем, что ле вые час ти уравнений содержат хi (элементы матрицы X(k +1)) с воз рас та ю щи ми но ме рами, т.е. сначала х1, затем х2, х3, ..., хn. Тог да ре ше ние системы уравнений Х = Вх + b на оче ред ной (k +1) итерации будет иметь вид
, (1.32)
т.е. каждое очередное найденное приближение хi(k +1) сразу же используется для определения . Ус ло вия схо ди мос ти для итерационного ме тода Зей де ля дают те же тео ре мы, что и в ме то де простых ите ра ций.
Вторая форма метода Зейделя требует пред ва ри тель но го приведения системы (1.22) к ви ду, ког да все ди а го наль ные элементы отличны от ну ля. Ес ли раз ло жить мат ри цу А на сумму двух мат риц R + С, где R - ниж няя ди а го наль ная мат ри ца, а С - верх няя с нулями на главной ди а го на ли, то систему (1.22) можно за пи сать как
Ax = (R + C)x = R x(k +1) + C x(k) = B
или x(k +1) = R -1B - R -1C x(k)
и тог да становится ясно, что метод Зейделя в ка но ни чес кой форме равносилен методу простой ите ра ции, при ме не н ному к системе
X = R -1B - R -1C X.
Для решения системы линейных уравнений ме то дом Зейделя можно использовать процедуру из п. 3.5, слег ка ее модифицировав для слу чая сис те мы n урав не ний:
procedure ITER (n, ik :integer; eps : real;
a : mas1; b : mas; var x : mas);
var x1 : mas; s : real; i, j, k : integer;
begin x1 := b;
x := x1;
k := 0;
repeat s := 0.0;
Inc(k);
for i := 1 to n do
begin
for j := 1 to n do
x[i] := a[i,j]*x1[j] + b[j];
s := s + abs (x[i]-x1[i]);
x1 := x;
end;
s := s / n; x1 := x;
until (s<eps) and (k>ik);
end.
Формальные параметры процедуры такие же, как и в п. 3.5. Для сравнения ре шим такую же сис те му ли ней ных уравнений, как в ме то де ите ра ций (см. п. 3.5). Ре зуль тат решения сис темы ме то дом Зей де ля (табл. 1.17) мож но срав нить с результатами табл. 1.16.
Таблица 1.17
|
x[1] |
x[2] |
x[3] |
x[4] |
№ итерации |
|
2.150000 |
0.440000 |
0.830000 |
1.160000 |
iter = 0 |
|
1.252800 |
1.484800 |
1.229600 |
1.299200 |
iter = 1 |
|
1.263936 |
1.523776 |
1.237952 |
1.315904 |
iter = 2 |
|
1.265272 |
1.528453 |
1.238954 |
1.317908 |
iter = 3 |
|
1.265433 |
1.529014 |
1.239075 |
1.318149 |
iter = 4 |
|
1.265452 |
1.529082 |
1.239089 |
1.318178 |
iter = 5 |
|
1.265452 1.529082 1.239089 1.318178 :РЕШЕНИЕ |
Решить системы методом Гаусса; методом простой итерации и методом Зейделя с точностью , сравнить итерационные методы по числу итераций; сравнить по эффективности (трудность реализации метода, объем памяти, общие затраты времени) итерационные методы и метод Гаусса
а)
|
|
Исходные данные |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
б)
|
|
Исходные данные |
|
|
|
|
|
Исходные данные |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|