d6r6.cpp

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

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

static double r[98];
static int iff,ix1,ix2,ix3;
double ran1(int& idum)
{
    double rm1,rm2,t;
	int m1,m2,m3,ia1,ia2,ia3,ic1,ic2,ic3,j;
	 
    m1 = 259200, ia1 = 7141, ic1 = 54773, rm1 = 0.0000038580247,
    m2 = 134456, ia2 = 8121, ic2 = 28411, rm2 = 0.0000074373773,
    m3 = 243000, ia3= 4561, ic3 = 51349;
	
    if ((idum < 0) || (iff == 0))
	{
        iff = 1;
        ix1 =  (ic1 - idum) % m1;
        ix1 = (ia1 * ix1 + ic1) % m1;
        ix2 = ix1 % m2;
        ix1 = (ia1 * ix1 + ic1) % m1;
        ix3 = ix1 % m3;
        for (j = 1; j<=97; j++)
		{
            ix1 = (ia1 * ix1 + ic1)% m1;
            ix2 = (ia2 * ix2 + ic2)% m2;
            r[j] = (double(ix1) + double(ix2) * rm2) * rm1;
        }
        idum = 1;
	}
      ix1 = (ia1 * ix1 + ic1) % m1;
      ix2 = (ia2 * ix2 + ic2) % m2;
      ix3 = (ia3 * ix3 + ic3) % m3;
      j = 1 + int((97 * ix3) / m3);
      if ((j > 97) || (j < 1)) 
	  {
		  cout<< "abnormal exit"<<endl;
	      return 1;
	  }
      t = r[j];
      r[j] = (float(ix1) + float(ix2) * rm2) * rm1;
	  return t;
}

double gasdev(int& idum)
{
	static int iset;
	static double gset;
	double t,v1,v2,fac,r;
    if (iset == 0.0)
	{
      do
	  {
        v1 = 2.0 * ran1(idum) - 1.0;
        v2 = 2.0 * ran1(idum) - 1.0;
        r = v1 *v1 + v2 *v2;
	  }
      while ((r >= 1.0) || (r == 0.0));
      fac = sqrt(-2.0 * log(r) / r);
      gset = v1 * fac;
      t = v2 * fac;
	  return t;
      iset = 1;
	}
    else
	{
      t = gset;
	  return t;
      iset = 0;
    }
}

int cint(double x)
{
	int temp;
	double iprt;
	if (x>0)
	{
	x=modf(x,&iprt);
		if(fabs(x)<0.5)
			temp=int(iprt);
		else
			temp=int(iprt+1);
	}
	else if(x==0)
		temp=0;
	else
	{
		x=modf(x,&iprt);
		if(fabs(x)<0.5)
			temp=int(iprt);
		else
			temp=int(iprt)-1;
	}
		return temp;
}

void main()
{
     //program d6r6
     //driver for routine gasdev
     double dist[22],x;
     char text[51] ;
	 int n,np1,nover2,npts,iscal,llen,idum,j,k,klim,i;
     n = 20;
     np1 = n + 1;
     nover2 = n / 2;
     npts = 10000;
     iscal = 400;
     llen = 50;
     idum = -13;
     for (j = 0; j<=np1; j++)
         dist[j] = 0.0;
     for (i = 1; i<=npts; i++)
	 {
		 j = cint(0.25 * n * gasdev(idum)) + nover2 + 1;
         if ((j >= 1) && (j <= np1)) dist[j] = dist[j] + 1;
     }
     cout<<endl;
     cout<<"Normally distributed deviate of  "<<npts<< " points"<<endl;
     cout<<endl;
     cout<<  "   x      p[x]     graph:"<<endl;
     for (j = 1; j<=np1; j++)
	 {
         dist[j] = dist[j] / npts;
         for(k = 1; k<=50; k++)
             text[k] = ' ';
         klim = cint(iscal * dist[j]);
         if (klim > llen ) klim = llen;
         for (k = 1; k<=klim; k++)
             text[k] = '*';
		 x = j / (0.25 * n);
         cout<<setw(5)<<x;
         cout<<setw(10)<<dist[j];
         cout<<"    ";
         for (k = 1; k<=50; k++)
         cout<< text[k];
         cout<<endl;
	 }
}