A Modular Arbitrary-Order Ocean-Atmosphere Model -- Stochastic implementation
Functions/Subroutines
util Module Reference

Utility module. More...

Functions/Subroutines

character(len=20) function, public str (k)
 Convert an integer to string. More...
 
character(len=40) function, public rstr (x, fm)
 Convert a real to string with a given format. More...
 
integer function, dimension(size(s)), public isin (c, s)
 Determine if a character is in a string and where. More...
 
subroutine, public init_random_seed ()
 Random generator initialization routine. More...
 
subroutine, public piksrt (k, arr, par)
 Simple card player sorting function. More...
 
subroutine, public init_one (A)
 Initialize a square matrix A as a unit matrix. More...
 
real(kind=8) function, public mat_trace (A)
 
real(kind=8) function, public mat_contract (A, B)
 
subroutine, public choldc (a, p)
 
subroutine, public printmat (A)
 
subroutine, public cprintmat (A)
 
real(kind=8) function, dimension(size(a, 1), size(a, 2)), public invmat (A)
 
subroutine, public triu (A, T)
 
subroutine, public diag (A, d)
 
subroutine, public cdiag (A, d)
 
integer function, public floordiv (i, j)
 
subroutine, public reduce (A, Ared, n, ind, rind)
 
subroutine, public ireduce (A, Ared, n, ind, rind)
 
subroutine, public vector_outer (u, v, A)
 

Detailed Description

Utility module.

Function/Subroutine Documentation

subroutine, public util::cdiag ( complex(kind=16), dimension(:,:), intent(in)  A,
complex(kind=16), dimension(:), intent(out)  d 
)

Definition at line 269 of file util.f90.

269  COMPLEX(KIND=16), DIMENSION(:,:), INTENT(IN) :: a
270  COMPLEX(KIND=16), DIMENSION(:), INTENT(OUT) :: d
271  INTEGER :: i
272 
273  DO i=1,SIZE(a,1)
274  d(i)=a(i,i)
275  END DO
subroutine, public util::choldc ( real(kind=8), dimension(:,:)  a,
real(kind=8), dimension(:)  p 
)

Definition at line 176 of file util.f90.

176  REAL(KIND=8), DIMENSION(:,:) :: a
177  REAL(KIND=8), DIMENSION(:) :: p
178  INTEGER :: n
179  INTEGER :: i,j,k
180  REAL(KIND=8) :: sum
181  n=size(a,1)
182  DO i=1,n
183  DO j=i,n
184  sum=a(i,j)
185  DO k=i-1,1,-1
186  sum=sum-a(i,k)*a(j,k)
187  END DO
188  IF (i.eq.j) THEN
189  IF (sum.le.0.) stop 'choldc failed'
190  p(i)=sqrt(sum)
191  ELSE
192  a(j,i)=sum/p(i)
193  ENDIF
194  END DO
195  END DO
196  RETURN
subroutine, public util::cprintmat ( complex(kind=16), dimension(:,:), intent(in)  A)

Definition at line 209 of file util.f90.

209  COMPLEX(KIND=16), DIMENSION(:,:), INTENT(IN) :: a
210  INTEGER :: i
211 
212  DO i=1,SIZE(a,1)
213  print*, a(i,:)
214  END DO
subroutine, public util::diag ( real(kind=8), dimension(:,:), intent(in)  A,
real(kind=8), dimension(:), intent(out)  d 
)

Definition at line 259 of file util.f90.

259  REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: a
260  REAL(KIND=8), DIMENSION(:), INTENT(OUT) :: d
261  INTEGER :: i
262 
263  DO i=1,SIZE(a,1)
264  d(i)=a(i,i)
265  END DO
integer function, public util::floordiv ( integer  i,
integer  j 
)

Definition at line 280 of file util.f90.

280  INTEGER :: i,j,floordiv
281  floordiv=int(floor(real(i)/real(j)))
282  RETURN
subroutine, public util::init_one ( real(kind=8), dimension(:,:), intent(inout)  A)

Initialize a square matrix A as a unit matrix.

Definition at line 139 of file util.f90.

