常微分方程的RK4解法

#include<stdio.h>
#include<stdlib.h>
#include<math.h>
#define H 0.1
int main()
{
	float x,y,r=3;
	float k1,k2,k3,k4;
	float fun1(float,float);
	int i,n;
	printf("请输入x,y的初值：\n");
	scanf("%f%f",&x,&y);
	printf("请输入运算次数：\n");
	scanf("%d",&n);
	printf("运算结果为：\n"); 
    for(i=0;i<n;i++)
	{
		printf("x%-2d=%.2f,\ty%-2d=%f,\ty(%-2d)=%f,\tr(%-2d)=%f",i,x,i,y,i,r,i,r-y);
		printf("\n"); 
	    k1=fun1(x,y);
	    k2=fun1(x+H/2,y+H*k1/2);
	    k3=fun1(x+H/2,y+H*k2/2);
	    k4=fun1(x+H,y+H*k3);
	    x=x+H;
	    y=y+H*(k1+2*k2+2*k3+k4)/6;
	    r=3/(1+x*x*x);
	}
	system("pause");
	return 0;
}
	float fun1(float x,float y)
	{
		float f;
		f=-x*x*y*y;
		return f;
	}