/[MITgcm]/MITgcm/model/src/port_rand.F
ViewVC logotype

Diff of /MITgcm/model/src/port_rand.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph | View Patch Patch

revision 1.1 by adcroft, Thu Jan 25 18:44:07 2001 UTC revision 1.2 by adcroft, Fri Feb 2 21:04:48 2001 UTC
# Line 0  Line 1 
1    C $Header$
2    C $Name$
3    
4    C Portable random number generator
5    
6    #undef _USE_INTEGERS
7    
8    ! ==============================================================================
9          real*8 function port_rand()
10          implicit none
11          integer nff
12          parameter(nff=55)
13    #ifdef _USE_INTEGERS
14          integer mbig,mseed,mZ
15    #else
16          real*8 mbig,mseed,mz
17    #endif
18          real*8 fac
19    #ifdef _USE_INTEGERS
20          parameter (mbig=1000000000,mseed=161803398,mz=0,fac=1./mbig)
21    #else
22          parameter (mbig=4000000.,mseed=1618033.,mz=0.,fac=1./mbig)
23    #endif
24          integer i,ii,inext,inextp,k,idum
25          parameter(idum=-2)
26    #ifdef _USE_INTEGERS
27          integer mj,mk,ma(nff)
28    #else
29          real*8 mj,mk,ma(nff)
30    #endif
31          logical firstCall
32          save firstCall,inext,inextp,ma
33          data firstCall /.true./
34    ! ------------------------------------------------------------------------------
35          if(firstCall)then
36            firstCall=.false.
37            mj=mseed-iabs(idum)
38            mj=mod(mj,mbig)
39            ma(nff)=mj
40            mk=1
41            do i=1,nff-1
42              ii=mod(21*i,nff)
43              ma(ii)=mk
44              mk=mj-mk
45              if(mk.lt.mz)mk=mk+mbig
46              mj=ma(ii)
47            enddo
48            do k=1,4
49              do i=1,nff
50                ma(i)=ma(i)-ma(1+mod(i+30,nff))
51                if(ma(i).lt.MZ)ma(i)=ma(i)+mbig
52              enddo
53            enddo
54            inext=0
55            inextp=31
56          endif
57          inext=mod(inext,nff)+1
58          inextp=mod(inextp,nff)+1
59          mj=ma(inext)-ma(inextp)
60          if(mj.lt.MZ)mj=mj+mbig
61          ma(inext)=mj
62          port_rand=mj*fac
63          return
64    ! ------------------------------------------------------------------------------
65          end
66    ! ==============================================================================
67    
68    ! ==============================================================================
69          subroutine port_ranarr(n,arr)
70          implicit none
71          integer n,i
72          real arr(n)
73          real*8 port_rand
74    ! ------------------------------------------------------------------------------
75          do i=1,n
76           arr(i)=port_rand()
77          enddo
78    
79          return
80    ! ------------------------------------------------------------------------------
81          end
82    ! ==============================================================================

Legend:
Removed from v.1.1  
changed lines
  Added in v.1.2

  ViewVC Help
Powered by ViewVC 1.1.22