<p><b>duda</b> 2009-08-07 18:30:27 -0600 (Fri, 07 Aug 2009)</p><p>Use a non-recursive implementation of quicksort rather than <br>
a recursive mergesort. Benefits include in-place sorting <br>
and no recursion; this implementation of quicksort is also<br>
notably faster than the existing mergesort implementation.<br>
<br>
This quicksort is not a stable sort, although this shouldn't <br>
be a problem for the arrays we currently sort with it, since<br>
none of them have duplicate keys. Mergesort is still available<br>
if needed.<br>
<br>
Note: There's currently a hardwired limit to the number of<br>
iterations (1000), which should easily be sufficient for all<br>
but the linear-time worst cases. <br>
<br>
<br>
M    module_io_input.F<br>
M    module_sort.F<br>
M    module_dmpar.F<br>
</p><hr noshade><pre><font color="gray">Modified: trunk/swmodel/module_dmpar.F
===================================================================
--- trunk/swmodel/module_dmpar.F        2009-08-03 18:33:16 UTC (rev 15)
+++ trunk/swmodel/module_dmpar.F        2009-08-08 00:30:27 UTC (rev 16)
@@ -632,8 +632,8 @@
          neededListSorted(2,i) = i
       end do
 
-      call mergesort(ownedListSorted, 2, 1, nOwnedList)
-      call mergesort(neededListSorted, 2, 1, nNeededList)
+      call quicksort(nOwnedList, ownedListSorted)
+      call quicksort(nNeededList, neededListSorted)
 
       ownerList(:) = -1
 

Modified: trunk/swmodel/module_io_input.F
===================================================================
--- trunk/swmodel/module_io_input.F        2009-08-03 18:33:16 UTC (rev 15)
+++ trunk/swmodel/module_io_input.F        2009-08-08 00:30:27 UTC (rev 16)
@@ -476,19 +476,19 @@
          cellIDSorted(1,i) = domain % blocklist % mesh % indexToCellID % array(i)
          cellIDSorted(2,i) = i
       end do
-      call mergesort(cellIDSorted, 2, 1, block_graph_2Halo % nVerticesTotal)
+      call quicksort(block_graph_2Halo % nVerticesTotal, cellIDSorted)
 
       do i=1,domain % blocklist % mesh % nEdges
          edgeIDSorted(1,i) = domain % blocklist % mesh % indexToEdgeID % array(i)
          edgeIDSorted(2,i) = i
       end do
-      call mergesort(edgeIDSorted, 2, 1, nlocal_edges)
+      call quicksort(nlocal_edges, edgeIDSorted)
 
       do i=1,domain % blocklist % mesh % nVertices
          vertexIDSorted(1,i) = domain % blocklist % mesh % indexToVertexID % array(i)
          vertexIDSorted(2,i) = i
       end do
-      call mergesort(vertexIDSorted, 2, 1, nlocal_vertices)
+      call quicksort(nlocal_vertices, vertexIDSorted)
 
 
       do i=1,domain % blocklist % mesh % nCells

Modified: trunk/swmodel/module_sort.F
===================================================================
--- trunk/swmodel/module_sort.F        2009-08-03 18:33:16 UTC (rev 15)
+++ trunk/swmodel/module_sort.F        2009-08-08 00:30:27 UTC (rev 16)
@@ -67,6 +67,67 @@
    end subroutine mergesort
 
 
+   subroutine quicksort(nArray, array)
+
+      implicit none
+
+      integer, intent(in) :: nArray
+      integer, dimension(2,nArray), intent(inout) :: array
+
+      integer :: i, j, top, l, r, pivot, s
+      integer :: pivot_value
+      integer, dimension(2) :: temp
+      integer, dimension(1000) :: lstack, rstack
+
+      top = 1
+      lstack(top) = 1
+      rstack(top) = nArray
+
+      do while (top &gt; 0)
+
+         l = lstack(top)
+         r = rstack(top)
+         top = top - 1
+
+         pivot = (l+r)/2
+
+         pivot_value = array(1,pivot)
+         temp(:) = array(:,pivot)
+         array(:,pivot) = array(:,r)
+         array(:,r) = temp(:)
+
+         s = l
+         do i=l,r-1
+            if (array(1,i) &lt;= pivot_value) then
+               temp(:) = array(:,s)
+               array(:,s) = array(:,i)
+               array(:,i) = temp(:)
+               s = s + 1
+            end if
+         end do
+
+         temp(:) = array(:,s)
+         array(:,s) = array(:,r)
+         array(:,r) = temp(:)
+
+         if (s-1 &gt; l) then
+            top = top + 1
+if (top &gt; 1000) write(0,*) 'Error: Quicksort exhausted its stack.'
+            lstack(top) = l
+            rstack(top) = s-1
+         end if
+
+         if (r &gt; s+1) then
+            top = top + 1
+if (top &gt; 1000) write(0,*) 'Error: Quicksort exhausted its stack.'
+            lstack(top) = s+1
+            rstack(top) = r
+         end if
+      end do
+
+   end subroutine quicksort
+
+
    integer function binary_search(array, d1, n1, n2, key)
 
       implicit none

</font>
</pre>