d2r7.cpp

Visual C++ 常用数值算法集 源代码

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

void polcoe(double x[], double y[], int n, double cof[])
{
    double phi,ff,b,s[15];
	int i,j,k;
    for(i = 1; i<=n; i++)
	{
		s[i] = 0.0;
		cof[i] = 0.0;
	}
    s[n] = -x[1];
    for (i = 2; i<=n; i++)
	{
		for (j = n + 1 - i;j<=n-1;j++)
		{
			s[j] = s[j] - x[i] * s[j + 1];
		}
		s[n] = s[n] - x[i];
	}
    for (j = 1; j<=n; j++)
	{
		phi = n;
		for (k = n - 1; k>=1; k--)
		{
			phi = k * s[k + 1] + x[j] * phi;
		}
		ff = y[j] / phi;
		b = 1.0;
		for (k = n; k>=1; k--)
		{
			cof[k] = cof[k] + b * ff;
			b = s[k] + x[j] * b;
		}
	}
}

void main()
{
    //program d2r7
    //driver for routine polcoe
    int nfunc,i,j,np = 5;
    double x,f,sum,xa[6], ya[6], coeff[6],pi = 3.141593;
    for (nfunc = 1; nfunc<=2; nfunc++)
	{
        if (nfunc == 1)
		{
            cout<<endl;
            cout<< "Sine function from 0 to pi"<<endl;
            for (i = 1; i<=np; i++)
			{
                xa[i] = i * pi /double( np);
                ya[i] = sin(xa[i]);
            }
		}
        else 
		{
			if (nfunc ==2)
			{
            cout<<endl;
            cout<< "Exponential function from 0 to 1"<<endl;
            for (i = 1; i<=np; i++)
				{
					xa[i] = 1.0 * i / double(np);
					ya[i] = exp(xa[i]);
				}
			}
            else
            return;
        }
        polcoe(xa, ya, np, coeff);
        cout<<"coefficients"<<endl;
        for (i = 1; i<=np; i++)
            cout<<setw(12)<<coeff[i];
		cout<<endl;
        cout<< "         x            f[x]       polynomial"<<endl;
		cout<<setprecision(6)<<setiosflags(ios::fixed);
        for (i = 1; i<=10; i++) 
		{
            if (nfunc == 1 )
			{
                x = (-0.05 + i / 10.0) * pi;
                f = sin(x);
			}
            else 
			{
				if (nfunc == 2)
				{
					x = (-0.05 + double(i) / 10.0);
					f = exp(x);
				}
            }
            sum = coeff[np];
            for (j = np - 1; j>=1; j--)
				{
					sum = coeff[j] + sum * x;
				}
            cout<<setw(14)<<x;
            cout<<setw(14)<<f;
            cout<<setw(14)<<sum<<endl;
		}
        cout<<  "**********************************"<<endl;
    }
}