乱数

#include <stdio.h>
#include <time.h>

#define MBIG 1000000000
#define MSEED 161803398
#define MZ 0
#define FAC (1.0 / MBIG)

float ran3(long *idum){
    
    static int inext, inextp;
    static long ma[56];
    static int iff = 0;
    long mj, mk;
    int i, ii, j, k;
    
    if(((*idum) < 0) || (iff == 0)){
        iff = 1;
        if(*idum < 0){
            mj = MSEED - (-(*idum));
        }else{
            mj = MSEED - (*idum);
        }
        mj %= MBIG;
        ma[55] = mj;
        
        mk = 1;
        for(i = 1; i <= 54; i += 1){
            ii = (21 * i) % 55;
            ma[ii] = mk;
            mk = mj - mk;
            if(mk < MZ){
                mk += MBIG;
            }
            mj = ma[ii];
        }
        
        for(k = 1; k <= 4; k += 1){
            for(i = 1; i <= 55; i += 1){
                j = 1 + (i + 30) % 55;
                ma[i] -= ma[j];
                if(ma[i] < MZ){
                    ma[i] += MBIG;
                }
            }
        }
        inext = 0;
        inextp = 31;
        *idum = 1;
    }
    
    inext += 1;
    if(inext == 56){
        inext = 1;
    }
    inextp += 1;
    if(inextp == 56){
        inextp = 1;
    }
    mj = ma[inext] - ma[inextp];
    if(mj < MZ){
        mj += MBIG;
    }
    ma[inext] = mj;
    return (mj * FAC);
}

int main(int argc, const char * argv[]) {
    
    int i;
    long idum;
    float ran;
    
    time(&idum);
    for(i = 1; i <= 5; i += 1){
        ran = ran3(&idum);
        printf("%.3f\n", ran);
    }
    
    return 0;
}