Question

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

Was it helpful?

Solution

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 :)

Licensed under: CC-BY-SA with attribution
Not affiliated with StackOverflow
scroll top