Исходный код

Теория.

Курсовая работа

С.И.Самыгин

ИСТОРИЯ

Ответственный

за выпуск:

Редактор:

Корректор:

Художник:

Верстка:

Кузнецов В. Федоров И. Тарасенко Н. Лойкова И. Патулова А.

Сдано в набор 3.10.2006 г. Подписано в печать 25.10.2006 г.

Формат 84x108 ¹/₃₂. Бумага типографская.

Гарнитура Школьная.

Тираж 5 000. Заказ № 606

Издательство «Феникс»

344082, г. Ростов-на-Дону,

пер. Халтуринский, 80

Отпечатано с готовых диапозитивов в ЗАО «Книга».

344019, г. Ростов-на-Дону, ул. Советская, 57. Качество печати соответствует предоставленным диапозитивам.

«Сравнение точных методов»

по дисциплине
«Вычислительная математика»

студентки группы УВА-211

Булгаковой Ирины Владимировны

Проверил:

проф. каф. ПМ-1

__________

Волосов Константин Александрович

Москва 2013

1) Метод Гаусса.

#include "stdafx.h"

#include "conio.h"

#include "iostream"

#include "stdio.h"

#include "clocale"

#include "math.h"

#include "complex"

using namespace std;

const int fact=10;

complex <float> normA(0,0);

complex <float> normA1(0,0);

complex <float> x[fact];

short n;

int N=0;

complex <float> Det2(complex <float> r[fact][fact])

{

complex <float> det2=r[0][0]*r[1][1]-r[0][1]*r[1][0];

N++;

return det2;

}

complex <float> Det(complex <float> h[fact][fact],short m)

{

complex <float> det(0,0);

complex <float> s[fact][fact];

complex <float> M(0,0);

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

{

int l,g;

l=0; g=0;

for (int d=0;d<m;d++)

for(int k=0;k<m;k++)

if ((d==i)||(k==0)) {}

else

{

s[l][g]=h[d][k]; //Формируем минор по данному элементу

if (g<m-2) g++;

else

{

g=0;

l++;

}

}

if ((m-1)==2) M=Det2(s);

else M=Det(s,m-1);

short st=pow(-1.00,(double)i+2);

complex <float> q(st,0);

det=det+q*h[i][0]*M;

}

return det;

}

complex <float> Obratnaya(complex <float> ob[fact][fact], short m, complex <float> D)

{

complex <float> s[fact][fact];

complex <float> c[fact][fact];

complex <float> Dg;

int l,g;

if (m==2)

{

c[0][0]=ob[1][1]; c[0][1]=-ob[0][1];

c[1][0]=-ob[1][0]; c[1][1]=ob[0][0];

}

else

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

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

{

l=0; g=0;

for (int d=0;d<m;d++)

for(int k=0;k<m;k++)

if ((d==i)||(k==j)) {}

else

{

s[l][g]=ob[d][k]; //Формируем минор по данному элементу

if (g<m-2) g++;

else

{

g=0;

l++;

}

}

if ((m-1)==2) Dg=Det2(s);

else Dg=Det(s,m-1);

short st=pow(-1.00,(double)i+j+2);

complex <float> q(st,0);

c[j][i]=q*Dg; //Записываем посчитанный элемент в транспонированную матрицу

}

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

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

c[i][j]=c[i][j]/D; //Делим каждый элемент транспонированной матрицы на определитель

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

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

normA1=normA1+pow(c[i][j],2);

normA1=sqrt(normA1); //Вычисляем норму обратной матрицы

return c[fact][fact];

}

complex <float> Gauss(complex <float> t[fact][fact])

