File: /nfs/home/0/users/jenkins/mfix.git/model/des/randomno_mod.f

1     ! vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2     !                                                                      C
3     !  Module name: Random Number Generation Utilities                     C
4     !  Purpose: Removed from interpolation mod and added built-in random
5     !           number routines instead of Pope's
6     !                                                                      C
7     !                                                                      C
8     !  Author: Sreekanth Pannala and Rahul Garg           Date: 23-Oct-08  C
9     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
10         MODULE randomno
11     
12         USE constant
13     
14         IMPLICIT NONE
15     
16         PRIVATE
17         PUBLIC :: uni_rno, nor_rno
18     
19     
20         CONTAINS
21     
22     
23     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv
24           SUBROUTINE UNI_RNO(Y)
25     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
26           IMPLICIT NONE
27     
28     !-----------------------------------------------
29     ! Local variables
30     !-----------------------------------------------
31           double precision, intent(out), dimension(:) :: y
32           double precision rmean, variance, sigma
33           integer i, nsize
34     !-----------------------------------------------
35     
36           nsize = size(y(:))
37     
38           call init_random_seed
39           call random_number(y)
40     
41           rmean = sum(y(:))/nsize
42     
43     !      write(*,*) 'Generating Uniform Random Variables for size', nsize
44     !      write(*,*) 'mean', rmean
45     
46           variance = 0.0
47           do i = 1, nsize
48     !         write(20,*) i, y(i)
49              variance = variance + (y(i)-rmean)**2
50           end do
51     
52           close(20)
53     
54           variance = variance/nsize
55           sigma = sqrt(variance)
56     
57     !      write(*,*) 'sigma', sigma
58     
59           RETURN
60           END SUBROUTINE UNI_RNO
61     
62     
63     
64     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv
65           SUBROUTINE NOR_RNO(Y, mean, sigma)
66     
67     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
68           IMPLICIT NONE
69     
70     !-----------------------------------------------
71     ! Local variables
72     !-----------------------------------------------
73           double precision, intent(out), dimension(:) :: y
74           double precision mean, sigma
75     
76           double precision lmean, lvariance, lsigma
77           double precision x(2), w
78           integer i, nsize, n
79     ! no. of times this routine has been called
80           integer, save :: COUNTER = 0
81     ! so all components are written
82           integer fileunit
83     !-----------------------------------------------
84           COUNTER = COUNTER + 1
85           fileunit = 20 + COUNTER
86     
87           nsize = size(y(:))
88     
89           call init_random_seed
90     
91           do i = 1, ceiling(real(nsize/2))
92              do n = 1,100000
93                 call random_number(x)
94                 x = 2.0 * x - 1.0
95                 w = x(1)**2 + x(2)**2
96                 if(w.lt.1.0) exit
97              end do
98     
99              w = sqrt( (-2.0 * log( w ) ) / w )
100              y(2*i-1) = x(1) * w * sigma + mean
101              if(2*i.lt.nsize) y(2*i) = x(2) * w * sigma + mean
102           end do
103     
104           lmean = sum(y(:))/nsize
105     
106           !write(*,'(7X,A)') 'Generating Normal Random Variables'
107           !write(*,'(7X,A,2X,ES15.5)') 'specified mean =', mean
108           !write(*,'(7X,A,2X,ES15.5)') 'computed mean =', lmean
109     
110           !write(fileunit,'(A)') 'FROM NOR_RNO'
111     ! specific to the call from init_particles_jn
112           !write(fileunit,'(A,I5,A)') 'FOR DIRECTION = ', &
113           !   COUNTER, ' where (1=X,2=Y,3=Z)'
114           !write(fileunit,'(5X,A,5X,A)') 'particle no.', 'velocity component'
115     
116           lvariance = 0.0
117           do i = 1, nsize
118              !write(fileunit,'(I10,5X,ES15.5)') i, y(i)
119              lvariance = lvariance + (y(i)-lmean)**2
120           end do
121     
122           !close(fileunit)
123     
124           lvariance = lvariance/nsize
125           lsigma = sqrt(lvariance)
126     
127           !write(*,'(7X,A,2X,ES15.5)') 'specified sigma =', sigma
128           !write(*,'(7X,A,2X,ES15.5)') 'computed sigma =', lsigma
129     
130           RETURN
131           END SUBROUTINE NOR_RNO
132     
133     
134     
135     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv
136           SUBROUTINE init_random_seed
137     
138     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
139           IMPLICIT NONE
140     
141     !-----------------------------------------------
142     ! Local variables
143     !-----------------------------------------------
144           INTEGER              :: isize,idate(8)
145           INTEGER,ALLOCATABLE  :: iseed(:)
146     !-----------------------------------------------
147     
148           CALL DATE_AND_TIME(VALUES=idate)
149           CALL RANDOM_SEED(SIZE=isize)
150           ALLOCATE( iseed(isize) )
151           CALL RANDOM_SEED(GET=iseed)
152           iseed = iseed * (idate(8)-500) ! idate(8) contains millisecond
153           CALL RANDOM_SEED(PUT=iseed)
154     
155           DEALLOCATE( iseed )
156     
157           END SUBROUTINE init_random_seed
158     
159     
160         END MODULE randomno
161     
162