Лабораторная работа № 2 Решение систем линейных уравнений методом простых итераций

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

Дана система линейных уравнений:

	..........

Необходимо решить систему линейных уравнений методом простых итераций.

Итерационный метод – это циклический процесс для получения последовательных приближений к решению.

Для реализации метода необходимо решить следующие задачи:

1) Проверка условия сходимости:

Преобладание диагональных элементов (диагональный элемент должен быть больше (или равен) суммы(е) остальных элементов в строке)

2) Получить итерационную зависимость (нахождение следующего приближения через предыдущее):

	делим элементы строки на
	диагональный элемент
	. . . . . . . . . . со знаком (-)

	. . . . . . . . . .

Иначе систему можно записать так: X_n=R_nxn*X_n+G_n

R_nxn – двумерный массив

X_n и G_n – одномерные массивы

Итерационная зависимость:

3) Получение начального приближения:

1) 0 – одномерный нулевой массив

2) G_n

4) Погрешность:

Eps > 0 допустимая погрешность вычислений

Если Eps достигнуто, расчёты заканчиваются.

Условия завершения итерационного процесса:

разность двух соседних итераций меньше Eps

Исходные данные:

Алгоритм решения:

1) Проверка условия сходимости:

2) Ввод исходных данных:

A_nxn, В_n, Eps

3) Сформировать матрицы R_nxn и G_n:

н. ц. i = 1, n

н. ц. j = 1, n

если i=j, то r_ij=0

иначе

к. ц. j

к. ц. i

4) Задание начального приближения:

н. ц. i = 1, n

(или 0)

к. ц. i

5) Итерационный процесс:

kol=0

do

н. ц. i = 1, n

XK1_i=0

н. ц. j = 1, n

к. ц. j

XK1_i=XK1_i+g_i

к. ц. i

н. ц. i = 1, n

d_i=|XK1_i-XK_i| получение вектора d

к. ц. i

max из d_n

если max > Eps, тогда

н. ц. i=1, n

XK_i= XK1_i

к. ц. i

kol++

while (max > Eps)

6) Вывод результата – XK1, kol (количество итераций)

7) Проверка результата:

XK1 подстав. в исходную систему

Текст программы:

#include<conio.h>

#include<iostream.h>

#include<math.h>

float **A,*B,h,*X,**a,*r,**R,*G,*XK,*XK1,*D; int n;

void init(int q)

{

int i;

if(q)

{

A=new float *[n]; for(i=0;i<n;i++) A[i]=new float[n];

B=new float [n];

X=new float [n]; for(i=0;i<n;i++) X[i]=0;

a=new float *[n]; for(i=0;i<n;i++) a[i]=new float[n];

r=new float [n]; for(i=0;i<n;i++) r[i]=0;

R=new float *[n]; for(i=0;i<n;i++) R[i]=new float[n];

G=new float [n];

XK=new float [n];

XK1=new float [n];

D=new float [n];

}

else

{

for(i=0;i<n;i++) delete[] A[i]; delete[] A;

delete[] B;

delete[] X;

for(i=0;i<n;i++) delete[] a[i]; delete[] a;

delete[] r;

for(i=0;i<n;i++) delete[] R[i]; delete[] R;

delete[] G;

delete[] XK;

delete[] XK1;

delete[] D;

}

}

void cro(void)

{

init(1);

float A1[4][4]={-0.12,1.00,-0.32,0.18, //исходные данные

0.77,0.14,-0.06,0.12,

-0.25,-0.22,-0.14,1.00,

-0.08,0.12,0.77,-0.32};

float B1[4]={-0.72,1.21,0.56,-0.58};

cout<<"\nVvod A[4][4]:\n";

for(int i=0;i<n;i++)

{

for(int j=0;j<n;j++)

{A[i][j]=A1[i][j]; //копируем в А

cout<<A[i][j]<<" ";}

cout<<"\n";}

cout<<"\nVvod B[4]:\n";

for(i=0;i<n;i++)

{

B[i]=B1[i];

cout<<B[i]<<" ";

}

cout<<"\n";

}

