I am facing a problem with a code that is related to the status of a pointer after a call to the BLAS subroutine zcopy.
The very simple code is the following:
program test integer, parameter :: n = 10, m = 21 ! if m=11 the code works real(8) :: A(n,m) real(8), pointer :: C(:,:) ! initialize A A = 0.0 A(1,2) = 4.5 allocate(C(1:m,1:n)) C = transpose(A) ! the call to zcopy should reassign A to its transpose. call zcopy(n*m, C,1,A,1) ! deallocation here is not succesful deallocate(C) end program test
The compiler 17.04 fails to deallocate C, complaining that it is not allocated. If I change m from 21 to 11, the programs doesn't fail and the deallocation is succesful. I have a linux system with Ubuntu 16.04.
Could anyone help me to understand what rules am I breaking about the status of pointers?
Your arguments to the ZCOPY call are inconsistent in type (or you have the number of elements wrong) with what ZCOPY wants - you are providing REAL(8) - it expects COMPLEX*16. As a result, execution of ZCOPY trashes memory, all sorts of bad things may then happen.
Did you intend to call DCOPY?
The Fortran 90 RESHAPE intrinsic might be of interest.
The BLAS routines use implicit interfaces, and do not distinguish between pointer and non-pointer arguments. Originally, they were written in Fortran 77; even if they have been rewritten, whether in Fortran or another language, they have to behave in the same way as Fortran 77 subroutines and functions.
In the specific code that you showed, the error is the programmer's. The routine ZCOPY is for copying double precision complex arrays. Your call passes double precision real arrays. Instead of copying m.n double precision numbers, the routine copies 2.m.n double precision numbers, which clobbers memory that is beyond that allocated for A.
Why did you call ZCOPY rather than DCOPY? Unless you are writing Fortran 77 code, why call an external routine?
Sorry, it was zcopy, but with complex(8) not real(8).
You are right about the abuse of this subroutines, but this is an extract of a code that is not mine. I think the best option for me is just to rewrite part of the code instead of debugging the use of zcopy.
I cannot simply write A = C as they have differente sizes.
Sorry for this comment. That part of the code is working fine. The problem with zcopy is still there in my full code but it is quite difficult to define where the error come from.
No apology needed. You are probably at a fork in the road that many of us have faced. You can either ask the original author of the code for help/guidance, or try to craft your own remedy. I appreciate your attempt to create a short example to illustrate the problem with the big code. We have just seen that there is such a thing as oversimplification.
I suggest that you try to isolate the problem. First, use the -traceback option when compiling, and see if the failure to deallocate is repeatable and locate the part of the code where that happens.
Second, check the arguments that were passed to the subprogram. Add WRITE statements, dump to a file, or use a symbolic debugger.
Finally, consider writing a short driver that passes the same arguments to the problematic subprogram.