{

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

{

if (t[k][k].real()!=1)

for (int i=n;i>=k;i--)

{

t[k][i]=t[k][i]/t[k][k]; // По диагонали делаем единицы

N++;

}

for (int i=k+1;i<n;i++)

if (t[i][k].real()!=0)

for (int j=n;j>=k;j--)

{

t[i][j]=t[i][k]*t[k][j]-t[i][j]; // Делаем нули под каждой единицей на диагонали

N++;

}

}

//Обратный ход Гаусса

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

{

complex <float> q(t[i][n].real(),t[i][n].imag());

x[i]=q; //Задаем начальные значения для решений системы

N++;

}

for (int i=n-2;i>=0;i--)

for (int j=i+1;j<n;j++)

{

x[i]=x[i]-t[i][j]*x[j]; //Вычисляем решения

N++;

}

return x[fact];

}

void main()

{

setlocale(LC_CTYPE,"");

puts("Введите размерность матрицы <10 ");

cin>>n;

complex <float> a[fact][fact];

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

{

float re,im;

cout<<"Помните, что в строке "<<n+1<<" элемента(ов),\n включая свободные члены";

cout<<"\nВведите "<<i+1<<"-ю строку матрицы\n";

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

{

printf("a[%d][%d]:",i,j);

printf(" Real="); scanf("%f",&re);

printf("Imagine="); scanf("%f",&im);

complex <float> q(re,im);

a[i][j]=q;

}

}

cout<<"Вывод матрицы\n";

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

{

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

printf(" [%d][%d]:%5.2f%5.2f",i,j,a[i][j].real(),a[i][j].imag());

cout<<"\n";

}

//Вычисляем норму матрицы

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

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

normA=normA+pow(a[i][j],2.0);

normA=sqrt(normA);

complex <float> z[fact][fact];

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

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

{

z[i][j]=a[i][j];

}

// Считаем определитель

complex <float> D;

if (n==2) D=Det2(a);

else D=Det(a,n);

if (D.real()==0 && D.imag()==0) printf("Система вырождена и не имеет решений, определитель равен %d",D);

else

{

printf("D: Re=%5.2f Im=%5.2f\n",D.real(),D.imag());

Obratnaya(a,n,D);

Gauss(z);

cout<<"Вывод решений\n";

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

printf(" x[%d]:Re=%8.3fIm=%8.3f\n",i,x[i].real(),x[i].imag());

cout<<"\n";

//Находим норму Х

complex <float> normX(0,0);

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

normX=normX+pow(x[i],2);

normX=sqrt(normX);

complex <float> xst[fact];

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

xst[i]=x[i];

//Находим норму В

complex <float> normB(0,0);

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

normB=normB+pow(a[i][n],2);

normB=sqrt(normB);

printf("Количество операций %d\n",N);

printf("Введите погрешность ");

float re;

scanf("%f",&re);

complex <float> delta(re,0);

printf("delta=%6.3f %6.3f\n",delta.real(),delta.imag());

//Находим норму дельтa b

complex <float> delb(0,0);

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

delb=delb+pow(delta,2);

delb=sqrt(delb);

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

a[i][n]=a[i][n]+delta;

Gauss(a);

cout<<"Вывод решений с погрешностью\n";

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

printf(" x[%d]:Re=%8.3fIm=%8.3f\n",i,x[i].real(),x[i].imag());

cout<<"\n";

printf("Норма матрицы=%5.2f %5.2f\n",normA.real(),normA.imag());

printf("Норма обратной матрицы=%5.2f %5.2f\n",normA1.real(), normA1.imag());

complex <float> condA;

condA=normA*normA1; //Число обусловленности

printf("Число обусловлености=%5.2f %5.2f\n",condA.real(),condA.imag());

//Находим норму дельта х

complex <float> Delx[fact];

complex <float> normDelx(0,0);

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

{

Delx[i]=x[i]-xst[i];

normDelx=normDelx+pow(Delx[i],2);

}

normDelx=sqrt(normDelx);

printf("Норма дельта x=%8.6f %8.6f\n",normDelx.real(),normDelx.imag());

printf("Норма дельта b=%8.6f %8.6\n",delb.real(),delb.imag());

printf("Норма b=%8.6f %8.6f\n",normB.real(),normB.imag());

printf("Нормa X=%8.6f %8.6f\n",normX.real(),normX.imag());

complex <float> left(0,0);

left=normDelx/normX;

complex <float> right(0,0);

right=condA*delb/normB;

printf("Апостериорная ошибка= %8.6f %8.6f\n",left.real(),left.imag());

printf("Априорная ошибка= %8.6f %8.6f\n",right.real(),right.imag());

}

getch();

}

