romberg1.cpp

本人的计算方法课程作业

//romberg.cpp()--ROMBERG方法计算定积分
#include<iostream>
#include<cmath>
using namespace std;
const int N=20;
float f(float x);
float t(int n);
float s(int n);
float c(int n);
float r(int n);
int p(int n);
static float a,b,T[N],S[N],C[N],R[N];
void main()
{
	float e,r1,r2;
	int i=1;
	cout<<"输入积分下限与积分上限:"<<endl;
	    cin>>a>>b;
	cout<<"输入误差精度:"<<endl;
		cin>>e;
	do
	{
	    r1=r(i);
		r2=r(2*i);
		i*=2;
	}while(fabs(r2-r1)>e);    //若两次积分所得值的差大于误差精度，继续计算
    cout<<"---------------------------------------"<<endl;
	cout<<"k   T   S   C   R"<<endl;
	i=0;
	while(T[i]!=0.0)          //若复化梯形的值不为零，输出结果
	{  
		cout.precision(5);
		cout<<i<<" "<<T[i]<<" "<<S[i]<<" "<<C[i]<<" "<<R[i]<<endl;
		i++;
	}
}
//t()用复化梯形公式计算T[i]
float t(int n)
{
	int k;
	float tn,h,t1;
	if(n==1)
	{h=fabs(b-a);
	tn=(f(a)+f(b))*h/2;}
	else
	{t1=0;h=fabs(b-a)/n;
	for(k=0;k<=n-1;k++)
		t1+=f(a+(2*k+1)*h/2);
	tn=t(n/2)/2+t1*h/2;
	}
	return (T[p(n)]=tn);
}
//s()用复化梯形公式的值计算复化Simpson公式的值S[i]
float s(int n)
{  
	return S[p(n)+1]=t(2*n)+(t(2*n)-t(n))/3.0;
}
//c()用复化Simpson公式的值计算复化Cotes公式的值C[i]
float c(int n)
{   
	return C[p(n)+2]=s(2*n)+(s(2*n)-s(n))/15.0;
}
//r()用复化Cotes公式的值计算Romberg公式的值R[i]
float r(int n)
{
    return R[p(n)+3]=c(2*n)+(c(2*n)-c(n))/63.0;
}
//f()定义被积函数
float f(float x)
{
	float y;
	if(x==0.0) x=0.000001;
	y=(x*x*x+sin(x))/x;
	return y;
}
//p()确定各数组的下标
int p(int n)
{
	int i=0,t=1;
	while(t!=n)
	{t*=2;i++;}
	return i;
}