newtonback.cpp

数值计算方法中的,牛顿向前和向后算法. 以前花了很长时间编写的,希望对后来的师妹师弟有用.

// NewtonBack.cpp: implementation of the NewtonBack class.
#include "NewtonBack.h"
#include "iostream.h"

extern float f(float x[],float y[],int s,int t);       //差商

NewtonBack::NewtonBack()
{
 NewtonTable=new NBTable;
 NewtonTable->builde=false;
 NewtonTable->n=NULL;
 NewtonTable->thead=NULL;
 h=NULL;
 RecNum=NULL;
 xn=NULL;
 xstart=NULL;
 yn=NULL;
 xarray=NULL;
}
NewtonBack::~NewtonBack()
{
}
void NewtonBack::SetData()
{
 //输入初始点数据
 cout<<"please enter the points number"<<endl;
 cin>>RecNum;
 xarray=new float[RecNum];
 NewtonTable->n=RecNum;
 NewtonTable->thead=new NBNode[RecNum];
 NBNode* p=NewtonTable->thead;
 cout<<"please enter the points position below:"<<endl;
 cout<<"enter origin x0 and step h:"<<endl;
 cin>>xstart>>h;
 //生成差分表xi列
 for(int k=0;k<RecNum;k++)
 //for(int k=RecNum;k>0;k--)
 {
  xarray[k]=xstart-k*h;
 }
 //录入差分表yi列
 for(int i=RecNum;i>0;i--)
 //for(int i=0;i<RecNum;i++)
 {
  cout<<"please enter point_"<<i-1<<": Y"<<endl;
  cin>>p->value;
  //初始化yi链域
  p->rhand=NULL;
  p->next=NULL;
  if(i>1) 
  {
   p->next=NewtonTable->thead+i+1;
   p=p->next;
  }
 }
}
void NewtonBack::SetData(float x[],float y[],int n)
{
 //:外部接口函数
 //:为统一接口,将x[],y[]数组倒置
 RecNum=n;                   //sizeof(x);
 float* Lanse=new float[RecNum];
 for(int k=0;k<RecNum;k++) 
 {
  Lanse[k]=x[RecNum-k-1];
  x[RecNum-k-1]=y[k];
 }
 y=x;
 x=Lanse;
 xstart=x[0];
 if(RecNum) h=x[0]-x[1];
 xarray=new float[RecNum];
 NewtonTable->thead=new NBNode[RecNum];
 NBNode* p=NewtonTable->thead;
 for(int i=0;i<RecNum;i++)
 {
  xarray[i]=x[i];
  p->value=y[i];
  p->rhand=NULL;
  p->next=NULL;
  if(i<RecNum-1) 
  {
   p->next=NewtonTable->thead+i+1;
   p=p->next;
  }
 }
 delete[] Lanse;
 Lanse=NULL;
}
void NewtonBack::CreateTable()
{
 //:生成牛顿插值列表中其他各项
 NBNode* p=NewtonTable->thead;
 int i,k;
 for(i=1;i<RecNum;i++)
 {
  //:生成新列
  NBNode* Lanse=new NBNode[RecNum-i];
  for(k=0;k<RecNum-i;k++)
  {
   //:给新节点赋值
   Lanse->value=p->value-p->next->value;
   Lanse->rhand=NULL;
   Lanse->next=NULL;
   //:连接同列next链域
   if(k<RecNum-i-1) Lanse->next=Lanse+1;
   //:连接新节点到差分表
   p->rhand=Lanse;
   //:指向下一节点
   if(k<RecNum-i-1) Lanse=Lanse->next;
   p=p->next;
  }
  //:指向下一列头节点
  p=NewtonTable->thead;
  for(k=0;k<i;k++) p=p->rhand;
 }
 NewtonTable->builde=true;
}
void NewtonBack::Setyn(float Curx)
{
 //取得差分链头节点
 yn=NewtonTable->thead;
 xn=xarray[0];
 if(NewtonTable->builde==true) 
 for(int i=1;i<RecNum;i++)
 {
  if(Curx<xarray[i])  
  {
   yn=yn->next;
   xn=xarray[i];
  }
  else break;
 }
}
float NewtonBack::Cal(float Curx)
{
 //计算x=Curx处近似值
 if(!RecNum) SetData();
 if(!NewtonTable->builde) CreateTable();
 if(!yn) Setyn(Curx);
 float t;
 int Lanse;
 float Rita=yn->value;
 NBNode* p=NULL;
 int i=0;
 for(p=yn->rhand;p;p=p->rhand)
 {
  t=(Curx-xn)/h;
  Lanse=1;
  for(int k=0;k<i;k++)
  {
   t=t*(t+k+1);
   Lanse=Lanse*(Lanse+k+1);
  }
  Rita=Rita+t*p->value/Lanse;
  i++;
 }
 return Rita;
}

float NewtonBack::Cal_Error(float Curx)
{
	float *xarray_temp=new float[RecNum+1];
	float *yarray_temp=new float[RecNum+1];
	float e=0;
    int i;
	NBNode* p=NewtonTable->thead;
	for(i=0;i<RecNum;i++)
	{
		if(Curx == xarray[i])
			goto end;                   //误差为0
		xarray_temp[i]=xarray[i];
		yarray_temp[i]=p->value;
		if(i<RecNum-1) 
		{
			p->next=NewtonTable->thead+i+1;
			p=p->next;
		}
	}
	xarray_temp[i]=Curx;
	yarray_temp[i]=Cal(Curx);
	e=f(xarray_temp,yarray_temp,0,RecNum);
	for(i=0;i<RecNum;i++)
		e=e*(Curx-xarray_temp[i]);
end:
	delete xarray_temp;
	delete yarray_temp;
	return e;
}