static double ran2(int& idum)
{
	static int iy,iff,ir[97];
	double rm,ran;
    int j,m,ic,ia;
    m = 714025; ia = 1366; ic = 150889; 
    if (idum < 0 || iff == 0) 
	{
        iff = 1;
        idum = (ic - idum) % m;
        for( j = 1; j<=97; j++)
		{
            idum = (ia*idum + ic) % m;
            ir[j] = idum;
        }
        idum = (ia * idum + ic) % m;
        iy = idum;
    }
    j = 1 + int((97 * iy) / m);
    if(j > 97 || j < 1)
    {
		cout<<"abnormal exit"; 
		exit(1);
	}
    rm = 0.0000014005112;
    iy = ir[j];
    ran = iy * rm;
    idum = (ia * idum + ic) % m;
    ir[j] = idum;
	return ran;
}