Merge branch 'devel' into AFM
[unres.git] / source / unres / src_MD-M-newcorr / ran.f
diff --git a/source/unres/src_MD-M-newcorr/ran.f b/source/unres/src_MD-M-newcorr/ran.f
new file mode 100644 (file)
index 0000000..dd23252
--- /dev/null
@@ -0,0 +1,128 @@
+ccccccccccccccccccccccccccccccccccccccccccccccccc
+      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.lt.0) idum=idum+IM
+      ran0=AM*idum
+      idum=ieor(idum,MASK)
+      return
+      END
+C  (C) Copr. 1986-92 Numerical Recipes Software *11915
+ccccccccccccccccccccccccccccccccccccccccccccccccc
+      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.le.0.or.iy.eq.0) then
+        idum=max(-idum,1)
+        do 11 j=NTAB+8,1,-1
+          k=idum/IQ
+          idum=IA*(idum-k*IQ)-IR*k
+          if (idum.lt.0) idum=idum+IM
+          if (j.le.NTAB) iv(j)=idum
+11      continue
+        iy=iv(1)
+      endif
+      k=idum/IQ
+      idum=IA*(idum-k*IQ)-IR*k
+      if (idum.lt.0) idum=idum+IM
+      j=1+iy/NDIV
+      iy=iv(j)
+      iv(j)=idum
+      ran1=min(AM*iy,RNMX)
+      return
+      END
+C  (C) Copr. 1986-92 Numerical Recipes Software *11915
+ccccccccccccccccccccccccccccccccccccccccccccccccc
+      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.le.0) then
+        idum=max(-idum,1)
+        idum2=idum
+        do 11 j=NTAB+8,1,-1
+          k=idum/IQ1
+          idum=IA1*(idum-k*IQ1)-k*IR1
+          if (idum.lt.0) idum=idum+IM1
+          if (j.le.NTAB) iv(j)=idum
+11      continue
+        iy=iv(1)
+      endif
+      k=idum/IQ1
+      idum=IA1*(idum-k*IQ1)-k*IR1
+      if (idum.lt.0) idum=idum+IM1
+      k=idum2/IQ2
+      idum2=IA2*(idum2-k*IQ2)-k*IR2
+      if (idum2.lt.0) idum2=idum2+IM2
+      j=1+iy/NDIV
+      iy=iv(j)-idum2
+      iv(j)=idum
+      if(iy.lt.1)iy=iy+IMM1
+      ran2=min(AM*iy,RNMX)
+      return
+      END
+C  (C) Copr. 1986-92 Numerical Recipes Software *11915
+ccccccccccccccccccccccccccccccccccccccccccccccccc
+      FUNCTION ran3(idum)
+      INTEGER idum
+      INTEGER MBIG,MSEED,MZ
+C     REAL MBIG,MSEED,MZ
+      REAL ran3,FAC
+      PARAMETER (MBIG=1000000000,MSEED=161803398,MZ=0,FAC=1./MBIG)
+C     PARAMETER (MBIG=4000000.,MSEED=1618033.,MZ=0.,FAC=1./MBIG)
+      INTEGER i,iff,ii,inext,inextp,k
+      INTEGER mj,mk,ma(55)
+C     REAL mj,mk,ma(55)
+      SAVE iff,inext,inextp,ma
+      DATA iff /0/
+      if(idum.lt.0.or.iff.eq.0)then
+        iff=1
+        mj=MSEED-iabs(idum)
+        mj=mod(mj,MBIG)
+        ma(55)=mj
+        mk=1
+        do 11 i=1,54
+          ii=mod(21*i,55)
+          ma(ii)=mk
+          mk=mj-mk
+          if(mk.lt.MZ)mk=mk+MBIG
+          mj=ma(ii)
+11      continue
+        do 13 k=1,4
+          do 12 i=1,55
+            ma(i)=ma(i)-ma(1+mod(i+30,55))
+            if(ma(i).lt.MZ)ma(i)=ma(i)+MBIG
+12        continue
+13      continue
+        inext=0
+        inextp=31
+        idum=1
+      endif
+      inext=inext+1
+      if(inext.eq.56)inext=1
+      inextp=inextp+1
+      if(inextp.eq.56)inextp=1
+      mj=ma(inext)-ma(inextp)
+      if(mj.lt.MZ)mj=mj+MBIG
+      ma(inext)=mj
+      ran3=mj*FAC
+      return
+      END
+C  (C) Copr. 1986-92 Numerical Recipes Software *11915
+ccccccccccccccccccccccccccccccccccccccccccccccccc