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

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

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

..........

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

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

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

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.81, 0.07, -0.38, 0.21, //исходные данные

-0.33, 0.41, 0.00, 1.00,

-0.51, 0.07, 0.91, 0.11,

0.22, 0.92, -0.11, -0.33};

float B1[4]={-0.81, -1.21, 1.71, -0.64};

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, тем больше итераций требуется для достижения результата, и тем точнее он получается.