d13r1.cpp

使用VC++编写的大量数学算法的源代码

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

int round(double x)
{
	int temp;
	double temp1;
	temp=int(x);
	temp1=x-temp;
	if(temp1>0.5)
	{
		return temp=temp+1;
	}
	else
	{
		return temp;
	}
}

void moment(double data[], int n, double& ave, double& adev, double& sdev, double& var, double& skew, double& curt)
{
	int j;
    if (n <= 1)
	{
		cout<<" 'n must be at least 2'";
		return;
    }
    double p,s = 0.0;
    for (j = 1; j<=n; j++)
	{
        s = s + data[j];
    }
    ave = s / n;
    adev = 0;
    var = 0;
    skew = 0;
    curt = 0;
    for (j = 1; j<=n; j++)
	{
        s = data[j] - ave;
        adev = adev + fabs(s);
        p = s * s;
        var = var + p;
        p = p * s;
        skew = skew + p;
        p = p * s;
        curt = curt + p;
    }
    adev = adev / n;
    var = var / (n - 1);
    sdev = sqrt(var);
    if (var != 0.0) 
	{
        skew = skew / (n * pow(sdev , 3));
        curt = curt / (n * var * var) - 3.0;
	}
    else
	{
        cout<<" no skew or kurtosis when zero variance"<<endl;
    }
}

void main()
{
    //program d13r1
    //driver for routine moment
    const double pi = 3.14159265;
    int npts = 10000;
    int nbin = 100;
    int ndat = npts + nbin;
    double data[100101];
    int nlim,j,k,i = 1;
	double x;
    for (j = 1; j<=nbin; j++)
	{
		double temp=0;
        x = pi * j / nbin;
        nlim = round(sin(x) * pi / 2.0 * (double)npts / (double)nbin);
        for (k = 1; k<=nlim; k++)
		{
            data[i] = x;
            i = i + 1;
        }
    }
	data[npts]=0.0;
	cout.setf(ios::fixed|ios::left);
	cout.precision(4);
    cout<<endl;
    cout<<"Moments of a sinusoidal distribution"<<endl;
	double adev,sdev,var,ave,skew,curt;
    moment(data, npts, ave, adev, sdev, var, skew, curt);
    cout<<"                     "<<"Calculated     Expected"<<endl;
    cout<<"Mean :                 ";
	cout.width(14);	   cout<<ave;
	cout.width(14);    cout<<pi / 2.0<<endl;
    cout<<"Average Deviation :    ";
	cout.width(14);    cout<<adev;
	cout.width(14);    cout<<0.5708<<endl;
    cout<<"Standard deviation :   ";
	cout.width(14);    cout<<sdev;
	cout.width(14);    cout<<0.6837<<endl;
    cout<<"Variance :             ";
	cout.width(14);    cout<<var;
	cout.width(14);    cout<<0.4674<<endl;
    cout<<"Skewness :             ";
	cout.width(14);    cout<<skew;
	cout.width(14);    cout<<0<<endl;
    cout<<"Kurtosis :             ";
	cout.width(14);    cout<<curt;
	cout.width(14);    cout<<-0.8062;
	cout<<endl;
}