Числові методи

Контрольная работа - Компьютеры, программирование

Другие контрольные работы по предмету Компьютеры, программирование

°триці А ;

d - вектор вільних членів;

n- степінь матриці А;

Вихідні дані:

М- вектор моментів кубічного сплайну.*/

{float k[nMax],fi[nMax];

int i;

// прямий хід

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

{k[i] = (i==0)? -c[i]/b[i] : -c[i]/(b[i]+a[i]*k[i-1]);

fi[i] = (i==0)? d[i]/b[i] : (-a[i]*fi[i-1]+d[i])/(b[i]+a[i]*k[i-1]);}

//зворотній хід

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

M[i] = (i==n)? fi[i-1] : k[i-1]*M[i+1]+fi[i-1];}

void fSplain( int n,float h,float Y[nMax+1],float M[nMax+1],float Ak[nMax][4])

/* Функція обчислює коефіцієнти кубічного сплайну

Вхідні дані:

n- кількість відрізків розбиття;

H - крок розбиття відрізку [X0; Xn]]

М- вектор моментів кубічного сплайну.

Y- вектор значень функції f(x,y) в точках x[0],x[1],...x[n].

Вихідні дані:

Ak- матриця коефіцієнтів кубічного сплайну.*/

{int i;

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

{Ak[i][0] = Y[i];

Ak[i][1] = (Y[i+1]-Y[i])/h-h/6*(2.*M[i]+M[i+1]);

Ak[i][2] = M[i]/2;

Ak[i][3] = (M[i+1]-M[i])/6*h;}}

void main()

{float Y[nMax+1],d[nMax],M[nMax+1],Ak[nMax][4];

int n,i,j;

n=nMax;

M[0]=0; M[n]=0; //крайові умови

char *strError="\n Error of file !";

FILE *FileIn,*FileOut,*FileOut2;

FileIn=fopen("data2_2in.txt","r"); // відкриваємо файл для читання

if (FileIn==NULL)

{ cout << " \"Data2_2in.txt\": Error open file or file not found !!!\n";

goto exit; }

FileOut=fopen("data2_2ou.txt","w"); // відкриваємо файл для запису

if (FileOut==NULL)

{ cout << " \"Data2_2ou.txt\": Error open file !!!\n";

goto exit; }

FileOut2=fopen("data2_3in.txt","w"); // відкриваємо файл для запису

if (FileOut2==NULL)

{ cout << " \"Data2_3in.txt\": Error open file !!!\n";

goto exit; }

// читаємо вектор Y

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

if(fscanf(FileIn,"%f,",&(Y[i]))==NULL)

{ cout << strError; goto exit;};

// обчислюємо вектор d

for (i=1; i<n; i++) d[i-1]=3/(h*h)*(Y[i+1]-2*Y[i]+Y[i-1]);

//fMetodProgonku(n-1,a,b,c,d,M);

fSplain( n,h,Y,M,Ak);

// Вивід результатів в тому числі і для наступного завдання

fprintf(FileOut2,"%d\n",n); // n - кількість відрізків

// координати точок сітки по Х

for(float xi=x0,i=0; i<n; i++) fprintf(FileOut2,"%2.2f ",xi+h*i);

fprintf(FileOut2,"\n");

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

{for (j=0; j<4; j++)

{printf("a[%d,%d]= %4.4f ",i,j,Ak[i][j]);

fprintf(FileOut,"a[%d,%d]= %4.4f ",i,j,Ak[i][j]);

fprintf(FileOut2,"%4.4f ",Ak[i][j]);}

cout << endl;

fprintf(FileOut,"\n");

fprintf(FileOut2,"\n");}

fclose(FileIn);

fclose(FileOut);

exit: cout << "\n Press any key ...";

getch();}

Результат роботи програми (" data2_2uo.txt"):

a[0,0]= 1.0000 a[0,1]= 0.5104 a[0,2]= 0.0000 a[0,3]= -0.0104

a[1,0]= 1.0500 a[1,1]= 0.4793 a[1,2]= -0.3107 a[1,3]= 0.0118

a[2,0]= 1.0960 a[2,1]= 0.4525 a[2,2]= 0.0429 a[2,3]= -0.0068

a[3,0]= 1.1410 a[3,1]= 0.4407 a[3,2]= -0.1607 a[3,3]= 0.0054

