d5r2.cpp

C语言的编程合集！是我自己编的！很有用的！大家可以下载！有意见麻烦大家提出来以便本人及时修改！

#include "iostream.h"
#include "math.h"

double gammln(double xx)
{
	double cof[7];
	cof[1] = 76.18009173;
	cof[2] = -86.50532033;
	cof[3] = 24.01409822;
	cof[4] = -1.231739516;
	cof[5] = 0.00120858003;
	cof[6] = -0.00000536382;
	double stp = 2.50662827465;
	double half = 0.5;
	double one = 1.0;
	double x,tmp,ser,fpf = 5.5;
	x = xx - one;
	tmp = x + fpf;
	tmp = (x + half) * log(tmp) - tmp;
	ser = one;
	for (int j = 1; j<=6; j++)
	{
		x = x + one;
		ser = ser + cof[j] / x;
	}
	return tmp + log(stp * ser);
}

double factrl(int n)
{
    double a[34];
    int ntop = 0;
    a[1] = 1.0;
    if (n < 0)
	{
        cout<< "negative factorial"<<endl;
        return -1;
	}
    else
	{
		if (n <= ntop)
		{
			return a[n + 1];
		}
		else
		{
			if (n <= 32)
			{
				for (int j = ntop + 1; j<=n; j++)
				{
					a[j + 1] = j * a[j];
				}
				ntop = n;
				return a[n + 1];
			}
		    else
			{
				return exp(gammln(n + 1.0));
			}
		}
	}
}

void ddpoly(double c[], int nc, double x, double pd[], int nd)
{
	int i,j,nnd;
    pd[1] = c[nc];
    for (j = 2; j<=nd; j++)
	{
        pd[j] = 0.0;
    }
    for (i = nc - 1; i>=1; i--)
	{
        if (nc + 1 - i < nd)
		{
           nnd = nc + 1 - i;
		}
        else
		{
           nnd = nd;
        }
        for (j = nnd; j>=2; j--)
		{
            pd[j] = pd[j] * x + pd[j - 1];
        }
        pd[1] = pd[1] * x + c[i];
    }
    double const1 = 2.0;
    for (i = 3; i<=nd; i++)
	{
        pd[i] = const1 * pd[i];
        const1 = const1 * i;
    }
}

void main()
{
    //program d5r2
    //driver for routine ddpoly
    //polynomial (x-1)^5
    char *a[6];
    int i,j,nc = 6;
    int ncm1 = 5;
    int np = 20;
	double aaa,x;
    double c[7], pd[6], d[6][21];
    a[1] = "polynomial:";
    a[2] = "first deriv:";
    a[3] = "second deriv:";
    a[4] = "third deriv:";
    a[5] = "fourth deriv:";
    c[1] = -1.0;
    c[2] = 5.0;
    c[3] = -10.0;
    c[4] = 10.0;
    c[5] = -5.0;
    c[6] = 1.0;
    for (i = 1; i<=np; i++)
	{
        x = 0.1 * i;
        ddpoly(c, nc, x, pd, nc - 1);
        for (j = 1; j<=nc - 1; j++)
		{
            d[j][i] = pd[j];
        }
    }
	cout.setf(ios::fixed|ios::right);
	cout.precision(6);
    for (i = 1; i<=nc-1; i++)
	{
        cout<<a[i]<<endl;
        cout<<"          x          ddpoly       actual"<<endl;
        for (j = 1; j<=np; j++)
		{
            x = 0.1 * j;
			cout.width(14);
            cout<<x;
			cout.width(14);
            cout<<d[i][j];
            aaa = factrl(nc - 1) / factrl(nc - i) * (pow((x - 1.0) , (nc - i)));
			cout.width(14);
            cout<<aaa<<endl;
        }
        cout<<endl;
    }
}