四个随机数程序

程序一：
FUNCTION ran0(idum)
INTEGER idum,IA,IM,IQ,IR,MASK
REAL ran0,AM
PARAMETER (IA=16807,IM=2147483647,AM=1./IM,&
                   IQ=127773,IR=2836,MASK=123459876)
INTEGER k
idum=ieor(idum,MASK)
k=idum/IQ
idum=IA*(idum-k*IQ)-IR*k
if (idum<0) idum=idum+IM
ran0=AM*idum
idum=ieor(idum,MASK)
END FUNCTION ran0

程序二：
FUNCTION ran1(idum)
INTEGER idum,ia,im,iq,ir,ntab,ndiv
REAL ran1,am,eps,rnmx
PARAMETER (ia=16807,im=2147483647,am=1./im,&
     iq=127773,ir=2836,ntab=32,ndiv=1+(im-1)/ntab,&
  eps=1.2e-7,rnmx=1.-eps)
INTEGER j,k,iv(ntab),iy
SAVE iv,iy
DATA iv /ntab*0/, iy /0/
if (idum<=0.or.iy==0) then
  idum=max(-idum,1)
  do j=ntab+8,1,-1
    k=idum/iq
    idum=ia*(idum-k*iq)-ir*k
    if (idum<0) idum=idum+im
    if (j<=ntab) iv(j)=idum
  end do
  iy=iv(1)
endif
k=idum/iq
idum=ia*(idum-k*iq)-ir*k
if (idum<0) idum=idum+im
j=1+iy/ndiv
iy=iv(j)
iv(j)=idum
ran1=min(am*iy,rnmx)
END FUNCTION ran1

程序三：
FUNCTION ran2(idum)
INTEGER idum,IM1,IM2,IMM1,IA1,IA2,IQ1,IQ2,&
        IR1,IR2,NTAB,NDIV
REAL ran2,AM,EPS,RNMX
PARAMETER (IM1=2147483563,IM2=2147483399,AM=1./IM1,&
           IMM1=IM1-1,IA1=40014,IA2=40692,IQ1=53668,&
     IQ2=52774,IR1=12211,IR2=3791,NTAB=32,&
     NDIV=1+IMM1/NTAB,EPS=1.2e-7,RNMX=1.-EPS)
INTEGER idum2,j,k,iv(NTAB),iy
SAVE iv,iy,idum2
DATA idum2/123456789/, iv/NTAB*0/, iy/0/
if (idum<=0) then
  idum=max(-idum,1)
  idum2=idum
  do j=NTAB+8,1,-1
    k=idum/IQ1
    idum=IA1*(idum-k*IQ1)-k*IR1
    if (idum<0) idum=idum+IM1
    if (j<=NTAB) iv(j)=idum
  end do
  iy=iv(1)
endif
k=idum/IQ1
idum=IA1*(idum-k*IQ1)-k*IR1
if (idum<0) idum=idum+IM1
k=idum2/IQ2
idum2=IA2*(idum2-k*IQ2)-k*IR2
if (idum2<0) idum2=idum2+IM2
j=1+iy/NDIV
iy=iv(j)-idum2
iv(j)=idum
if(iy<1) iy=iy+IMM1
ran2=min(AM*iy,RNMX)
END FUNCTION ran2

程序四：
FUNCTION ran3(idum)
INTEGER idum
INTEGER MBIG,MSEED,MZ
!REAL MBIG,MSEED,MZ
REAL ran3,FAC
PARAMETER (MBIG=1000000000,MSEED=161803398,MZ=0,FAC=1./MBIG)
!PARAMETER (MBIG=4000000.,MSEED=1618033.,MZ=0.,FAC=1./MBIG)
INTEGER i,iff,ii,inext,inextp,k
INTEGER mj,mk,ma(55)
!REAL mj,mk,ma(55)
SAVE iff,inext,inextp,ma
DATA iff /0/
if(idum<0.or.iff==0) then
  iff=1
  mj=MSEED-iabs(idum)
  mj=mod(mj,MBIG)
  ma(55)=mj
  mk=1
  do i=1,54
    ii=mod(21*i,55)
    ma(ii)=mk
    mk=mj-mk
    if(mk    mj=ma(ii)
  end do
  do k=1,4
    do i=1,55
      ma(i)=ma(i)-ma(1+mod(i+30,55))
      if(ma(i)    end do
  end do
  inext=0
  inextp=31
  idum=1
endif
inext=inext+1
if(inext==56) inext=1
inextp=inextp+1
if(inextp==56) inextp=1
mj=ma(inext)-ma(inextp)
if(mjma(inext)=mj
ran3=mj*FAC
END FUNCTION ran3