139  REAL(KIND=8), DIMENSION(:,:),INTENT(INOUT) :: a
140  INTEGER :: i,n
141  n=size(a,1)
142  a=0.0d0
143  DO i=1,n
144  a(i,i)=1.0d0
145  END DO
146 
subroutine, public util::init_random_seed ( )

Random generator initialization routine.

Definition at line 64 of file util.f90.

real(kind=8) function, dimension(size(a,1),size(a,2)), public util::invmat ( real(kind=8), dimension(:,:), intent(in)  A)

Definition at line 218 of file util.f90.

218  REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: a
219  REAL(KIND=8), DIMENSION(SIZE(A,1),SIZE(A,2)) :: ainv
220 
221  REAL(KIND=8), DIMENSION(SIZE(A,1)) :: work ! work array for LAPACK
222  INTEGER, DIMENSION(SIZE(A,1)) :: ipiv ! pivot indices
223  INTEGER :: n, info
224 
225  ! Store A in Ainv to prevent it from being overwritten by LAPACK
226  ainv = a
227  n = size(a,1)
228 
229  ! DGETRF computes an LU factorization of a general M-by-N matrix A
230  ! using partial pivoting with row interchanges.
231  CALL dgetrf(n, n, ainv, n, ipiv, info)
232 
233  IF (info /= 0) THEN
234  stop 'Matrix is numerically singular!'
235  ENDIF
236 
237  ! DGETRI computes the inverse of a matrix using the LU factorization
238  ! computed by DGETRF.
239  CALL dgetri(n, ainv, n, ipiv, work, n, info)
240 
241  IF (info /= 0) THEN
242  stop 'Matrix inversion failed!'
243  ENDIF
subroutine, public util::ireduce ( real(kind=8), dimension(:,:), intent(out)  A,
real(kind=8), dimension(:,:), intent(in)  Ared,
integer, intent(in)  n,
integer, dimension(:), intent(in)  ind,
integer, dimension(:), intent(in)  rind 
)

Definition at line 314 of file util.f90.

314  REAL(KIND=8), DIMENSION(:,:), INTENT(OUT) :: a
315  REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: ared
316  INTEGER, INTENT(IN) :: n
317  INTEGER, DIMENSION(:), INTENT(IN) :: ind,rind
318  INTEGER :: i,j
319  a=0.d0
320  DO i=1,n
321  DO j=1,n
322  a(ind(i),ind(j))=ared(i,j)
323  END DO
324  END DO
integer function, dimension(size(s)), public util::isin ( character, intent(in)  c,
character, dimension(:), intent(in)  s 
)

Determine if a character is in a string and where.

Remarks
: return positions in a vector if found and 0 vector if not found

Definition at line 47 of file util.f90.

47  CHARACTER, INTENT(IN) :: c
48  CHARACTER, DIMENSION(:), INTENT(IN) :: s
49  INTEGER, DIMENSION(size(s)) :: isin
50  INTEGER :: i,j
51 
52  isin=0
53  j=0
54  DO i=size(s),1,-1
55  IF (c==s(i)) THEN
56  j=j+1
57  isin(j)=i
58  END IF
59  END DO
real(kind=8) function, public util::mat_contract ( real(kind=8), dimension(:,:)  A,
real(kind=8), dimension(:,:)  B 
)

Definition at line 162 of file util.f90.

162  REAL(KIND=8), DIMENSION(:,:) :: a,b
163  REAL(KIND=8) :: mat_contract
164  INTEGER :: i,j,n
165  n=size(a,1)
166  mat_contract=0.d0
167  DO i=1,n
168  DO j=1,n
169  mat_contract=mat_contract+a(i,j)*b(i,j)
170  END DO
171  ENDDO
172  RETURN
real(kind=8) function, public util::mat_trace ( real(kind=8), dimension(:,:)  A)

Definition at line 150 of file util.f90.

150  REAL(KIND=8), DIMENSION(:,:) :: a
151  REAL(KIND=8) :: mat_trace
152  INTEGER :: i,n
153  n=size(a,1)
154  mat_trace=0.d0
155  DO i=1,n
156  mat_trace=mat_trace+a(i,i)
157  END DO
158  RETURN
subroutine, public util::piksrt ( integer, intent(in)  k,
integer, dimension(k), intent(inout)  arr,
integer, intent(out)  par 
)

