d14r8.cpp

这个是常用的数值算法用VC编写的。相信对大家有用哦。

#include<math.h>
#include<iomanip.h>
#include<iostream.h>
#include<process.h>

void pzextr(int iest, double xest,double yest[], double yz[],double dy[],int nv,int nuse)
{
    static double x[12], qcol[11][ 8];
	int j,ncol,nmax,m1,k1;
	double delta,f1,f2,q;
    ncol = 7;
    nmax = 10;
    double d[11];
    x[iest] = xest;
    for (j = 1; j<=nv; j++)
	{
        dy[j] = yest[j];
        yz[j] = yest[j];
    }
    if (iest == 1)
	{
        for (j = 1; j<=nv; j++)
		{
            qcol[j][1] = yest[j];
        }
	}
    else
	{
        m1 = iest;
        if (iest > nuse) m1 = nuse;
        for (j = 1; j<=nv; j++)
		{
            d[j] = yest[j];
		}
        for (k1 = 1; k1<=m1 - 1; k1++)
		{
            delta = 1.0 / (x[iest - k1] - xest);
            f1 = xest * delta;
            f2 = x[iest - k1] * delta;
            for (j = 1; j<=nv; j++)
			{
                q = qcol[j][k1];
                qcol[j][k1] = dy[j];
                delta = d[j] - q;
                dy[j] = f1 * delta;
                d[j] = f2 * delta;
                yz[j] = yz[j] + dy[j];
            }
        }
        for (j = 1; j<=nv; j++)
		{
            qcol[j][m1] = dy[j];
        }
    }
    for (j = 1; j<=10; j++)
	{
    d[j]=0.0;
	}
}

void main()
{
    //program d14r8
    //driver for routine pzextr
    //feed values from a rational function
    //fn(x)=(1-x+x**3)/(x+1)**n
	int nv,nuse,i,j,iest;
	double xest,dum;
    nv = 4;
    nuse = 5;
    double  yest[5], yz[5], dy[5];
    for (i = 1; i<=10; i++)
	{
        iest = i;
        xest = 1.0 / float(i);
        dum = 1.0 - xest + pow(xest , 3);
        for (j = 1; j<= nv; j++)
		{
            dum = dum / (xest + 1.0);
            yest[j] = dum;
        }
        pzextr(iest, xest, yest, yz, dy, nv, nuse);
        cout<<endl;
        cout<<"I= "<< i<<endl;
        cout<<"Extrap. function: ";
		cout<<setprecision(6)<<setiosflags(ios::fixed);
        for (j = 1; j<= nv; j++)
		{
            cout<<setw(13)<<yz[j];
        }
		cout<<endl;
        cout<<"Estimated error:  ";
        for (j = 1; j<= nv; j++)
		{
            cout<<setw(13)<<dy[j];
        }
		cout<<endl;
    }
	cout<<endl;
    cout<<"actual values:    "<<setw(13)<<"1.000000";
    cout<<setw(13)<<"1.000000"<<setw(13)<<"1.000000";
	cout<<setw(13)<<"1.000000"<<endl;
}