four.cpp

这是我数值分析课程的实习题目的源代码(并附有文档),文档包括对问题进行了描述,并对结果进行分析,希望对大家有所帮助.

/*微分方程初值问题
y'=1/(1+x*x)-2*y*y;  0=<x<=4;  y(0)=0;
解析式:y=x/(1+x*x);
*/
#include<iostream.h>
#include<stdio.h>
#include<math.h>
#include<iomanip.h>
double Funtion (double x,double y);
main()
{
	double a=0,b=4,A0=0;
	int n=20;

	//尤拉方法
	double h;
	h=(b-a)/n;
	double x=a,y=A0;
	double You_lai[30],xi[30];
	You_lai[0]=A0;
	xi[0]=a;
	
    for(int i=1;i<=n;i++)
	{
		
		y=y+h*Funtion(x,y);
		x=x+h;
		You_lai[i]=y;
		xi[i]=x;
		
	}
	cout<<endl;
	//龙库方法
	 x=a,y=A0;
	 double Long_ku[30];
     Long_ku[0]=A0;
	 double A_da[30];
     A_da[0]=A0;
	double X[50],Y[50];
	X[0]=a;
	Y[0]=A0;
	for(int j=1;j<=n;j++)
	{
	
		double K1=h*Funtion(x,y);
		double K2=h*Funtion(x+h/2,y+K1/2);
		double K3=h*Funtion(x+h/2,y+K2/2);
		double K4=h*Funtion(x+h,y+K3);
		y=y+(K1+2*K2+2*K3+K4)/6;
		x=x+h;
		Long_ku[j]=y;
		if(j<=3)
		{
			X[j]=x;
			Y[j]=y;
			 A_da[j]=y;
   
		}
	}
	cout<<endl;
	//阿达姆斯方法
     x=a,y=A0;
	for(int k=4;k<=n;k++)
	{
		x=X[k-1]+h;
		y=Y[k-1]+(55*Funtion(X[k-1],Y[k-1])
			-59*Funtion(X[k-2],Y[k-2])
			+37*Funtion(X[k-3],Y[k-3])
			-9*Funtion(X[k-4],Y[k-4]))*h/24;
		A_da[k]=y;
        X[k]= x;
		Y[k]=y;		
	}
	cout<<endl;
	//解析解
	x=a,y=A0;
	double Jie_xi[30];
     Jie_xi[0]=A0;
	for(int g=1;g<=n;g++)
	{
		x=x+h;
		y=x/(1+x*x);
		//y=exp(-x)+x;
		Jie_xi[g]=y;

	}
	cout<<" xi    解析解           尤拉         龙库         阿达姆斯"<<endl;
	for(g=0;g<=n;g++)
	{
		printf("%3.3f  %9.9f   %9.9f  %9.9f    %9.9f"
			,xi[g],Jie_xi[g],You_lai[g],Long_ku[g],A_da[g]);
		cout<<endl;

	}


	return 0;

}
double Funtion (double x,double y)
{
	return (double)1/(1+x*x)-(double)2*y*y;
	//return -y+x+1;
}