float shdm(int t)

{

float tp=0;

for(int j=0;j<n;j++)

if(j!=t) tp+=A[t][j];

return tp;

}

void swap(int a, int b) //меняет местами две строки, когда метод не сходится

{

float tp,rt;

for(int i=0;i<n;i++)

{

tp=A[a][i];

A[a][i]=A[b][i];

A[b][i]=tp;

}

rt=B[a]; B[a]=B[b]; B[b]=rt;

}

float maxd(void)

{

float mx;

mx=D[0];

for(int i=1;i<n;i++) if(D[i]>mx) mx=D[i];

return mx;

}

void main(void)

{

int i,j,k,col=0;

long double max,eps;

textmode(2); clrscr();

cout<<"Vvod n (0)"; gotoxy(9,1); cin>>n; //выводит текущий пример по 0

if(n<=0) или задаём размер и вводим

{ матрицы

n=4;

cro();

}

else

{

init(1);

cout<<"Vvod A["<<n<<"]["<<n<<"]:\n";

for(i=0;i<n;i++)

for(j=0;j<n;j++)

cin>>A[i][j];

cout<<"\n";

cout<<"Vvod B["<<n<<"]:\n";

for(i=0;i<n;i++)

cin>>B[i];}

int *y=new int[3];

y[2]=0;

for(j=0;j<n;j++)

for(i=0;i<n;i++)

if(A[i][i]<shdm(i))

{

y[y[2]]=i; y[2]++;

if(y[2]>=2)

{

swap(y[0],y[1]); y[2]=0;

}

}

delete[] y;

for(i=0;i<n;i++)

for(j=0;j<n;j++)

a[i][j]=A[i][j];

cout<<endl;

for(i=0;i<n;i++)

{

for(j=0;j<n;j++)

{

if(i!=j) R[i][j]=(-A[i][j]/A[i][i]); else R[i][j]=0;

}

G[i]=B[i]/A[i][i];

}

for(i=0;i<n;i++)

{

XK[i]=G[i];

}

cout<<"\nVvod E: "; cin>>eps;

do

{

for(i=0;i<n;i++)

{

XK1[i]=0;

for(j=0;j<n;j++)

{

XK1[i]+=(R[i][j]*XK[j]);

}

XK1[i]=XK1[i]+G[i];

}

for(i=0;i<n;i++)

D[i]=fabs(XK1[i]-XK[i]);

max=maxd();

if(max>=eps)

for(i=0;i<n;i++)

XK[i]=XK1[i];

col++;

}

while(max>=eps);

cout<<"\nIteraziy: "<<col<<"\n";

cout<<"\nResenie:\n";

for(i=0;i<n;i++)

cout<<XK1[i]<<" ";

cout<<"\n";

cout<<"\nProverka:\n";

for(i=0;i<n;i++)

{

for(int j=0;j<n;j++)

r[i]+=(a[i][j]*XK1[j]);

cout<<r[i]<<" ";

}

getch();

init(0);

}

Скриншот результата программы (при Eps=0,01 5 итераций):

Результаты работы программы (проверка говорит о том, что решение правильное):

x₁ = 1, 570428

x₂ = - 0, 721282

x₃ = - 0, 15574

x₄ = 0, 770323

Если сравнить результаты решения систем линейных уравнений методом Гаусса и методом простых итераций, то результаты почти совпадают.

Скриншот результата программы (при Eps=0, 0000001 17 итераций):

Результаты работы программы (проверка говорит о том, что решение правильное):

x₁ = 1, 569947

x₂ = - 0, 720794

x₃ = - 0, 157001

x₄ = 0, 771932

Если сравнить результаты решения систем линейных уравнений методом Гаусса и методом простых итераций, то результаты почти совпадают.

Чем меньше Eps, тем больше итераций требуется для достижения результата, и тем точнее он получается.