2) Метод Крамера.

#include "stdafx.h"

#include "conio.h"

#include "iostream"

#include "stdio.h"

#include "clocale"

#include "math.h"

#include "complex"

using namespace std;

const int fact=10;

complex <float> normA(0,0);

complex <float> normA1(0,0);

complex <float> b[fact];

complex <float> answ[fact];

short int n;

int N=0;

complex <float> Det2(complex <float> r[fact][fact])

{

complex <float> det2=r[0][0]*r[1][1]-r[0][1]*r[1][0];

N++;

return det2;

}

complex <float> Det(complex <float> a[fact][fact],short m)

{

complex <float> det(0,0);

complex <float> s[fact][fact];

complex <float> M;

int l,g;

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

{

l=0; g=0;

for (int d=0;d<m;d++)

for(int k=0;k<m;k++)

if ((d==i)||(k==0)) {}

else

{

s[l][g]=a[d][k]; //Формируем минор по данному элементу

N++;

if (g<m-2) g++;

else

{

g=0;

l++;

}

N++;

}

if ((m-1)==2) M=Det2(s);

else M=Det(s,m-1);

short st=pow(-1.00,(double)i+2);

complex <float> q(st,0);

det=det+q*a[i][0]*M;

}

return det;

}

complex <float> Obratnaya(complex <float> ob[fact][fact], short m, complex <float> D)

{

complex <float> s[fact][fact];

complex <float> c[fact][fact];

complex <float> Dg;

int l,g;

if (m==2)

{

c[0][0]=ob[1][1]; c[0][1]=-ob[0][1];

c[1][0]=-ob[1][0]; c[1][1]=ob[0][0];

}

else

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

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

{

l=0; g=0;

for (int d=0;d<m;d++)

for(int k=0;k<m;k++)

if ((d==i)||(k==j)) {}

else

{

s[l][g]=ob[d][k]; //Формируем минор по данному элементу

if (g<m-2) g++;

else

{

g=0;

l++;

}

}

if ((m-1)==2) Dg=Det2(s);

else Dg=Det(s,m-1);

short st=pow(-1.00,(double)i+j+2);

complex <float> q(st,0);

c[j][i]=q*Dg; //Записываем посчитанный элемент в транспонированную матрицу

}

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

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

c[i][j]=c[i][j]/D; //Делим каждый элемент транспонированной матрицы на определитель

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

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

normA1=normA1+pow(c[i][j],2);

normA1=sqrt(normA1); //Вычисляем норму обратной матрицы

return c[fact][fact];

}

complex <float> Kramer(complex <float> a[fact][fact],complex <float> b[fact],complex <float> D)

{

complex <float> detxy[fact][fact];

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

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

detxy[i][j]=a[i][j]; //Формируем матрицу для подсчета определителя каждой переменной

N++;

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

{

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

{

detxy[i][j]=b[i]; //Замещаем один столбец матрицы на свободные члены

N++;

}

if (n==2) answ[j]=Det2(detxy)/D;

else answ[j]=Det(detxy,n)/D; // Находим значения переменных

N++; N++;

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

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

detxy[i][j]=a[i][j]; //Возвращаем матрице значение исходной матрицы

}

return answ[fact];

}

void main()

