romberg.cpp

数值计算方法中几个重要的算法用VC++实现

//龙贝格算法
#include<iostream>
#include<cmath>
#include<iomanip>
using namespace std;
double f(double x)
{
	if(x<0) return 1;
	return sqrt(x);
}
void main()
{
	int k,i,j,k0,flag=0;
	double a,b,e,sum;
	double **T;
	cout<<"Input a,b,e,k0"<<endl;
	cin>>a>>b>>e>>k0;
	T=new double*[k0+1];
	for(i=0;i<=k0;i++)
		T[i]=new double[i+1];
	T[0][0]=(b-a)/2.0*(f(a)+f(b));
	try{
		for(k=1;k<=k0;k++)
		{    
			sum=0;
			for(i=1;i<=(int)pow(2.0,(double)(k-1));i++)
				sum+=f(a+(2*i-1)*(b-a)/pow(2.0,(double)k));
			T[0][k]=1/2.0*T[0][k-1]+(b-a)/pow(2.0,(double)k)*sum;

			for(j=1;j<=k;j++)
				T[j][k-j]=(pow(4.0,(double)j)*T[j-1][k-j+1]-T[j-1][k-j])/(pow(4.0,(double)j)-1);
			if(fabs(T[k][0]-T[k-1][0])<e){flag=1;break;}
		}
	}
	catch(...)
	{
		cout<<"exception!"<<endl;
	}
	if(flag)
	{
		cout.precision (7);//设置精度
		for(j=0;j<=k;j++)
		{
			for(i=0;i<=j;i++)
				cout<<"T["<<i<<"]["<<(j-i)<<"]="<<T[i][j-i]<<"  ";
		    cout<<endl;
		}
	}
	else 
		cout<<"Failure!"<<endl;

}