[ previous ] [ Contents ] [ 1 ] [ 2 ] [ 3 ] [ 4 ] [ 5 ] [ 6 ] [ 7 ] [ 8 ] [ 9 ] [ 10 ] [ 11 ] [ next ]


Some Mini-Howtos of Interest
Chapter 6 - Fortran Tips and Templates


6.1 Compile and install the LAPACK95 interface driver routine.

Updated on February 22nd, 2018.

Updated on September 16th, 2011.

This howto was checked initially in a Debian Squeeze box, with a gfortran 4.4.5 compiler and in a Ubuntu Lucid Lynx with a compiler gfortran 4.4.3. Last version has been checked in a Debian Stretch box, with a gfortran 6.3.0 compiler.

In this howto we install from scratch the Fortran 77 libraries BLAS and LAPACK, and then compile the FORTRAN 95 extension to LAPACK. Both Debian and Ubuntu systems provide BLAS and LAPACK libraries as native packages and this first step could be skipped and you can ignore next two paragraphs. In any case it may be interesting to compile all the libraries with optimized flags for the system in question.

First we should upload the LAPACK and LAPACK 95 libraries, for example from the Netlib website (LAPACK tgz from Netlib, (LAPACK95 tgz from Netlib ). The LAPACK version we will install is 3.3.1.

      
     $ wget  -c http://www.netlib.org/lapack/lapack-3.3.1.tgz
     --2011-09-16 13:49:37--  http://www.netlib.org/lapack/lapack-3.3.1.tgz
     Resolving www.netlib.org... 160.36.58.108
     .
     .
     .
     2011-09-16 13:49:48 (439 KB/s) - `lapack-3.3.1.tgz' saved [4945204/4945204]
     
     $ wget http://www.netlib.org/lapack95/lapack95.tgz
     --2011-09-16 13:49:54--  http://www.netlib.org/lapack95/lapack95.tgz
     .
     .
     .
     2011-09-16 13:49:59 (368 KB/s) - `lapack95.tgz' saved [1579613/1579613]

Untar the LAPACK library tarball.

     $ tar xzf lapack-3.3.1.tgz 
     $ cd lapack-3.3.1/
     lapack-3.3.1$

Edit the file make.inc to a file conveniently tuned for your system. A working example for gfortran is shown in make.inc for LAPACK compilation, Section 6.1.1. The library can now be compiled, starting with the compilation of the included BLAS library

     lapack-3.3.1$ cd BLAS/SRC
     lapack-3.3.1/BLAS/SRC$ make
     .
     .
     .
      zhemm.o zherk.o zher2k.o lsame.o xerbla.o xerbla_array.o
     ranlib ../../blas_linux_gfortran.a
     lapack-3.3.1/BLAS/SRC$

The next step is the compilation of the BLAS library

     lapack-3.3.1/BLAS/SRC$ cd ../../
     lapack-3.3.1$ make
     .
     .
     .
     make[1]: Leaving directory `/home/curro/1/lapack-3.3.1/BLAS/TESTING'
     ( cd BLAS; ./xblat3s < sblat3.in     ; \
     	           ./xblat3d < dblat3.in     ; \
     	           ./xblat3c < cblat3.in     ; \
     	           ./xblat3z < zblat3.in     ) 
     lapack-3.3.1$ ls *.a
     blas_linux_gfortran.a  lapack_linux_gfortran.a  tmglib_linux_gfortran.a

The next step is to copy the compiled libraries to their final destination

     lapack-3.3.1$ sudo mkdir /usr/local/lapack-3.3.1
     [sudo] password for curro: 
     lapack-3.3.1$ sudo cp *.a /usr/local/lapack-3.3.1/

We now proceed to compile and install the LAPACK 95 library. This is where you should start if you are using the BLAS and LAPACK provided by the distribution. First we unpack the LAPACK 95 tarball and edit the make.inc file. A working example for gfortran is shown in make.inc for LAPACK95 compilation, Section 6.1.2.

     $ tar xzf lapack95.tgz 
     $ cd LAPACK95/
     LAPACK95$ cd SRC
     LAPACK95/SRC$ make single_double_complex_dcomplex
     .
     .
     .
     ranlib ../lapack95.a
     rm -fr ../lapack95_modules
     mkdir ../lapack95_modules
     'cp' *.mod ../lapack95_modules/
     rm -f f77_lapack.* f95_lapack.*
     rm -f *_lapack_single_double_complex_dcomplex.o
     LAPACK95/SRC$

The final step is to copy the libfile and modules to a convenient location.

     LAPACK95$ sudo mkdir /usr/local/lib/lapack95
     LAPACK95$ sudo cp -r lapack95.a lapack95_modules /usr/local/lapack95

We include a program template invoking this library in Example of program using LAPACK95, Section 6.1.3 and a makefile that can be used to compile this program in makefile for compiling programs with calls to LAPACK95, Section 6.1.4.


6.1.1 make.inc for LAPACK compilation

     # -*- Makefile -*-
     ####################################################################
     #  LAPACK make include file.                                       #
     #  LAPACK, Version 3.3.1                                           #
     #  April 2011                                                      #
     ####################################################################
     #
     # See the INSTALL/ directory for more examples.
     #
     SHELL = /bin/sh
     #
     #  The machine (platform) identifier to append to the library names
     #
     PLAT = _linux_gfortran
     #  
     #  Modify the FORTRAN and OPTS definitions to refer to the
     #  compiler and desired compiler options for your machine.  NOOPT
     #  refers to the compiler options desired when NO OPTIMIZATION is
     #  selected.  Define LOADER and LOADOPTS to refer to the loader
     #  and desired load options for your machine.
     #
     FORTRAN  = gfortran -O2 -m32
     OPTS     =
     DRVOPTS  = $(OPTS)
     NOOPT    = -g -O0
     LOADER   = gfortran -g
     LOADOPTS =
     #
     # Timer for the SECOND and DSECND routines
     #
     # Default : SECOND and DSECND will use a call to the EXTERNAL FUNCTION ETIME
     # TIMER    = EXT_ETIME
     # For RS6K : SECOND and DSECND will use a call to the EXTERNAL FUNCTION ETIME_
     # TIMER    = EXT_ETIME_
     # For gfortran compiler: SECOND and DSECND will use a call to the INTERNAL FUNCTION ETIME
     TIMER    = INT_ETIME
     # If your Fortran compiler does not provide etime (like Nag Fortran Compiler, etc...)
     # SECOND and DSECND will use a call to the Fortran standard INTERNAL FUNCTION CPU_TIME 
     # TIMER    = INT_CPU_TIME
     # If neither of this works...you can use the NONE value... In that case, SECOND and DSECND will always return 0
     # TIMER     = NONE
     #
     #  The archiver and the flag(s) to use when building archive (library)
     #  If you system has no ranlib, set RANLIB = echo.
     #
     ARCH     = ar
     ARCHFLAGS= cr
     RANLIB   = ranlib
     #
     #  The location of BLAS library for linking the testing programs.
     #  The target's machine-specific, optimized BLAS library should be
     #  used whenever possible.
     #
     BLASLIB      = ../../blas$(PLAT).a
     #
     #  Location of the extended-precision BLAS (XBLAS) Fortran library
     #  used for building and testing extended-precision routines.  The
     #  relevant routines will be compiled and XBLAS will be linked only if
     #  USEXBLAS is defined.
     #
     # USEXBLAS    = Yes
     XBLASLIB     =
     # XBLASLIB    = -lxblas
     #
     #  Names of generated libraries.
     #
     LAPACKLIB    = lapack$(PLAT).a
     TMGLIB       = tmglib$(PLAT).a
     EIGSRCLIB    = eigsrc$(PLAT).a
     LINSRCLIB    = linsrc$(PLAT).a

6.1.2 make.inc for LAPACK95 compilation

     #
     #  -- LAPACK95 interface driver routine (version 2.0) --
     #     UNI-C, Denmark; Univ. of Tennessee, USA; NAG Ltd., UK
     #     August 5, 2000
     #
     FC	 = gfortran 
     FC1      = gfortran  
     # -dcfuns  Enable recognition of non-standard double
     #          precision  complex intrinsic functions
     # -dusty   Allows the compilation and execution of "legacy"
     #          software by downgrading the category  of  common
     #          errors  found  in  such software from "Error" to
     # -ieee=full enables all IEEE arithmetic facilities
     #          including non-stop arithmetic.
     OPTS0    = -O3
     MODLIB   = -I./../lapack95_modules
     OPTS1    = -c $(OPTS0)
     OPTS3    = $(OPTS1) $(MODLIB)
     OPTL     = -o 
     OPTLIB   =
     
     #LAPACK_PATH = /usr/local/lib/LAPACK3/
     
     LAPACK95 = ../lapack95.a
     LAPACK77 = -llapack
     BLAS     = -lblas
     
     LIBS     = $(LAPACK95) $(LAPACK77) $(BLAS)
     SUF      = f90
     
     XX = 'rm' -f $@; \
             'rm' -f $@.res; \
     	$(FC) $(OPTS0) -o $@ $(MODLIB) $@.$(SUF) $(OPTLIB) $(LIBS); \
             $@ < $@.dat > $@.res; \
             'rm' -f $@
     
     YY = $(FC) $(OPTS0) -o $@ $(MODLIB) $@.$(SUF) $(OPTLIB) $(LIBS)
     
     .SUFFIXES: .f90 .f .o
     
     .$(SUF).o: 
     	$(FC) $(OPTS3) $<
     
     .f.o:
     	$(FC1) $(OPTS3) $<

6.1.3 Example of program using LAPACK95

This source should be saved in a file called ejemplo_la_spsv.f90 to use the Makefile provided in n makefile for compiling programs with calls to LAPACK95, Section 6.1.4.

     PROGRAM LA_SSPSV_EXAMPLE
     
       !  -- LAPACK95 EXAMPLE DRIVER ROUTINE (VERSION 1.0) --
       !
       !  .Use Statements
       USE LA_PRECISION, ONLY: WP => SP
       USE F95_LAPACK, ONLY: LA_SPSV
       !  Implicit Statement
       IMPLICIT NONE
       !  Local Scalars
       INTEGER :: I, N, NN, NRHS
       !  Local Arrays
       INTEGER, ALLOCATABLE :: IPIV(:)
       REAL(WP), ALLOCATABLE :: B(:,:), AP(:)
       !  Executable Statements
       WRITE (*,*) 'SSPSV Example Program Results.'
       N = 5; NRHS = 1
       WRITE(*,'(5H N = , I4, 9H; NRHS = , I4)') N, NRHS
       NN = N*(N+1)/2
       ALLOCATE ( AP(NN), B(N,NRHS), IPIV(N) )
       !
       CALL  RANDOM_NUMBER(AP)
       AP = 10*AP
       !
       WRITE(*,*)'Matrix AP :'
       DO I=1,NN; WRITE(*,"(15(I3,1X,1X),I3,1X)") INT(AP(I));
       ENDDO
       !
       WRITE(*,*)'Matrix B :'
       CALL  RANDOM_NUMBER(B)
       B = 10*B
       DO I=1,N; WRITE(*,"(10(I3,1X,1X),I3,1X)") INT(B(I,1));
       ENDDO
       !
       WRITE(*,*)" CALL LA_SPSV( AP, B, 'L', IPIV )"
       !
       CALL LA_SPSV( AP, B, 'L', IPIV )
       !
       WRITE(*,*)'AP on exit: '
       DO I=1,NN; WRITE(*,"(15(E13.5))") AP(I); 
       ENDDO
       !
       WRITE(*,*)'Matrix B on exit :'
       DO I=1,N; WRITE(*,"(F9.5)") B(I,1);
       ENDDO
       !
       WRITE(*,*)'IPIV = ', IPIV
       !
     END PROGRAM LA_SSPSV_EXAMPLE

6.1.4 makefile for compiling programs with calls to LAPACK95

     #
     #  -- LAPACK95 makefile (version 1.0) --
     #
     FC	 = gfortran
     #
     MODLIB   = -I/usr/local/lib/lapack95_modules
     OPTS1    = -c
     OPTS3    = $(OPTS1) $(MODLIB)
     OPTL     = -o 
     OPTLIB   = -lblas -llapack
     
     LAPACK_PATH = /usr/local/lib
     LAPACK95_PATH = /usr/local/lib
     
     LAPACK95 = $(LAPACK95_PATH)/lapack95.a
     
     LIBS     = $(LAPACK95)
     SUF      = f90
     
     YY = $(FC) -o $@ $(MODLIB) $@.$(SUF) $(OPTLIB) $(LIBS)
     
     .SUFFIXES: .f90 .f .o
     
     .$(SUF).o: 
     	$(FC) $(OPTS3) $<
     
     ejemplo_la_spsv: 
     	$(YY)
     
     clean:
     	'rm' -f *.o *.mod core

6.2 Compile and link statically with NAG and LAPACK

It is important to be able to compile and link statically programs when libraries are not available in all nodes. This is the case with the NAG library which is not compatible with gfortran, the only Fortran compiler in Debian Lenny. The program statically linked in one node (where Etch is installed and g77 is available can the be executed in any other node.

An example of compilation is the following

      
     g77 -static -o infsq_box_1D  infsq_box_1D.f ../Potentials/wsaxon_Box_pot.f -L /usr/lib/atlas -L/usr/local2/NAG -lnag  -llapack -lblas-3

Some important points:

  1. It is necessary to include both Lapack and Blas libraries.

  1. The Blas library should be blas-3.

  1. The use of standard Lapack and Blas libraries give an error due to the different sizes of object files. Something like

          
         (xerbla.o): In function `xerbla_': multiple definition of `xerbla_'
         /usr/lib/liblapack.a(xerbla.o): first defined here
         /usr/bin/ld: Warning: size of symbol `xerbla_' changed from 86
         in /usr/lib/liblapack.a(xerbla.o) to 38
         in /usr/lib/libblas.a(xerbla.o)
         collect2: ld returned 1 exit status
    

    This is a known bug and can be solved using the libraries provided with the Atlas packages and adding the corresponding path to the compilation: -L /usr/lib/atlas.

  1. The order of the libraries is not irrelevant. In particular I found that lapack has to be invoked prior to blas-3.


6.3 Creating a Makefile to use cpp with gfortran

Added on November 09th, 2015.

There is a way of using the C language preprocessor, cpp, with gfortran. In this way we can communicate via Makefile with our program. Let's assume that we have a Fortran 90 program to compute two normally distributed set of random numbers that use the Box-Muller transform. The program is called box_muller.f90 and we also compute the mean, median and standard variance of the resulting distributions using a subroutine called stats.f90[13]. The FORTRAN 90 codes are provided in Fortran Codes, Section 6.3.1.

The Makefile used for compilation is

      
     # Fortran Compiler
     FC=gfortran
     #
     # C Preprocessor
     CPP=cpp
     #
     TYPE_DEF=-DPREC=DP 
     ifdef USE_SINGLE
     TYPE_DEF = -DPREC=SP 
     endif
     #
     # Optimization Flags
     OPTFLAGS = -O3 -funroll-loops -march=native
     OPTFLAGS_P = "\"-O3 -funroll-loops -march=native\""
     #
     all: box_muller
     #
     .FORCE:
     #
     box_muller: box_muller.f90 stats.f90 .FORCE
     	$(CPP)  -std=c89 $(TYPE_DEF) -DFLAGS=$(OPTFLAGS_P) $< > /tmp/$<
     	$(CPP)  -std=c89 $(TYPE_DEF)  stats.f90 > /tmp/stats.f90
     	$(FC) $(OPTFLAGS) -o $@ /tmp/stats.f90 /tmp/$< 
     	rm /tmp/stats.f90 /tmp/$<
     #
     clean:
     	rm -f *.o box_muller

Note how cpp is called twice, once for each file and the PREC string are replaced either by DP or SP, the last case if the option USE_SINGLE is set. In the main program the string FLAGS is replaced by the contents of the variable OPTFLAGS_P to print the compilation flags used. Thus to run the program in single precision the user proceed as follows

      
     $ make all USE_SINGLE=1
     cpp  -std=c89 -DPREC=SP  -DFLAGS="\"-O3 -funroll-loops -march=native\"" box_muller.f90 > /tmp/box_muller.f90
     cpp  -std=c89 -DPREC=SP   stats.f90 > /tmp/stats.f90
     gfortran -O3 -funroll-loops -march=native -o box_muller /tmp/stats.f90 /tmp/box_muller.f90 
     rm /tmp/stats.f90 /tmp/box_muller.f90
     $ ./box_muller 
     10000
      Time Box Muller subroutine :   4.00000019E-03
       MEAN =   -7.31508713E-03
       STANDARD DEVIATION =   0.995381057    
       MEDIAN IS =    1.01202750E-04
      Time STATS subroutine 1:  0.100005999    
       MEAN =    2.96602229E-04
       STANDARD DEVIATION =   0.999957144    
       MEDIAN IS =    1.29099786E-02
      Time STATS subroutine 2:  0.100006007    
      KIND = SINGLE ;  Optimization Flags: -O3 -funroll-loops -march=native

6.3.1 Fortran Codes

Main program: box_muller.f90

      
     PROGRAM Box_Muller_Prog
       !
       IMPLICIT NONE
       !
       ! Single and double precision kind parameters
       INTEGER, PARAMETER :: SP = KIND(1.0)
       INTEGER, PARAMETER :: DP = KIND(1.0D0)
       !
       INTEGER :: I, IERR
       REAL(KIND = PREC), DIMENSION(:), ALLOCATABLE :: X, Y
       REAL(KIND = PREC) :: M, SD, MEDIAN
       REAL(KIND = SP) :: time_end, time_start
       !
       !
       CHARACTER(LEN = 6) :: PRECISION 
       !
       ! interface block   
       INTERFACE
          SUBROUTINE STATS(VECTOR,N,MEAN,STD_DEV,MEDIAN)
            IMPLICIT NONE
            ! Single and double precision kind parameters
            INTEGER, PARAMETER :: SP = KIND(1.0)
            INTEGER, PARAMETER :: DP = KIND(1.0D0)
            INTEGER , INTENT(IN)                    ::  N
            REAL(KIND = PREC)      , INTENT(IN) , DIMENSION(:)   :: VECTOR  
            REAL(KIND = PREC)      , INTENT(OUT)                 :: MEAN
            REAL(KIND = PREC)      , INTENT(OUT)                 :: STD_DEV
            REAL(KIND = PREC)      , INTENT(OUT)                 :: MEDIAN
          END SUBROUTINE STATS
       END INTERFACE
       !
       ! Set Precision for output
       IF (PREC == DP) THEN
          PRECISION = "DOUBLE"
       ELSE
          PRECISION = "SINGLE"
       ENDIF
       !
       ! Input number of data
       READ*, I  
       !
       ALLOCATE(X(1:I), STAT = IERR)    
       IF (IERR /= 0) THEN
          PRINT*, "X allocation request denied."
          STOP
       ENDIF
       !
       ALLOCATE(Y(1:I), STAT = IERR)    
       IF (IERR /= 0) THEN
          PRINT*, "Y allocation request denied."
          STOP
       ENDIF
       !
       ! Set time start
       CALL CPU_TIME(time_start)
       !
       CALL BOX_MULLER(I)
       !
       ! Set time end
       CALL CPU_TIME(time_end)
       !
       !PRINT*, X
       !PRINT*, Y
       !
       PRINT*, "Time Box Muller subroutine :", time_end - time_start
       time_start = time_end
       !
       CALL STATS(X,I,M,SD,MEDIAN)
       !
       PRINT *,' MEAN = ',M
       PRINT *,' STANDARD DEVIATION = ',SD
       PRINT *,' MEDIAN IS = ',MEDIAN
       !
       ! Set time end
       CALL CPU_TIME(time_end)
       !
       PRINT*, "Time STATS subroutine 1:", time_end - time_start
       time_start = time_end
       !
       IF (ALLOCATED(X)) DEALLOCATE(X, STAT = IERR) 
       IF (IERR /= 0) THEN
          PRINT*, "X NON DEALLOCATED!"
          STOP
       ENDIF
       !
       CALL STATS(Y,I,M,SD,MEDIAN)
       !
       !
       PRINT *,' MEAN = ',M
       PRINT *,' STANDARD DEVIATION = ',SD
       PRINT *,' MEDIAN IS = ',MEDIAN
       !
       ! Set time end
       CALL CPU_TIME(time_end)
       !
       PRINT*, "Time STATS subroutine 2:", time_end - time_start
       time_start = time_end
       !
       IF (ALLOCATED(Y)) DEALLOCATE(Y, STAT = IERR)   
       IF (IERR /= 0) THEN
          PRINT*, "Y NON DEALLOCATED!"
          STOP
       ENDIF
       !
       PRINT*, "KIND = ", PRECISION, " ;  Optimization Flags: ", FLAGS 
       !
     CONTAINS
       !
       SUBROUTINE BOX_MULLER(dim)
         ! 
         ! Uses the Box-Muller method to create two normally distributed vectors
         !
         INTEGER, INTENT(IN) :: dim
         !
         REAL(KIND = PREC), PARAMETER :: PI = ACOS(-1.0)
         REAL(KIND = PREC), DIMENSION(dim) :: RANDOM_u, RANDOM_v ! Automatic arrays
         !
         CALL RANDOM_NUMBER(RANDOM_u)
         CALL RANDOM_NUMBER(RANDOM_v)
         !
         X = SQRT(-2.0_ PREC*LOG(RANDOM_u))
         Y = X*SIN(2*PI*RANDOM_v)
         X = X*COS(2*PI*RANDOM_v)
         !
       END SUBROUTINE BOX_MULLER
       !
     END PROGRAM Box_Muller_Prog

Subroutine: stats.f90

      
     SUBROUTINE STATS(VECTOR,N,MEAN,STD_DEV,MEDIAN)
       IMPLICIT NONE
       ! Single and double precision kind parameters
       INTEGER, PARAMETER :: SP = KIND(1.0)
       INTEGER, PARAMETER :: DP = KIND(1.0D0)
       ! Argument definition
       INTEGER , INTENT(IN)                    ::  N
       REAL(KIND = PREC)      , INTENT(IN) , DIMENSION(:)    ::  VECTOR   
       REAL(KIND = PREC)      , INTENT(OUT)                  ::  MEAN
       REAL(KIND = PREC)      , INTENT(OUT)                  ::  STD_DEV
       REAL(KIND = PREC)      , INTENT(OUT)                  ::  MEDIAN
       ! Local variables
       REAL(KIND = PREC) :: VARIANCE = 0.0
       REAL(KIND = PREC)      :: SUMXI = 0.0, SUMXI2 = 0.0
       REAL(KIND = PREC)      , DIMENSION(1:N)              ::  Y
       !
       SUMXI=SUM(VECTOR)      
       SUMXI2=SUM(VECTOR*VECTOR) 
       MEAN=SUMXI/N       
       VARIANCE=(SUMXI2-SUMXI*SUMXI/N)/(N-1)
       STD_DEV = SQRT(VARIANCE)
       Y=VECTOR
       ! Sort values
       CALL SELECTION
       IF (MOD(N,2) == 0) THEN
          MEDIAN=(Y(N/2)+Y((N/2)+1))/2
       ELSE
          MEDIAN=Y((N/2)+1)
       ENDIF
     CONTAINS   
       SUBROUTINE SELECTION
         IMPLICIT NONE
         INTEGER :: I,J,K
         REAL :: MINIMUM
         DO I=1,N-1
            K=I
            MINIMUM=Y(I)
            DO J=I+1,N
               IF (Y(J) < MINIMUM) THEN
                  K=J
                  MINIMUM=Y(K)
               END IF
            END DO
            Y(K)=Y(I)
            Y(I)=MINIMUM
         END DO
       END SUBROUTINE SELECTION
     END SUBROUTINE STATS

6.3.2 References

  1. StackOverflow


6.4 Obtaining environment variables values from Fortran

Created on March 9th, 2018

With the gfortran compiler the intrinsic procedure GET_ENVIRONMENT_VARIABLE allows the reading of environment variables such as hostname, username, or . The syntax for this procedure is

           CALL GET_ENVIRONMENT_VARIABLE(NAME[, VALUE, LENGTH, STATUS, TRIM_NAME])

This procedure replaces the now deprecated GETENV procedure. In the next program we use this subroutine to get values of the environment variables HOSTNAME, USER, and HOME. Note that unless the user has exported this variable HOSTNAME reading may fail. We also use the intrinsic procedure GETPID to obtain the process PID number.

     PROGRAM GET_ENV
       !                                           
       IMPLICIT NONE
       !
       INTEGER :: ivar = 0, len, stat, pid
       CHARACTER(LEN=256) :: home, host, user
       !
       ! Hostname variable
       CALL GET_ENVIRONMENT_VARIABLE('HOSTNAME', host, len, stat)
       IF (stat == 0) THEN
          PRINT*, "Hostname read: ", TRIM(host)
       ELSE
          PRINT*, "Hostname read failed: stat = ", stat
       ENDIF
       !
       ! Username variable
       CALL GET_ENVIRONMENT_VARIABLE('USER', user, len, stat)
       IF (stat == 0) print*, "Username read: ", TRIM(user)
       !
       ! Homedir variable
       CALL GET_ENVIRONMENT_VARIABLE('HOME', home, len, stat)
       IF (stat == 0) THEN
          PRINT*, "Homedir read: ", TRIM(home)
       ELSE
          PRINT*, "Homedir read failed: stat = ", stat
       ENDIF
       !
       ! Process ID
       pid = GETPID()
       PRINT*, "Process ID: ", pid
       !
     END PROGRAM GET_ENV

A possible output of the previous code is

     $ ./get_env_gfortran 
      Hostname read failed: stat =            1
      Username read: bilbo
      Homedir read: /home/bilbo
      Process ID:         5614

6.5 Creating a vector of random data

Created on March 8th, 2018

The best way to deal with random numbers with gfortran is making use of the intrinsic procedures RANDOM_NUMBER and RANDOM_SEED. Note that RAN and RAND procedures have been deprecated. From the gfortran documentation:

           Returns a single pseudorandom number or an array of pseudorandom numbers from the uniform distribution over the range 0 \leq x \le 1.

Therefore the creation of a vector with, e.g. 100 random integers in the 0-99 range should be accomplished with a simple program as this one

     PROGRAM RANDOM_VECTOR
       IMPLICIT NONE
       REAL, DIMENSION(1:100) :: REAL_RAND_VALUES
       INTEGER, DIMENSION(1:100) :: INT_RAND_VALUES
       INTEGER :: INDEX
       !
       CALL RANDOM_NUMBER(REAL_RAND_VALUES)
       !
       INT_RAND_VALUES = NINT(REAL_RAND_VALUES*100)
       !
       DO INDEX = 1, 100
          PRINT*, INDEX, INT_RAND_VALUES(INDEX)
       ENDDO
       !
     END PROGRAM RANDOM_VECTOR

And everything seems pretty okay until we run several times the program and observe that our random number vector always have the same values... Surprise...

                1         100
                2          57
                3          97
               ..............
               98          11
               99          22
              100          10

The problem is that RANDOM_NUMBER provides a sequence of pseudorandom values. If it starts the sequence from the same point, it gives you always the same sequence. We need to chose a different starting point every time that we use the procedure. This can be done using RANDOM_SEED and providing a particular value to this procedure, e.g. the process PID. This is done in the code that follows

     PROGRAM RANDOM_VECTOR
       IMPLICIT NONE
       REAL, DIMENSION(1:100) :: REAL_RAND_VALUES
       INTEGER, DIMENSION(1:100) :: INT_RAND_VALUES
       INTEGER, DIMENSION(:), ALLOCATABLE :: pid
       INTEGER :: index, size=1
       !
       ! Get minimum size of the PUT array for RANDOM_SEED
       CALL RANDOM_SEED(SIZE=size)
       !
       ! Allocate and give values to pid array
       ALLOCATE(pid(1:size))
       pid = GETPID() 
       !
       ! Restart random seed
       CALL RANDOM_SEED(PUT=pid)
       !
       ! Generate random vector
       CALL RANDOM_NUMBER(REAL_RAND_VALUES)
       !
       ! Transform to integers
       INT_RAND_VALUES = NINT(REAL_RAND_VALUES*100)
       !
       DO index = 1, 100
          PRINT*, index, INT_RAND_VALUES(index)
       ENDDO
       !
     END PROGRAM RANDOM_VECTOR

6.6 Obtaining the values and positions of the larger values in an array

Created on July 30th, 2018

We use the random vector building method of Creating a vector of random data, Section 6.5 to build an array of random values and then show how to obtain the indexes and values of the largest elements in each row or column of the array. In order to do so we use the intrinsic procedures MAXVAL and MAXLOC.

Therefore we first create a two dimensional array with 100 random integers in the 0-1 range, display the array and then calculate the maximum values and its index in a column-wise and a row-wise way.

     PROGRAM MAXIMA_VALUES_POSITIONS
       IMPLICIT NONE
       INTEGER, PARAMETER :: dim = 10
       REAL, DIMENSION(1:dim,1:dim) :: ARRAY_VALUES
       REAL, DIMENSION(1:dim) :: MAX_VALUES
       INTEGER, DIMENSION(1:dim) :: MAX_INDEX
       INTEGER, DIMENSION(:), ALLOCATABLE :: pid
       INTEGER :: index, size=1
       !
       ! Get minimum size of the PUT array for RANDOM_SEED
       CALL RANDOM_SEED(SIZE=size)
       !
       ! Allocate and give values to pid array
       ALLOCATE(pid(1:size))
       pid = GETPID() 
       !
       ! Restart random seed
       CALL RANDOM_SEED(PUT=pid)
       !
       ! Generate random vector
       CALL RANDOM_NUMBER(ARRAY_VALUES)
       !
       ! Display array
       DO index = 1, dim
          PRINT 10, index, ARRAY_VALUES(index,:)
       ENDDO
       !
       !
       PRINT*
       !
       !
       PRINT*, 'MAXIMUM COLUMN-WISE VALUES'
       MAX_VALUES = MAXVAL(ARRAY_VALUES, 1)
       MAX_INDEX = MAXLOC(ARRAY_VALUES,1)
       DO index = 1, dim
          PRINT 20, index, MAX_INDEX(index), MAX_VALUES(index)
       ENDDO
       !
       PRINT*
       !
       PRINT*, 'MAXIMUM ROW-WISE VALUES'
       MAX_VALUES = MAXVAL(ARRAY_VALUES, 2)
       MAX_INDEX = MAXLOC(ARRAY_VALUES, 2)
       DO index = 1, dim
          PRINT 20, index, MAX_INDEX(index), MAX_VALUES(index)
       ENDDO
       !
     10 FORMAT(1X, I2, 1X, *(F6.4, 1X)) ! Note the * in the format statement for variable number of fields
       !
     20 FORMAT(2X, I2, 1X, I2, 2X, F6.4)
     END PROGRAM MAXIMA_VALUES_POSITIONS

Notice that if the second argument to the procedures is one the array is sliced column-wise and if it is two it is sliced row wise. There are analogous procedures for minimum values with the intrinsics MINVAL and MINLOC.


6.7 Tips for building a Fortran Makefile

Created on March 20th, 2018

I enclose a template of a Fortran 90 Makefile that solves some common problems that appear when using the make application with this programming language. Let's assume that we have a program file called prog_main.f90 that depends on two module files, module_a.f90 and module_b.f90. And we want to use the preprocessor and to build single thread and multithread execs making use of openmp. The template Makefile is as follows:

     BINPATH = ../bin/
     ###################################################
     prog_SRC = module_a.f90 module_b.f90 prog_main.f90 
     #
     prog_OBJ = $(prog_SRC:.f90=.o)
     #
     prog_OPENMP_OBJ = $(prog_SRC:.f90=.oOMP)
     ###################################################
     #
     .SUFFIXES:
     .SUFFIXES: .o .f90 .oOMP
     ###
     FC = gfortran
     FOPT	= -O3 -Wall
     PREP = -cpp
     #
     MODLIB   = -I/usr/local/include/lapack95_modules
     LAPACK95 = -L/usr/local/lib -llapack95
     LAPACK77 = -llapack
     BLAS     = -lblas
     # 
     LIBS     = $(LAPACK95)  $(LAPACK77) $(BLAS) 
     ###################################################
     all: prog prog_openmp 
     ###################################################
     .PHONY : all
     ###################################################
     .f90.o:
     	$(info )
     	$(info Compiling single thread object file:)
     	$(FC) -c $(PREP) $(MODLIB)  "$<"
     ###################################################
     .f90.oOMP:
     	$(info )
     	$(info Compiling multi thread object file:)
     	$(FC) -c -fopenmp $(PREP) $(MODLIB)  -o "$@" "$<"
     ###################################################
     prog: $(prog_OBJ) Makefile 
     	$(info )
     	$(info Linking prog single thread executable:)
     	$(FC) $(FOPT) $(MODLIB) -o $(BINPATH)/$@_$(FC) $(prog_OBJ) $(LIBS) 
     #######################
     prog_openmp: $(prog_OPENMP_OBJ) Makefile 
     	$(info )
     	$(info Linking prog multithread executable:)
     	$(FC) -fopenmp $(FOPT) $(MODLIB) -o $(BINPATH)/$@_$(FC) $(prog_OPENMP_OBJ) $(LIBS) 
     #######################
     clean:
     	$(info Cleaning object, module, and executable files.)
     	@rm -f  *.o *.oOMP *.mod $(BINPATH)/prog_$(FC) $(BINPATH)/prog_openmp_$(FC)
     #######################

Notice, for example, the empty SUFFIXES line, as it resets the predefined suffixes to avoid dependency problems associated with the .mod module files. The executable files are called prog_gfortran and prog_openmp_gfortran and are stored in the folder defined by the BINPATH variable.


[ previous ] [ Contents ] [ 1 ] [ 2 ] [ 3 ] [ 4 ] [ 5 ] [ 6 ] [ 7 ] [ 8 ] [ 9 ] [ 10 ] [ 11 ] [ next ]


Some Mini-Howtos of Interest

Curro Perez-Bernal mailto:francisco.perez@dfaie.uhu.es