blas compilation issue

Accelerated Computing HPC Compilers Legacy PGI Compilers
#1

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:

program zgemm_wrap
implicit none
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.

final remarks:
(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!

#2

Your Zgemm call is attempting to do the operation

B = AA + betaB

The first error is that the multiplier beta is required to be double complex; you are passing the double precision real constant 1d0

The second error is that when beta is not zero, B has to be set before the call. You have left B undefined.

x and y are not used in the calculation of B, so their values are immaterial.

#3

Wow! Thank you so much. Of course, you are correct. I was worried my issue was in calling BLAS, which is confirmed. I was just so confused because initializing y seemed to solve the problem, as did other random things.
Interesting the issue didn’t arise on my university cluster, or when I invoked the -C option on the pgi compiler. I naively thought the variables would be initialized to zero; lesson learned the hard way…
Again, thank you, I’m quite new to programming and have found it difficult. I’ll have to be more careful.

#4

It is quite straightforward to check if you have the arguments correctly arranged and initialized, although it does take a bit of reading and thinking during the first consultation.

For Zgemm (and, similarly, any BLAS routine), the documentation at

http://www.netlib.org/lapack/explore-html/d7/d76/zgemm_8f.html

lists the names of the arguments, their types, whether they are IN, OUT, or INOUT arguments. and further details if relevant (such as when B has to be set if Beta is not zero).