adamsout.h

常微分方程数值解法

#ifndef Adamsout_h
#define Adamsout_h
#include "Differential.h"
#define Adamsoutf "adamsoinitial.txt"
template<class Type> class Adamsout:public Differential<Type>{
	private:
		Type b[6][6];
	protected:
		Type *u1,*e1;
	public:
		Adamsout(Type (*pf)(Type,Type)):Differential<Type>
		(pf){u1=new Type[N+1];e1=new Type[N+1];}/**/
		Adamsout(Type (*pu)(Type),Type (*pf)(Type,Type));
		~Adamsout(){delete []u1;delete []e1;}
		void Calculate();
		void Solution();
		void Error();
};
template<class Type> Adamsout<Type>::Adamsout(Type (*pu)(Type),Type (*pf)(Type,Type))
:Differential<Type>(pf){
	int i,j;
	u1=new Type[N+1];e1=new Type[N+1];	u=pu;
	cout<<"Input k:";cin>>k;
    n=N-k+1;
	ifstream fin(Adamsoutf,ios::in|ios::nocreate);
	if(!fin){
		cout<<Adamsoutf<<" cannot be openned!"<<endl;exit(1);
	}
	for(i=0;i<6;i++)
		for(j=0;j<i+1;j++)
		{fin>>b[i][j];}
		for(j=0;j<6;j++){
			b[1][j]=b[1][j]/2;
			b[2][j]=b[2][j]/12;
			b[3][j]=b[3][j]/24;
			b[4][j]=b[4][j]/720;
			b[5][j]=b[5][j]/1440;
		}
		fin.close();
}
template<class Type> void Adamsout<Type>::Calculate(){
	int i,l;Type temp=0;
    cout<<"Input u("<<low<<"):";cin>>u1[0];
	for( i=1;i<=k;i++)
		u1[i]=u1[i-1]+h*F(t[i-1],u1[i-1],f);
	for(i=k+1;i<N+1;i++){
		for(l=0;l<k+1;l++)
		{temp=temp+b[k][l]*F(t[i-l-1],u1[i-l-1],f);
		}//cout<<temp<<endl;
		u1[i]=u1[i-1]+h*temp;//cout<<u1[i]<<endl;
		temp=0;
	}
}
template<class Type> void Adamsout<Type>::Error(){
	Type x=0;
	Exact();
    cout<<setw(6)<<"t"<<setw(14)<<"Adamsout"<<setw(14)<<"Precision"
		<<setw(14)<<"Error"<<endl;
	for(int i=0;i<N+1;i++)
	{
		x=fabs(u0[i]-u1[i]);
		cout<<setw(6)<<t[i]<<setw(14)<<u1[i]<<setw(14)<<u0[i]<<
			setw(14)<<x<<endl;
	}
}
#endif