<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 > 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) <= 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 > l) then
+ top = top + 1
+if (top > 1000) write(0,*) 'Error: Quicksort exhausted its stack.'
+ lstack(top) = l
+ rstack(top) = s-1
+ end if
+
+ if (r > s+1) then
+ top = top + 1
+if (top > 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>