// A simple random number generator and an example of a program which
// calls it.
//
float randpy(int *seed)
{
    int  IMAX = 2147483647;
    float INV_IMAX = 0.46566128730e-9 ;
    *seed *= 16807;	// multiply seed by the "magic" number 16807
    *seed &= IMAX;	// take bitwise "and" with IMAX (= 2147483647) 
    			// (this removes a possible minus sign)
    return *seed * INV_IMAX; // convert to a float between 0 and 1
}

#include <math.h>
#include <stdio.h>

main()
{
    float randpy(int *seed);
    double sum1, sum2, x, sigma;
    int i, seed, N;

    printf ("seed = ? N = ? " );
    scanf("%d %d", &seed, &N);	  // read in the seed, and no. of random numbers
    				  // to be averaged over
    printf ("\nseed = %d, N = %d \n\n", seed, N) ;

    sum1 = 0; sum2 = 0;
    for (i = 0; i < N; i++) 		
    {
	x = randpy(&seed);
	sum1 += x;
	sum2 += x*x;
    }
    sum1 /= N;
    sum2 /= N;
    sigma = sqrt(sum2 - sum1*sum1); // print mean and standard deviation
    printf ("mean = %10.4f, sigma = %10.4f \n\n", sum1, sigma) ;
    printf ("mean should equal 1/2, sigma should equal 1/sqrt(12) = 0.2887\n\n");
}