{

setlocale(LC_CTYPE,"");

complex <float> a[fact][fact];

puts("Введите размерность матрицы <=10 ");

cin>>n;

float re,im;

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

{

cout<<"Помните, что в строке "<<n<<" элемента(ов)\n";

cout<<"Введите "<<i+1<<"-ю строку матрицы\n";

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

{

printf("a[%d][%d]:",i,j);

printf(" Real="); scanf("%f",&re);

printf("Imagine="); scanf("%f",&im);

complex <float> q(re,im);

a[i][j]=q;

}

}

cout<<"Введите свободные члены\n";

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

{

printf("b[%d]:",i);

printf(" Real="); scanf("%f",&re);

printf("Imagine="); scanf("%f",&im);

complex <float> q(re,im);

b[i]=q;

}

cout<<"Вывод матрицы\n";

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

{

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

printf(" a[%d][%d]:%8.6f %8.6f",i,j,a[i][j].real(),a[i][j].imag());

cout<<"\n";

}

cout<<"Вывод свободных членов\n";

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

printf("b[%d]:%8.6f %8.6f ",i,b[i].real(),b[i].imag());

cout<<"\n";

//Вычисляем норму матрицы

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

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

normA=normA+pow(a[i][j],2.0);

normA=sqrt(normA);

// Считаем определитель

complex <float> D;

if (n==2) D=Det2(a);

else D=Det(a,n);

if (D.real()==0 && D.imag()==0) printf("Система не имеет решений, определитель равен %d",D);

else

{

printf("D:%8.6f %8.6f\n",D.real(),D.imag());

Obratnaya(a,n,D);

Kramer(a,b,D);

printf("Вывод решений\n");

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

printf("answ[%d]:%8.6f %8.6f ",i,answ[i].real(),answ[i].imag());

cout<<"\n";

//Находим норму Х

complex <float> normX(0,0);

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

normX=normX+pow(answ[i],2);

normX=sqrt(normX);

complex <float> answst[fact];

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

answst[i]=answ[i];

//Находим норму В

complex <float> normB(0,0);

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

normB=normB+pow(b[i],2);

normB=sqrt(normB);

printf("Количество операций %d\n",N);

printf("Введите погрешность ");

float re;

scanf("%f",&re);

complex <float> delta(re,0);

printf("delta=%6.3f %6.3f\n",delta.real(),delta.imag());

//Находим норму дельтa b

complex <float> delb(0,0);

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

delb=delb+pow(delta,2);

delb=sqrt(delb);

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

b[i]=b[i]+delta;

Kramer(a,b,D);

cout<<"Вывод решений с погрешностью\n";

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

printf(" answ[%d]:Re=%8.6fIm=%8.6f\n",i,answ[i].real(),answ[i].imag());

cout<<"\n";

printf("Норма матрицы=%8.6f %8.6f\n",normA.real(),normA.imag());

printf("Норма обратной матрицы=%8.6f %8.6f\n",normA1.real(), normA1.imag());

complex <float> condA;

condA=normA*normA1; //Число обусловленности

printf("Число обусловлености=%8.6f %8.6f\n",condA.real(),condA.imag());

//Находим норму дельта х

complex <float> Delx[fact];

complex <float> normDelx(0,0);

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

{

Delx[i]=answ[i]-answst[i];

normDelx=normDelx+pow(Delx[i],2);

}

normDelx=sqrt(normDelx);

printf("Норма дельта x=%8.6f %8.6f\n",normDelx.real(),normDelx.imag());

printf("Норма дельта b=%8.6f %8.6\n",delb.real(),delb.imag());

printf("Норма b=%8.6f %8.6f\n",normB.real(),normB.imag());

printf("Нормa X=%8.6f %8.6f\n",normX.real(),normX.imag());

complex <float> left(0,0);

left=normDelx/normX;

complex <float> right(0,0);

right=condA*delb/normB;

printf("Апостериорная ошибка= %8.6f %8.6f\n",left.real(),left.imag());

printf("Априорная ошибка= %8.6f %8.6f\n",right.real(),right.imag());

}

getch();

}

Дата добавления: 2014-11-29