Simple card player sorting function.

Definition at line 118 of file util.f90.

118  INTEGER, INTENT(IN) :: k
119  INTEGER, DIMENSION(k), INTENT(INOUT) :: arr
120  INTEGER, INTENT(OUT) :: par
121  INTEGER :: i,j,a,b
122 
123  par=1
124 
125  DO j=2,k
126  a=arr(j)
127  DO i=j-1,1,-1
128  if(arr(i).le.a) EXIT
129  arr(i+1)=arr(i)
130  par=-par
131  END DO
132  arr(i+1)=a
133  ENDDO
134  RETURN
subroutine, public util::printmat ( real(kind=8), dimension(:,:), intent(in)  A)

Definition at line 200 of file util.f90.

200  REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: a
201  INTEGER :: i
202 
203  DO i=1,SIZE(a,1)
204  print*, a(i,:)
205  END DO
subroutine, public util::reduce ( real(kind=8), dimension(:,:), intent(in)  A,
real(kind=8), dimension(:,:), intent(out)  Ared,
integer, intent(out)  n,
integer, dimension(:), intent(out)  ind,
integer, dimension(:), intent(out)  rind 
)

Definition at line 286 of file util.f90.

286  REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: a
287  REAL(KIND=8), DIMENSION(:,:), INTENT(OUT) :: ared
288  INTEGER, INTENT(OUT) :: n
289  INTEGER, DIMENSION(:), INTENT(OUT) :: ind,rind
290  LOGICAL, DIMENSION(SIZE(A,1)) :: sel
291  INTEGER :: i,j
292 
293  ind=0
294  rind=0
295  sel=.false.
296  n=0
297  DO i=1,SIZE(a,1)
298  IF (any(a(i,:)/=0)) THEN
299  n=n+1
300  sel(i)=.true.
301  ind(n)=i
302  rind(i)=n
303  ENDIF
304  END DO
305  ared=0.d0
306  DO i=1,SIZE(a,1)
307  DO j=1,SIZE(a,1)
308  IF (sel(i).and.sel(j)) ared(rind(i),rind(j))=a(i,j)
309  ENDDO
310  ENDDO
character(len=40) function, public util::rstr ( real(kind=8), intent(in)  x,
character(len=20), intent(in)  fm 
)

Convert a real to string with a given format.

Definition at line 38 of file util.f90.

38  REAL(KIND=8), INTENT(IN) :: x
39  CHARACTER(len=20), INTENT(IN) :: fm
40  WRITE (rstr, trim(adjustl(fm))) x
41  rstr = adjustl(rstr)
character(len=20) function, public util::str ( integer, intent(in)  k)

Convert an integer to string.

Definition at line 31 of file util.f90.

31  INTEGER, INTENT(IN) :: k
32  WRITE (str, *) k
33  str = adjustl(str)
character(len=20) function, public str(k)
Convert an integer to string.
Definition: util.f90:31
subroutine, public util::triu ( real(kind=8), dimension(:,:), intent(in)  A,
real(kind=8), dimension(:,:), intent(out)  T 
)

Definition at line 247 of file util.f90.

247  REAL(KIND=8), DIMENSION(:,:), INTENT(IN) :: a
248  REAL(KIND=8), DIMENSION(:,:), INTENT(OUT) :: t
249  INTEGER i,j
250  t=0.d0
251  DO i=1,SIZE(a,1)
252  DO j=i,SIZE(a,1)
253  t(i,j)=a(i,j)
254  END DO
255  END DO
subroutine, public util::vector_outer ( real(kind=8), dimension(:), intent(in)  u,
real(kind=8), dimension(:), intent(in)  v,
real(kind=8), dimension(:,:), intent(out)  A 
)

Definition at line 328 of file util.f90.

328  REAL(KIND=8), DIMENSION(:), INTENT(IN) :: u,v
329  REAL(KIND=8), DIMENSION(:,:), INTENT(OUT) :: a
330  INTEGER :: i,j
331 
332  a=0.d0
333  DO i=1,SIZE(u)
334  DO j=1,SIZE(v)
335  a(i,j)=u(i)*v(j)
336  ENDDO
337  ENDDO