#include<stdio.h>
#include<math.h>
#define  f(x,y)  -(x*x*y*y)//導數
#define  y(x)    2.0/(1.0+x*x*x)//準確解
void main()
{
	int i,j;
	long double h[4]={0.1,0.05,0.025,0.0125};//步長
	int m[4]={15,30,60,120};//不同步長時迭代到1.5的循環次數
	long double x0,y0,xn,yn,yn1,yn2;
	long double k1,k2,k3,k4;
	yn2=(long double)y(1.5);//計算準確值
	x0=0.0;y0=3.0;
    printf("h        數值結果y(1.5)          精確解           數值誤差  \n");
	for(i=0;i<4;i++)
	{   xn=x0;yn=y0;//初始化
		for(j=0;j<m[i];j++)
		{
			k1=f(xn,yn);
			k2=f((xn+h[i]/2.0),(yn+h[i]*k1/2.0));
			k3=f((xn+h[i]/2.0),(yn+h[i]*k2/2.0));
			k4=f((xn+h[i]),(yn+h[i]*k3));
			yn1=yn+h[i]/6.0*(k1+2*k2+2*k3+k4);//龍格-庫塔四階公式
			xn=xn+h[i];
			yn=yn1;
		}
			printf("%.12lf       %.12lf         %.12lf          %.12lf    \n",h[i],yn1,yn2,fabs(yn2-yn1));//輸出
	}
}