when to use iso_Fortran_env ,selected_int_kind,real(8),or -fdefault-real-8 for writing or compiling fortran code?

Question

I have this simple code which uses DGEMM routine for matrix multiplication

program check

  implicit none
  real(8),dimension(2,2)::A,B,C
  A(1,1)=4.5
  A(1,2)=4.5
  A(2,1)=4.5
  A(2,2)=4.5
  B(1,1)=2.5
  B(1,2)=2.5
  B(2,1)=2.5
  B(2,2)=2.5

  c=0.0

  call DGEMM('n','n',2,2,2,1.00,A,2,B,2,0.00,C,2)
  print *,C(1,1)    
  print *,C(1,2)  

  print *,C(2,1)    
  print *,C(2,2)
end program check

now when i compile this code with command

gfortran -o check check.f90   -lblas 

I get some random garbage values. But when I add

-fdefault-real-8

to the compiling options I get correct values.

But since it is not a good way of variable declaration in Fortran. So I used the iso_fortran_env intrinsic module and added two lines to the code

use iso_fortran_env
real(kind=real32),dimension(2,2)::A,B,C

and compiled with

 gfortran -o check check.f90 -lblas

Again I got wrong output . Where I'm erring in this code? I'm on 32bit linux and using GCC

Answer 1

DGEMM expects double precision values for ALPHA and BETA.

Without further options, you are feeding single precision floats to LAPACK - hence the garbage. Using -fdefault-real-8 you force every float specified to be double precision by default, and DGEMM is fed correctly.

In your case, the call should be:

call DGEMM('n','n',2,2,2,1.00_8,A,2,B,2,0.00_8,C,2)

which specifies the value for alpha to be 1 as a float of kind 8, and zero of kind 8 for beta.

If you want to perform the matrix-vector product in single precision, use SGEMM.

Note that this is highly compiler-specific, you should consider using REAL32/REAL64 from the ISO_Fortran_env module instead (also for the declaration of A, B, and C).