Автоматическим контролем шага

Вариант 1

№	Дифференциальное уравнение	y0	x0	xn	Точное решение ДУ

 

1. Метод Рунге - Кутта 4-ого порядка точности

Код программы для диф. уравнения 1-ого порядка:

 

#include<stdio.h>

#include<conio.h>

#include<math.h>

float f(float x,float y)

{

return (y*y -2*x)/(2*x*y);

}

int main()

{

clock_t start = clock();

 

float k1,k2,k3,k4,x,yi,h=0.1;

const float x0=1,xn=2;

float y0=3;

yi=y0;

for(x=x0;x<=xn;)

{

k1=h*f(x,yi);

k2=h*f(x+h/2,yi+k1/2);

k3=h*f(x+h/2,yi+k2/2);

k4=h*f(x+h,yi+k3);

yi=yi+(k1+2*k2+2*k3+k4)/6;

x+=h;

printf("\nx=%f",x);

printf("\nyi=%f",yi);

printf("\nyt=%f",sqrt(9*x -2*x*log(x)));

}

getch();

return 0;

}

 

printf("Work time = %.11Lf", (long double)((clock() -1.0L*start)/CLOCKS_PER_SEC));

system("pause");

return 0;

}

Результаты:

Work time = 0.6710

 

 

 

Метод Рунге - Кутта 4-ого порядка точности с

автоматическим контролем шага

Код программы для диф. уравнения 1-ого порядка:

#include <math.h>

#include <stdio.h>

#include <stdlib.h>

#include <conio.h>

#include <iostream>

#include <time.h>

 

using namespace std;

 

const int N=20;

float func(float x, float y)

{

return y*y-2*x/(2*x*y);

}

 

void auto_method()

{

float x0, xn, h1, eps, h[N], x[N], y[N], k1[N], k2[N], k3[N], k4[N], yy[N];

int i=0;

eps=0.01;

x0=1;

xn=2;

h1=0.1;

h[0]=2*h1;

x[0]=1;

y[0]=3;

for (i=0; i<N-1; i++)

{

do

{

h[i]=h[i]/2;

x[i+1]=x[i]+h[i];

k1[i]=h[i]*func(x[i],y[i]);

k2[i]=h[i]*func(x[i]+h[i]/2,y[i]+k1[i]/2);

k3[i]=h[i]*func(x[i]+h[i]/2,y[i]+k2[i]/2);

k4[i]=h[i]*func(x[i]+h[i],y[i]+k3[i]);

y[i+1]=y[i]+(k1[i]+2*k2[i]+2*k3[i]+k4[i])/6;

}

while(fabs(2*(k1[i]-k2[i]-k3[i]+k4[i])/3)>eps);

if(fabs(2*(k1[i]-k2[i]-k3[i]+k4[i])/3)<(eps/32))

h[i+1]=4*h[i];

else

h[i+1]=2*h[i];

}

for(i=0;i<N-1;i++)

{

yy[i]=sqrt(9*x[i] -2*x[i]*log(x[i])) ;

printf("x=%.4f ", x[i]);

printf("yi=%.4f ", y[i]);

printf("yt=%.4f h=%.4f ", yy[i], h[i]);

printf("|yi-yt|=%.4f ", fabs(y[i]-yy[i]));

printf("\n");

}

}

 

 

int main()

{

clock_t start = clock();

auto_method();

printf("Work time = %.11Lf", (long double)((clock() -1.0L*start)/CLOCKS_PER_SEC));

system("pause");

return 0;

}

Результаты: