Tips and Tricks for Programming in Fortran 95

Introduction

This documention contains some tips and tricks for programming in Fortran 95. More will be added from time to time as they are collected.

Random-Number Generation in Fortran 95
Many algorithms call for the generation of "random" numbers. (Mathematically, only pseudorandom numbers can be generated by any finite algorithm, but pseudorandomness is sufficient for nearly all purposes.) In the past, Fortran programmers were forced to call C routines such as rand. Since rand and other such routines were not standard to Fortran, the manner in which they were called and even their availability were dependent on the compiler. Fortran 90 introduced two new intrinsic subroutines, random_seed and random_number, which standardized the process. For most purposes, these intrinsics are sufficient for obtaining statistically valid pseudorandom numbers.

The random_seed intrinsic sets or returns the seeds, while random_number returns a uniformly-distributed random number normalized to lie between 0 and 1. Random_number can return a vector of random numbers as well as single numbers. The simplest usage would be: 


real, dimension(100) :: v
.....
call random_number(v)
This will use the default seeds. If you want to control the seeds more directly, use the intrinsic random_seed before calling random_number. 

call random_seed()
This sets the seed as determined by the compiler's method (often it is based on the system date or time). However, it is often more effective for the user to set the seed directly. For example, for debugging purposes it is often desirable to set the seed to some explicit value (often 0) in order to obtain a reproducible sequence. This is accomplished by code such as 

REAL                 :: rnd
INTEGER              :: isize
INTEGER,ALLOCATABLE  :: iseed(:)
......
  CALL RANDOM_SEED(SIZE=isize)
  ALLOCATE( iseed(isize) )
  CALL RANDOM_SEED(GET=iseed)
  iseed = 0
  CALL RANDOM_NUMBER(RND)
Compilers do not always set the seed in a sufficiently "random" manner. To ensure a randomized seed, you can use the date_and_time intrinsic to obtain a seed from the system date. An example follows; note that this example also shows how to set 100 random numbers simultaneously. 

REAL, DIMENSION(100) :: rnd
INTEGER              :: isize,idate(8)
INTEGER,ALLOCATABLE  :: iseed(:)

  CALL DATE_AND_TIME(VALUES=idate)
  CALL RANDOM_SEED(SIZE=isize)
  ALLOCATE( iseed(isize) )
  CALL RANDOM_SEED(GET=iseed)
  iseed = iseed * (idate(8)-500)      ! idate(8) contains millisecond
  CALL RANDOM_SEED(PUT=iseed)

  CALL RANDOM_NUMBER(RND)

  DEALLOCATE( iseed )
If you need state-of-the-art (for most purposes other than cryptography) random number generation, you can also try the Mersenne Twister. This is described at http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/emt.html. A Fortran 90 version of the Mersenne Twister algorithm was prepared by Josi Rui Faustino de Sousa from the original C code written by Takuji Nishimura and Makoto Matsumoto. We have placed a copy here. Here is an example of the usage of the module: 

PROGRAM testmt
USE mt19937
IMPLICIT NONE

REAL    :: rnd
INTEGER :: num = 10
INTEGER :: idate(8)
INTEGER :: i, iseed

  CALL DATE_AND_TIME(VALUES=idate)
  iseed = idate(7)*1000 + idate(8)
  CALL init_genrand(iseed)

  DO i=1,num
    rnd = genrand_real1()
    PRINT *,rnd
  END DO

  STOP
END PROGRAM testmt
Many commercial libraries such as IMSL and and the Intel Math Kernel Libraries (both licensed by ITC) also contain random-number generating subroutines. Some of these routines can return pseudorandom numbers distributed other than uniformly, such as normally or over a unit sphere. See the documention for IMSL or the Intel MKL for more information.

Easy Switching Between Single and Double Precision
Often, it is desirable to switch from one precision to another when compiling for different systems or compilers. Fortran 90/95 makes this fairly simple. Create a module; this example will call the module prec. 

module prec

  !Uncomment the appropriate line in order to select the desired precision
  !These kinds are tailored to IEEE 754 standard kinds.

  ! Single precision
  !integer, parameter   :: rk = selected_real_kind(6,37)

  ! Double precision
  integer, parameter   :: rk = selected_real_kind(15,307)

end module prec
Then in every program unit, 

use prec
Remember that the use statement must immediately follow the unit invocation (program, subroutine, function). Then declare all real variables using the form 

real (rk)  :: r, s, t
and so forth. It will be necessary to recompile every time you change from real to double or vice versa, but this module makes it much simpler.

Remember that in Fortran 90/95, all intrinsic functions with real arguments are automatically overloaded to handle single precision, double precision, and complex. It is no longer necessary to use e.g. dcos or ccos; cos(x) will behave appropriately regardless of whether x is real, double, or complex.

© 2009 by the Rector and Visitors of the University of Virginia.

The information contained on the University of Virginia’s Department of Information Technology and Communication (ITC) website is provided as a public service with the understanding that ITC makes no representations or warranties, either expressed or implied, concerning the accuracy, completeness, reliability or suitability of the information, including warrantees of title, non-infringement of copyright or patent rights of others. These pages are expected to represent the University of Virginia community and the State of Virginia in a professional manner in accordance with the University of Virginia’s Computing Policies.