в) Розібємо відрізок на частин.

Складова формула Сімпсона буде мати вигляд:

 

;

 

де - крок розбиття, значення функції в точках сітки.

Замінимо значеннями кубічних сплайнів із пункту б) цього завдання.

 

 

Для оцінки похибки використаємо правило Рунге. Для цього обчислимо наближені значення інтегралу з кроком (), а потім з кроком ().

За наближене значення інтегралу, обчисленого за формулою Сімпсона з поправкою по Рунге приймемо: .

//------------------------------------------------------------

// Work2_3.cpp

//------------------------------------------------------------

// "Числові методи"

// Завдання 2

// Обчислення інтегралу методом Сімпсона з використанням кубічного сплайну

#include

#include

#include

#include

// визначення сплайнового класу

class Tsplain

{public:

int kol; // кількість рівнянь (відрізків розбиття)

float ** Ak; // масив коефіцієнтів

float * Xi; // вектор початків відрізків

float vol(float x); // функція повертає значення сплайну в точці х

Tsplain(int k); // constructor};

Tsplain::Tsplain(int k)

{kol=k;

Xi=new float[kol];

Ak=new float*[kol];

for(int i=0; i<kol; i++) Ak[i]=new float[kol];}

float Tsplain::vol(float x)

{float s=0.;

int i,t;

// шукаємо відрізок t де знаходиться точка х

for (i=0; i=Xi[i]) { t=i; break; }

s=Ak[t][0];

for (i=1; i<kol; i++)

s+=Ak[t][i]*pow((x-Xi[t]),i);

return s;}

float fSimps(float down,float up, int n, Tsplain *spl)

/* Функція обчислює інтеграл методом Сімпсона з використанням кубічного сплайну

Вхідні дані:

down,up -границі інтегрування ;

n- число відрізків , на яке розбиваєтьться відрізок інтегрування ;

spl - вказівник на обєкт класу Tsplain ( кубічний сплайн );

Вихідні дані:

функція повертає знайдене значення інтегралу.*/

{float s=0;

float h=(up-down)/(float)n;

int i;

s=spl->vol(down)+spl->vol(up-h);

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

s+=2*(spl->vol(down+i*h));

for (i=1; i<n; i+=2)

s+=4*(spl->vol(down+i*h));

return s*h;}

void main()

{int kol; // кількість рівняннь кубічного сплайну

float down,up;

float I1,I2,I,eps;

int n,i,j;

char *strError="\n Error of file !";

FILE *FileIn,*FileOut;

FileIn=fopen("data2_3in.txt","r"); // відкриваємо файл для читання

if (FileIn==NULL)

{ cout << " \"Data2_3in.txt\": Error open file or file not found !!!\n";

goto exit; }

FileOut=fopen("data2_3ou.txt","w"); // відкриваємо файл для запису

if (FileOut==NULL)

{ cout << " \"Data2_3ou.txt\": Error open file !!!\n";

goto exit; }

// читаємо kol

if(fscanf(FileIn,"%d,",&kol)==NULL)

{ cout << strError; goto exit;};

Tsplain *sp;

sp=new Tsplain(kol);

// читаємо вектор Xi

for(i=0; iXi[i]));

// читаємо масив Ak

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

for (j=0; jAk[i][j]));

// читаємо n - кількість відрізків розбиття відрізку інтегрування

if(fscanf(FileIn,"%d,",&n)==NULL)

{ cout << strError; goto exit;};

down=sp->Xi[0];

up=sp->Xi[sp->kol-1]+(sp->Xi[sp->kol-1]-sp->Xi[sp->kol-2]);

I1=fSimps(down,up, n, sp);

I2=fSimps(down,up, 2*n, sp);

eps=(I2-I1)/15;

I=I2+eps;

// Вивід результатів

printf("I= %5.5f\n",I);

printf("EPS= %5.5f\n",eps);

fprintf(FileOut,"I= %5.5f\n",I);

fprintf(FileOut,"EPS= %5.5f\n",eps);

fclose(FileIn);

fclose(FileOut);

exit: cout << "\n Press any key ...";

getch();}

Результат роботи програми ("data2_3ou.txt")

I= 1.32213

EPS= 0.00004

 

Завдання 3

 

Знайти розвязок системи нелінійних рівнянь

 

,

 

Рішення.

Умову завдання