In the code below I am having the following issue which arose while trying to write a wrapper (which I’m still struggling with and could use some help) to zgemm in BLAS so that I can use it more easily:
integer, parameter :: N = 2
double complex, dimension(N, N) :: A, B
real, dimension(N, N) :: x
real, dimension(2) :: y
y = (/ 1.0, 2.0 /)
A = reshape((/ 1.0, 2.0, 3.0, 4.0 /), (/ 2, 2/))
call ZGEMM(‘N’, ‘N’, N, N, N, 1.0d0, A, N, A, N, 1.0d0, B, N)
write(*, *) B
write(*, *) matmul(A, A)
end program zgemm_wrap
If I comment out line 8 (i.e. the initialization of y) the results from zgemm are correct, namely (7, 10, 15, 22).
However, if I uncomment it, the results are incorrect; (7 + i14, 10 + i20, 15 + i30, 22 + i44). That is, the imaginary component becomes non-zero (which it shouldn’t) and curiously has a value twice that of the real component. I’ve had a similar problem in the past where the presence of a write statement fixed or destroyed the calculations
I’m using PGI 19.10 community edition on windows 10
I’m not sure what elso to include regarding my machine.
(1) Compiling with gcc (and pgi) on a linux machine at my university does not reproduce the error
(2) Running with the option -C for the pgi compiler seems to fix the error, so I suspect it has to do with array bounds.
(3) Due to the above I initially thought it may be an issue with BLAS on my machine, since I am just using what came with the pgi program.
(4) I just realized that writing out the variable y at the end of the program solved the issue, but I have no idea why. Curiously, it doesn’t fix the problem if I wrote it out at the beginning. I suspect the issue lies with y
(5) Just to be sure, I compile as: pgf90 zgemm_wrap.f90 -lblas, and then run the following executable.
thanks for any help!