Размещено на http://www.allbest.ru/

Министерство образования и науки Российской Федерации

Государственное образовательное учреждение

Высшего профессионального образования

"ОРЕНБУРГСКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ"

Факультет экономики и управления

Кафедра математических методов и моделей в экономике

Лабораторная работа

По курсу численные методы

по теме "Решение системы линейных алгебраических уравнений методом Зейделя"

Оренбург 2011

1. Постановка задачи

Найти корни системы линейных алгебраических уравнений, используя метод Зейделя.

2. Краткие теоретические сведения

Система линейных алгебраических уравнений (СЛАУ) имеет вид Ax=b, где

, ,

СЛАУ имеет решение, если ранг матрицы А равен n (числу переменных) и если матрица А не является вырожденной.

Метод Зейделя является модификацией метода простых итераций.

Суть метода простых итераций состоит в выборе начального приближения (обычно за начальное приближение берут столбец свободных членов) и последующим нахождением элементов вектора xⁱ⁺¹ по формуле . Модификация Зейделя ускоряет сходимость за счёт того, что при расчёте элементов xⁱ⁺¹ используются уже вычисленные элементы этого вектора, т.е. . Если взять за начальное приближение вектор свободных членов, то получаем следующую итерационную последовательность:



Однако, метод Зейделя сходится далеко не для любой произвольной матрицы. Для сходимости метода Зейделя необходимо либо диагональное преобладание в матрице А, либо нормальность системы Ах=b (т.е., А является симметричной и положительно определённой матрицей). В случае, если не выполняется ни одно из этих условий, а матрица А является невырожденной, то следует умножить матрицу и вектор b на матрицу А^Т слева. Тогда мы получим нормальную систему, которая может быть решена при любом начальном приближении х⁰. Перед оценкой погрешности необходимо провести некоторые преобразования. Пусть B=-D^-1(L+R), где D состоит из элементов главной диагонали матрицы А, L - из элементов, лежащих ниже и левее главной диагонали, R - из элементов, лежащих выше и правее главной диагонали. Тогда матрица B будет иметь вид , а матрица B_R будет состоять из элементов матрицы B, расположенных выше главной диагонали. Оценка погрешности производится следующим образом: . Исходя из этого, получаем - критерий останова.

3. Алгоритм

0) А, b, е, х⁰=b, i:=0.

1) Находим матрицу B=-D^-1(L+R), считаем .

2) Если , то к шагу 4. Иначе - к шагу 3.

3) Если , то А=А^ТА, b=A^Tb, к шагу 4). Иначе метод Зейделя не обладает сходимостью для данной матрицы.

4) .

5) Если , то на конец алгоритма, х^*=хⁱ. Иначе i:=i+1, к шагу 4.

4. Блок-схема

 

5. Код программы

unit Unit1;

interface

uses

Windows, Messages, SysUtils, Variants, Classes, Graphics, Controls, Forms,

Dialogs, XPMan, StdCtrls, Grids, ComCtrls;

type

TForm1 = class(TForm)

strngrd1: TStringGrid;

btn1: TButton;

strngrd2: TStringGrid;

lbl1: TLabel;

xpmnfst1: TXPManifest;

edt1: TEdit;

ud1: TUpDown;

lbl2: TLabel;

mmo1: TMemo;

procedure edt1Change(Sender: TObject);

procedure btn1Click(Sender: TObject);

private

{ Private declarations }

public

{ Public declarations }

end;

const ep=0.0001;

var n:integer;

eps:real;

x0,x,b: array [1..5] of real;

a: array [1..5,1..5] of Real;

Form1: TForm1;

implementation

{$R *.dfm}

procedure xnext(num:integer);

var j:integer;

t:Real;

begin

for j:=1 to n do x0[j]:=x[j];

t:=0;

for j:=1 to (num-1) do t:=t+a[num, j]*x[j];

for j:=(num+1) to n do t:=t+a[num, j]*x[j];

x[num]:=-(t-b[num])/a[num,num];

end;

function dxnorm:real;

var i:integer;

dx:real;

begin

dx:=0;

for i:=1 to n do dx:=dx+sqr(x[i]-x0[i]);

dxnorm:=Sqrt(dx);

end;

function shod:Boolean;

var i,j:integer;

t:Real;

begin

shod:=True;

for i:=1 to n do begin

t:=0;

for j:=1 to (i-1) do

t:=t+a[i,j];

for j:=(i+1) to n do

t:=t+a[i,j];

if a[i,i]<t then begin

shod:=False;

Break;

end;

end;

end;

procedure normalize;

var i,j,k:integer;

bnormed: array [1..5] of Real;

normed: array [1..5,1..5] of Real;

begin

for i:=1 to n do

for j:=1 to n do

normed[i,j]:=0;

for i:=1 to n do begin

bnormed[i]:=0;

for j:=1 to n do bnormed[i]:=bnormed[i]+b[j]*a[j,i];

end;

for i:=1 to n do b[i]:=bnormed[i];

for i:=1 to n do

for j:=1 to n do

for k:=1 to n do normed[i,j]:=normed[i,j]+a[k,i]*a[k,j];

for i:=1 to n do

for j:=1 to n do

a[i,j]:=normed[i,j];

end;

procedure ocenkapogr;

var i,j:integer;

bb: array [1..5, 1..5] of real;

bnorm,brnorm: Real;

procedure bfind;

var i,j:integer;

begin

for i:=1 to n do bb[i,i]:=0;

for i:=1 to n do

for j:=1 to n do

if i<>j then bb[i,j]:=-a[i,j]/a[i,i];

end;

begin

bfind;

for i:=1 to n do

for j:=(i+1) to n do brnorm:=brnorm+sqr(bb[i,j]);

bnorm:=brnorm;

for i:=1 to n do

for j:=1 to (i-1) do bnorm:=bnorm+sqr(bb[i,j]);

brnorm:=Sqrt(brnorm);

bnorm:=Sqrt(bnorm);

eps:=ep*Abs((1-bnorm)/brnorm);

end;

procedure TForm1.edt1Change(Sender: TObject);

begin

n:=StrToInt(Edt1.Text);

strngrd1.ColCount:=n;

strngrd1.RowCount:=n;

strngrd2.RowCount:=n;

end;

procedure TForm1.btn1Click(Sender: TObject);

var i,j:integer;

begin

for i:=0 to (n-1) do begin

for j:=0 to (n-1) do

a[i+1,j+1]:=StrToFloat(strngrd1.Cells[j,i]);

b[i+1]:=StrToFloat(strngrd2.Cells[0,i]);

end;

if shod=False then normalize;

ocenkapogr;

for i:=1 to n do x[i]:=b[i];

for i:=1 to n do xnext(i);

while dxnorm>eps do

for i:=1 to n do xnext(i);

mmo1.Text:='';

for i:=1 to n do

mmo1.Text:=mmo1.Text+'x'+IntTostr(i)+'='+FloatTostr(x[i])+'; ';

end;

end.

6. Тестовый пример и проверка в MathCad

зейдель линейный алгебраический уравнение

Размещено на Allbest.ru