newtonfor.cpp

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

// NewtonFor.cpp: implementation of the NewtonFor class.

#include "NewtonFor.h"
#include "iostream.h"

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

NewtonFor::NewtonFor()
{
 NewtonTable=new NFTable;
 NewtonTable->builde=false;
 NewtonTable->n=NULL;
 NewtonTable->thead=NULL;
 h=NULL;
 RecNum=NULL;
 x0=NULL;
 xstart=NULL;
 y0=NULL;
 xarray=NULL;
}

NewtonFor::~NewtonFor()
{
} 

void NewtonFor::SetData()
{
 //输入初始点数据:
 //RecNum,xstart(表中x列第一个数),h(步长),xarray,yi列所有数据
 cout<<"please enter the points number:"<<endl;
 cin>>RecNum;
 xarray=new float[RecNum];
 NewtonTable->n=RecNum;

 NewtonTable->thead=new NFNode[RecNum];
 NFNode* p=NewtonTable->thead;
 //NTNode* p=new NTNode[RecNum];
 //p=NewtonTable->thead;

 cout<<"please enter the points position below"<<endl;
 cout<<"please enter the origin x and h"<<endl;
 cin>>xstart>>h;
 for(int k=0;k<RecNum;k++)
 {
  xarray[k]=xstart+k*h;
 }
 for(int i=0;i<RecNum;i++)
 {
  cout<<"please enter point_"<<i<<": Y"<<endl;
  cin>>p->value;
  p->rhand=NULL;
  p->next=NULL;
  if(i<RecNum-1) 
  {
   p->next=NewtonTable->thead+i+1;
   p=p->next;
  }
 }
}

void NewtonFor::SetData(float x[],float y[],int n)
{
 //:外部接口函数
 RecNum=n;                   //sizeof(x);
 xstart=x[0];
 if(RecNum) h=x[1]-x[0];
 xarray=new float[RecNum];
 NewtonTable->thead=new NFNode[RecNum];
 NFNode* p=NewtonTable->thead;
 //xarray=x;
 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;
  }
 }
}

void NewtonFor::AddData()
{
 //添加数据结点
}

void NewtonFor::CreateTable()
{
 //生成牛顿插值列表
 NFNode* p=NewtonTable->thead;//已知差分列指针
 int i,k;
 for(i=1;i<RecNum;i++)
 {
  NFNode* Lanse=new NFNode[RecNum-i];
  for(k=0;k<RecNum-i;k++)
  {
   //:给新节点赋值
   Lanse->value=p->next->value-p->value;
   Lanse->rhand=NULL;
   Lanse->next=NULL;
   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;

  //delete[] Lanse;
  //Lanse=NULL;
 }
 NewtonTable->builde=true;
}

void NewtonFor::Sety0(float Curx)
{
 //取得差分链头节点
 y0=NewtonTable->thead;
 x0=xarray[0];
 if(NewtonTable->builde==true) 
 for(int i=1;i<RecNum;i++)
 {
  if(Curx>xarray[i])  
  {
   y0=y0->next;
   x0=xarray[i];
  }
  else break;
 }
}

float NewtonFor::Cal(float Curx)
{
 //计算x=Curx处近似值
 if(!RecNum) SetData();
 if(!NewtonTable->builde) CreateTable();
 if(!y0) Sety0(Curx);
 float t;
 int Lanse;
 float Rita=y0->value;
 NFNode* p=NULL;
 int i=0;
 for(p=y0->rhand;p;p=p->rhand)
 {
  t=(Curx-x0)/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 NewtonFor::Cal_Error(float Curx)
{
	float *xarray_temp=new float[RecNum+1];
	float *yarray_temp=new float[RecNum+1];
	float e=0;
    int i;
	NFNode* 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;
}