MergeSort vs. antialising rule

Viewed 116

I have the following implementation of the MergeSort algorithm in Fortran.

My question is about call merge(work(1 : half), A(half + 1:), A). Obviously I have overlapping memory, but from looking at the code in merge, this should be no problem, as long as the input arrays are sorted. (Which they are assumed to be anyway.)

On the other hand Fortran compilers may assume non aliased memory, so I always think "don't do this".

I have two questions now:

  1. When and how can I run into problems with my merge subroutine.
  2. If I cannot implement MergeSort like this, how do I do it without creating a temporary array.
!> Merge sorted arrays A and B into C while preversing order.
        subroutine merge(A, B, C)
        implicit none
        integer, intent(in) :: A(:), B(:)
        integer, intent(inout) :: C(:)

        integer :: i, j, k

        if (size(A) + size(B) > size(C)) abort

        i = 1; j = 1
        do k = 1, size(C)
          if (i <= size(A) .and. j <= size(B)) then
            if (A(i) <= B(j)) then
              C(k) = A(i)
              i = i + 1
            else
              C(k) = B(j)
              j = j + 1
            end if
          else if (i <= size(A)) then
            C(k) = A(i)
            i = i + 1
          else if (j <= size(B)) then
            C(k) = B(j)
            j = j + 1
          end if
        end do
      end subroutine merge

      recursive subroutine MergeSort(A, work)
        implicit none
        integer, intent(inout) :: A(:)
        integer, intent(inout) :: work(:)

        integer :: half
        half = (size(A) + 1) / 2
        if (size(A) < 2) then
          continue
        else if (size(A) == 2) then
          call naive_sort(A)
        else
          call MergeSort(A( : half), work)
          call MergeSort(A(half + 1 :), work)
          if (A(half) > A(half + 1)) then
            work(1 : half) = A(1 : half)
! TODO: Non aliasing rule.
            call merge(work(1 : half), A(half + 1:), A)
          endif
        end if
      end subroutine MergeSort

PS: As you perhaps notice, the array C in the merge subroutine is declared as an inout parameter, because it is later used with overlapping memory.

1 Answers

This use of aliasing in calling merge is erroneous.

With

call merge(work(1 : half), A(half + 1:), A)

the dummy argument B is associated with A(half+1:) and the dummy argument C with A which is the understood overlap.

This aliasing means that the elements of B may not be defined (which is additionally required by the intent) and that the last few elements of C may not be defined.

However, if we look at the main loop in merge we see that in general every element of C appears in a statement looking like C(k)=...: we expect at least one of those conditions inside to be true. This is therefore invalid.

To be clear: a statement like C(k)=B(j) would be an illegal definition even if the value of C(k) doesn't change as a result.

Fortunately, perhaps, there is an easy way to create a temporary array to avoid aliasing: give the dummy argument B the value attribute. You could even do the same to A and remove the workspace array.

Related