compilationfortrangfortranblasf2py

Blas function return 2 different result when compiled with relative path and absolute path


I have a very weird problem when I compile Fortran code using the zgerc subroutine from BLAS. Basically, this subroutine calculate the outer product of vector x with the conjugate of vector y. More about that function here. My simple code is as following:

program main
    implicit none
    integer :: i
    complex(8), dimension(10) :: a = [(i, i=0,9)]
    complex(8), dimension(10) :: b = [(i, i=0,9)]
    complex(8), dimension(10, 10) :: c
    c = 0
    CALL zgerc(10, 10, 1.D0, a, 1, b, 1, c, 10)
    WRITE(*, *) c
end program main

I have here 2 complex vectors, a and b, both goes from 0 to 9 and their imaginary part is 0.

Now for the weird part. If I compile that code with absolute path: gfortran -c /home/myUser/Fortran/tests/main.f90 -o main.o I get the correct result, but if I compile with gfortran -c main.f90 -o main.o (of course I'm in the right directory, and I've also tried ./main.f90) the result of the real part is correct, but for the imaginary part I get numbers like 1E+225 (and if I use ./main.f90 I get numbers like 1E+163).

I can't understand why the path to my code change the result of the imaginary part... I'll be glad for your help.

I use Ubuntu 20.04.2 with the default gfortran (9.3.0)

P.S, my final goal is to use this as part of a more complex subroutine in Python with f2py.

EDIT: my full commands:

#gfortran -c /home/myUser/Fortran/tests/main.f90 -o main.o
gfortran -c main.f90 -o main.o
gfortran -o test main.o /home/myUser/PycharmProjects/GSIE_2D/fortran_scripts/libblas.a /home/myUser/PycharmProjects/GSIE_2D/fortran_scripts/liblapack.a
rm ./main.o
./test

line 1 and 2 are the 2 cases, so I run only one of them each time.


Solution

  • You supply 1d0 which is a double precision literal whereas zgerc assumes a double complex value.

    call zgerc(10, 10, cmplx(1, kind=8), a, 1, b, 1, c, 10)
    

    By including explicit interfaces (via some kind of blas module) you would have gotten compile time errors when supplying arguments of wrong data type. Intel's mkl provides such explicit interfaces in its blas95 module as well as generic routines (gerc instead of {c,z}gerc). There is also this open-source module providing explicit interfaces to standard blas routines.

    Also make use of the portable types defined in iso_fortran_env.

    program main
      use blas95,          only: gerc
      use iso_fortran_env, only: real64
      implicit none
    
      integer, parameter :: n = 10
      integer            :: i
      complex(real64)    :: a(n) = [(i, i=0,n-1)], b(n) = [(i, i=0,n-1)], c(n,n)
    
      c = 0
      ! call zgerc(10, 10, cmplx(1, kind=8), a, 1, b, 1, c, 10) ! standard blas 
      call gerc(c, a, b, alpha=cmplx(1, kind=real64))           !  generic blas95 
      print *, c
    end program