I'm trying to understand a bug that I encountered in my code. When trying to allocate an array using the transpose of an array, I use an allocate statement with source=transpose(original_array) in my code. However, I don't get the expected result using this method. It seems that the indexing is off by one and the first row of the source array is skipped.
Example:
program testalloc
real*8, allocatable :: a(:, :)
real*8, allocatable :: b(:, :)
allocate(b(2, 3))
b(1, :) = [1, 2, 3]
b(2, :) = [4, 5, 6]
call printmat(b)
a = transpose(b)
call printmat(a) ! Good
deallocate(a)
allocate(a(3, 2), source=transpose(b))
call printmat(a) ! Bad
deallocate(a)
allocate(a(3, 2))
a = transpose(b)
call printmat(a) ! Good
contains
subroutine printmat(mat)
real*8, intent(in) :: mat(:, :)
integer :: i
write(*,*) 'print'
do i = 1, size(mat, 1)
write(*,*) mat(i, :)
end do
end subroutine
end program
which gives
print
5.0000000000000000 3.0000000000000000
6.0000000000000000 0.0000000000000000
3.2114266979681025E-322 5.0000000000000000
for the sourced allocation after compiling with gfortran (gcc version 7.3.0 (Ubuntu 7.3.0-27ubuntu1~18.04)) instead of the transposed original array. Am I doing something wrong here or is it a compiler bug?
The sourced allocation
allocate(a(3, 2), source=transpose(b))
is a valid one. The shape specified for a
is the same as that of the source expression transpose(b)
. As a result, a
takes the value of the expression given.
It is incorrect for the compiler to give a different result. The output routine is not to blame here.
gfortran 8 appears to give the expected output.
Interestingly, with gfortran 7 the expected result appears if b
is not itself allocatable.
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With