Программы с автоматическим и оптимальным выбором шага

// Sys_ode solution by Runge-Kutta with automatic choose of step

#include "stdafx.h"

double A,B,x0,xk,y1x0,y2x0,e;

double a21,c2,b1,b2;

double _x,_y1,_y2,_y11,_y21,fi,h;

void init()

{

A=1.0/20.0;

B=1.0/25.0;

x0=0.0;

xk=3.14;

y1x0=B*3.14;

y2x0=A*3.14;

}

void init_method()

{

e=1.0/19.0;

c2=e;

a21=c2;

b2=1/(2*c2);

b1=1-b2;

}

double Kutta_fun11(double h, double y1,double y2)

{

double temp;

temp=h*A*y2;

return temp;

}

double Kutta_fun21(double h, double y1,double y2)

{

double temp;

temp=-h*B*y1;

return temp;

}

double Kutta_fun12(double h, double y1,double y2)

{

double temp;

temp=h*A*(y2+Kutta_fun21(h,y1,y2));

return temp;

}

double Kutta_fun22(double h, double y1,double y2)

{

double temp;

temp=-h*B*(y1+Kutta_fun11(h,y1,y2));

return temp;

}

double abs_x(double x)

{

if (x<0) return (-x);

else return x;

}

double max_x_y(double x, double y)

{

if (x>y) return x;

else return y;

}

double norm_x_y(double x,double y)

{

return max_x_y(abs_x(x),abs_x(y));

}

void main()

{

init();

printf("x=%f   y1=%f   y2=%f\n",x0,y1x0,y2x0);

init_method();

_x=x0;

_y1=y1x0;

_y2=y2x0;

_y11=y1x0;

_y21=y2x0;

double eps=0.001,eps1=eps/8;

double buff1,buff2,buff11,buff21;

h=0.31415770;

while (_x<=xk)

{

buff1=_y1;

buff2=_y2;

buff11=_y11;

buff21=_y21;

_y1+=b1*Kutta_fun11(h,_y1,_y2)+b2*Kutta_fun12(h,_y1,_y2);

_y2+=b1*Kutta_fun21(h,_y1,_y2)+b2*Kutta_fun22(h,_y1,_y2);

_y11+=b1*Kutta_fun11((h/2),_y11,_y21)+b2*Kutta_fun12((h/2),_y11,_y21);

_y21+=b1*Kutta_fun21((h/2),_y11,_y21)+b2*Kutta_fun22((h/2),_y11,_y21);

_y11+=b1*Kutta_fun11((h/2),_y11,_y21)+b2*Kutta_fun12((h/2),_y11,_y21);

_y21+=b1*Kutta_fun21((h/2),_y11,_y21)+b2*Kutta_fun22((h/2),_y11,_y21);

fi=norm_x_y(_y1-_y11,_y2-_y21)/(1-0.5);

if (fi>eps)

{

printf("%f deJluM nonoJlam\n",h);

h=h/2;

_y1=buff1;

_y2=buff2;

_y11=buff11;

_y21=buff21;

}

if (fi<=eps1)

{

printf("%f yMHo*aeM Ha 2\n",h);

h=h*2;

_y1=buff1;

_y2=buff2;

_y11=buff11;

_y21=buff21;

}

if ((fi<=eps)&&(fi>eps1))

{

_x+=h;

printf("x=%f  kaka  y1=%f  y2=%f\n\n",_x,_y1,_y2);

}

}

}

// Sys_ode solution by Runge-Kutta with static optimal choose of step

#include "stdafx.h"

double A,B,x0,xk,y1x0,y2x0,e;

double a21,c2,b1,b2;

double x,y1,y2,y11,y21,fi,h;

void init()

{

A=1.0/20.0;

B=1.0/25.0;

x0=0.0;

xk=3.14;

y1x0=B*3.14;

y2x0=A*3.14;

}

void init_method()

{

e=1.0/19.0;

c2=e;

a21=c2;

b2=1/(2*c2);

b1=1-b2;

}

double Kutta_fun11(double h, double y1,double y2)

{

double temp;

temp=h*A*y2;

return temp;

}

double Kutta_fun21(double h, double y1,double y2)

{

double temp;

temp=-h*B*y1;

return temp;

}

double Kutta_fun12(double h, double y1,double y2)

{

double temp;

temp=h*A*(y2+Kutta_fun21(h,y1,y2));

return temp;

}

double Kutta_fun22(double h, double y1,double y2)

{

double temp;

temp=-h*B*(y1+Kutta_fun11(h,y1,y2));

return temp;

}

double abs_x(double x)

{

if (x<0) return (-x);

else return x;

}

double max_x_y(double x, double y)

{

if (x>y) return x;

else return y;

}

double norm_x_y(double x,double y)

{

return max_x_y(abs_x(x),abs_x(y));

}

void main()

{

init();

printf("x=%f   y1=%f   y2=%f\n",x0,y1x0,y2x0);

init_method();

x=x0;

y1=y1x0;

y2=y2x0;

y11=y1x0;

y21=y2x0;

h=0.089959852;

double eps=0.0001;

double delta;

while (x<=xk)

{

y1+=b1*Kutta_fun11(h,y1,y2)+b2*Kutta_fun12(h,y1,y2);

y2+=b1*Kutta_fun21(h,y1,y2)+b2*Kutta_fun22(h,y1,y2);

y11+=b1*Kutta_fun11((h/2),y11,y21)+b2*Kutta_fun12((h/2),y11,y21);

y21+=b1*Kutta_fun21((h/2),y11,y21)+b2*Kutta_fun22((h/2),y11,y21);

y11+=b1*Kutta_fun11((h/2),y11,y21)+b2*Kutta_fun12((h/2),y11,y21);

y21+=b1*Kutta_fun21((h/2),y11,y21)+b2*Kutta_fun22((h/2),y11,y21);

fi=norm_x_y(y1-y11,y2-y21)/(1-0.5);

x+=h;

printf("x = %f     y1 = %f     y2 = %f\n",x,y1,y2);

printf("norpewHocmb = %f\n\n",fi);

}

}

Вывод. С автоматическим шагом в Maple получается слишком много итераций цикла, поэтому работать нормально он не может. Однако, в Си, при автоматическом шаге, программа работает, правда времени тратит горазда больше, чем программа с постоянным оптимальным шагом. К сожалению Visual studio не позволяет брать корни n-ой степени, поэтому некоторые результаты приходилось переносить непосредственно из Maple. На графиках, построенных для программ с автоматическим выбором шага и с оптимальным, видно, что оптимальный шаг менее точен.

Сравнив две программы можно однозначно заключить, что наибыстрейшая – с оптимальным постоянным шагом, а наиболее точная – с автоматическим шагом.