I'm writing a Fortran 95 program to generate random correlated surfaces. The code worked in a "program" file, but I need to put it in a module. For some reason, the sorting algorithm I use goes out of bounds if I call it inside the while loop after I put it in a module. If I call it outside the while loop, or if I call it on a vector of the same size (xo in the code) it stays in bounds. I tried calling it on a copy of the vector I want to sort, but it failed. I really have no clue why, does anyone have an idea?
I am working in Ubuntu 13.10 with the gfortran 4.7.3 compiler
Code with variables and while loop:
subroutine generate_skewed_surface(x, system_length, corr_length, alpha)
real(wp), intent(inout), allocatable :: x(:)
real(wp), intent(in) :: system_length, corr_length, alpha
real(wp) :: deltat, deltao
real(wp) :: prefactor, factor
integer :: dim, i
integer, allocatable :: indx(:), indxp(:)
real(wp), allocatable :: omega_k(:), Sxx(:), xo(:), Rp(:)
complex(wp), allocatable :: Xf(:), Rk(:)
logical :: converged
if(.not. allocated(x)) then
!main should be alerted
return
end if
[Skipping initializations and intermediate steps. The only thing that happens to x is that it is filled with random numbers]
!store a copy
xo = x
!sort original signal
call QsortC(xo) !<- this call works
!create surface by reordering original signal
converged = .false.
do while (.not. converged)
Rk = cmplx(x, 0._wp, wp)
call fft(Rk)
Rp = atan2(real(-imu * Rk, wp), real(Rk, wp))
Rk = exp(imu*Rp)*abs(Xf)
call ifft(Rk)
x = real(Rk, wp)
indx = (/(i, i=1, dim)/)
call QsortC(x, indx) !<- this call produces the error during the first call
do i=1, dim
x(indx(i)) = xo(i)
end do
converged = .true.
do i=1, dim
if (indx(i) .ne. indxp(i)) then
converged = .false.
end if
end do
indxp = indx
end do
Sorting algorithm code (separate module file):
module module_qsort
use SFL_precision, only : DP, SP
implicit none
public :: QsortC
private :: Partition
integer, parameter, private :: wp = DP
contains
recursive subroutine QsortC(A, indices)
real(wp), intent(in out), dimension(:) :: A
integer, optional, intent(inout), dimension(:) :: indices
integer :: iq
if(size(A) > 1) then
if (present(indices)) then
call Partition(A, iq, indices)
call QsortC(A(:iq-1), indices(:iq-1))
call QsortC(A(iq:), indices(iq:))
else
call Partition(A, iq)
call QsortC(A(:iq-1))
call QsortC(A(iq:))
end if
endif
end subroutine QsortC
subroutine Partition(A, marker, indices)
real(wp), intent(in out), dimension(:) :: A
integer, intent(out) :: marker
integer, optional, intent(inout), dimension(:) :: indices
integer :: i, j, tempint
real(wp) :: temp
real(wp) :: x ! pivot point
x = A(1)
i= 0
j= size(A) + 1
do
j = j-1
do
if (A(j) <= x) exit !<- j goes to zero, giving A(0) below lower bound
j = j-1
end do
i = i+1
do
if (A(i) >= x) exit
i = i+1
end do
if (i < j) then
! exchange A(i) and A(j)
temp = A(i)
A(i) = A(j)
A(j) = temp
if (present(indices)) then
tempint = indices(i)
indices(i) = indices(j)
indices(j) = tempint
end if
elseif (i == j) then
marker = i+1
return
else
marker = i
return
endif
end do
end subroutine Partition
end module module_qsort
The specific error with the -fbounds-check flag is
At line 51 of file module_quicksort.f90 Fortran runtime error: Index '0' of dimension 1 of array 'a' below lower bound of 1
As I mentioned, the code worked in a "program" file, but not after I made a separate module file.
Cheers
It turns out I had the wrong syntax in initializing the argument variables so they were set to 0, which in turn meant I was sorting a vector of NaNs. Sorry to bother you guys with a rookie mistake :)