<p><b>dwj07@fsu.edu</b> 2012-07-19 14:41:47 -0600 (Thu, 19 Jul 2012)</p><p><br>
        -- BRANCH COMMIT --<br>
<br>
        Sorting cell ID's based on global ID order within each halo.<br>
        This allows bit reproducibility with the same number of blocks but different number of processors.<br>
<br>
        Renaming the quicksort interface to mpas_quicksort.<br>
        Implementing 1d versions of the quicksort subroutines.<br>
</p><hr noshade><pre><font color="gray">Modified: branches/omp_blocks/multiple_blocks/src/framework/mpas_block_creator.F
===================================================================
--- branches/omp_blocks/multiple_blocks/src/framework/mpas_block_creator.F        2012-07-19 20:10:27 UTC (rev 2036)
+++ branches/omp_blocks/multiple_blocks/src/framework/mpas_block_creator.F        2012-07-19 20:41:47 UTC (rev 2037)
@@ -92,7 +92,8 @@
! Allocate array, and copy indices into array
allocate(fieldCursor % array(fieldCursor % dimSizes(1)))
- fieldCursor % array(1:fieldCursor % dimSizes(1)) = cellList(blockStart(i)+1:blockStart(i)+blockCount(i))
+ fieldCursor % array(:) = cellList(blockStart(i)+1:blockStart(i)+blockCount(i))
+ call mpas_quicksort(fieldCursor % dimSizes(1), fieldCursor % array)
! Advance cursors, and create new blocks as needed
if(i < nBlocks) then
@@ -643,8 +644,8 @@
! Setup haloIndices
haloCursor % dimSizes(1) = blockGraphWithHalo % nVerticesTotal - blockGraphWithHalo % nVertices
allocate(haloCursor % array(haloCursor % dimSizes(1)))
- haloCursor % array(:) = -1
haloCursor % array(:) = blockGraphWithHalo % vertexID(blockGraphWithHalo % nVertices+1:blockGraphWithHalo % nVerticesTotal)
+ call mpas_quicksort(haloCursor % dimSizes(1), haloCursor % array)
haloCursor % sendList => indexCursor % sendList
haloCursor % recvList => indexCursor % recvList
haloCursor % copyList => indexCursor % copyList
@@ -700,7 +701,6 @@
field1dArrayHolder => indexCursor % array
indexCursor % dimSizes(1) = nCellsSolveCursor % array(iHalo+1)
allocate(indexCursor % array(indexCursor % dimSizes(1)))
- indexCursor % array = -1
indexCursor % array(1:nCellsInBlock) = field1dArrayHolder(:)
indexCursor % array(nCellsInBlock+1:nCellsSolveCursor % array(iHalo+1)) = haloCursor % array(1:nCellsInHalo)
deallocate(field1dArrayHolder)
@@ -1145,7 +1145,7 @@
!-----------------------------------------------------------------------
subroutine mpas_block_creator_reindex_block_fields(blocklist)!{{{
- type (block_type), pointer :: blocklist !< Input/Output: Linked list of blocks
+ type (block_type), pointer :: blocklist !< Input/Output: Linked list of blocks
type (block_type), pointer :: block_ptr
@@ -1166,19 +1166,19 @@
cellIDSorted(1,i) = block_ptr % mesh % indexToCellID % array(i)
cellIDSorted(2,i) = i
end do
- call quicksort(block_ptr % mesh % nCells, cellIDSorted)
+ call mpas_quicksort(block_ptr % mesh % nCells, cellIDSorted)
do i=1,block_ptr % mesh % nEdges
edgeIDSorted(1,i) = block_ptr % mesh % indexToEdgeID % array(i)
edgeIDSorted(2,i) = i
end do
- call quicksort(block_ptr % mesh % nEdges, edgeIDSorted)
+ call mpas_quicksort(block_ptr % mesh % nEdges, edgeIDSorted)
do i=1,block_ptr % mesh % nVertices
vertexIDSorted(1,i) = block_ptr % mesh % indexToVertexID % array(i)
vertexIDSorted(2,i) = i
end do
- call quicksort(block_ptr % mesh % nVertices, vertexIDSorted)
+ call mpas_quicksort(block_ptr % mesh % nVertices, vertexIDSorted)
do i=1,block_ptr % mesh % nCells
Modified: branches/omp_blocks/multiple_blocks/src/framework/mpas_block_decomp.F
===================================================================
--- branches/omp_blocks/multiple_blocks/src/framework/mpas_block_decomp.F        2012-07-19 20:10:27 UTC (rev 2036)
+++ branches/omp_blocks/multiple_blocks/src/framework/mpas_block_decomp.F        2012-07-19 20:41:47 UTC (rev 2037)
@@ -163,7 +163,7 @@
block_count(local_block_id+1) = block_count(local_block_id+1) + 1
end do
- call quicksort(local_nvertices(dminfo % my_proc_id + 1), sorted_local_cell_list)
+ call mpas_quicksort(local_nvertices(dminfo % my_proc_id + 1), sorted_local_cell_list)
do i = 1, local_nvertices(dminfo % my_proc_id+1)
local_cell_list(i) = sorted_local_cell_list(2, i)
Modified: branches/omp_blocks/multiple_blocks/src/framework/mpas_dmpar.F
===================================================================
--- branches/omp_blocks/multiple_blocks/src/framework/mpas_dmpar.F        2012-07-19 20:10:27 UTC (rev 2036)
+++ branches/omp_blocks/multiple_blocks/src/framework/mpas_dmpar.F        2012-07-19 20:41:47 UTC (rev 2037)
@@ -688,7 +688,7 @@
kk = mpas_hash_size(neededHash)
call mpas_hash_destroy(neededHash)
- call quicksort(nUniqueNeededList, uniqueSortedNeededList)
+ call mpas_quicksort(nUniqueNeededList, uniqueSortedNeededList)
!
! Get list of index offsets for all blocks
@@ -784,14 +784,14 @@
j = j + 1
k = k + 1
end do
- call quicksort(nOwnedList, ownedListSorted)
+ call mpas_quicksort(nOwnedList, ownedListSorted)
allocate(ownedBlockSorted(2,nOwnedList))
do i=1,nOwnedList
ownedBlockSorted(1,i) = ownedList(i)
ownedBlockSorted(2,i) = ownedBlock(i)
end do
- call quicksort(nOwnedList, ownedBlockSorted)
+ call mpas_quicksort(nOwnedList, ownedBlockSorted)
allocate(neededListIndex(nNeededList))
@@ -1038,7 +1038,7 @@
ownedListSorted(2, i) = i
end do
- call quicksort(nOwnedList, ownedListSorted)
+ call mpas_quicksort(nOwnedList, ownedListSorted)
fieldCursor2 => neededListField
do while(associated(fieldCursor2))
Modified: branches/omp_blocks/multiple_blocks/src/framework/mpas_sort.F
===================================================================
--- branches/omp_blocks/multiple_blocks/src/framework/mpas_sort.F        2012-07-19 20:10:27 UTC (rev 2036)
+++ branches/omp_blocks/multiple_blocks/src/framework/mpas_sort.F        2012-07-19 20:41:47 UTC (rev 2037)
@@ -2,16 +2,16 @@
use mpas_kind_types
- interface quicksort
- module procedure mpas_quicksort_int
- module procedure mpas_quicksort_real
+ interface mpas_quicksort
+ module procedure mpas_quicksort_1dint
+ module procedure mpas_quicksort_1dreal
+ module procedure mpas_quicksort_2dint
+ module procedure mpas_quicksort_2dreal
end interface
-
contains
-
- recursive subroutine mpas_mergesort(array, d1, n1, n2)
+ recursive subroutine mpas_mergesort(array, d1, n1, n2)!{{{
implicit none
@@ -71,14 +71,137 @@
array(1:d1,n1:n2) = temp(1:d1,1:k-1)
- end subroutine mpas_mergesort
+ end subroutine mpas_mergesort!}}}
+ subroutine mpas_quicksort_1dint(nArray, array)!{{{
- subroutine mpas_quicksort_int(nArray, array)
+ implicit none
+ integer, intent(in) :: nArray
+ integer, dimension(nArray), intent(inout) :: array
+
+ integer :: i, j, top, l, r, pivot, s
+ integer :: pivot_value
+ integer, dimension(1) :: temp
+ integer, dimension(1000) :: lstack, rstack
+
+ if (nArray < 1) return
+
+ 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(pivot)
+ temp(1) = array(pivot)
+ array(pivot) = array(r)
+ array(r) = temp(1)
+
+ s = l
+ do i=l,r-1
+ if (array(i) <= pivot_value) then
+ temp(1) = array(s)
+ array(s) = array(i)
+ array(i) = temp(1)
+ s = s + 1
+ end if
+ end do
+
+ temp(1) = array(s)
+ array(s) = array(r)
+ array(r) = temp(1)
+
+ 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 mpas_quicksort_1dint!}}}
+
+ subroutine mpas_quicksort_1dreal(nArray, array)!{{{
+
implicit none
integer, intent(in) :: nArray
+ real (kind=RKIND), dimension(nArray), intent(inout) :: array
+
+ integer :: i, j, top, l, r, pivot, s
+ real (kind=RKIND) :: pivot_value
+ real (kind=RKIND), dimension(1) :: temp
+ integer, dimension(1000) :: lstack, rstack
+
+ if (nArray < 1) return
+
+ 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(pivot)
+ temp(1) = array(pivot)
+ array(pivot) = array(r)
+ array(r) = temp(1)
+
+ s = l
+ do i=l,r-1
+ if (array(i) <= pivot_value) then
+ temp(1) = array(s)
+ array(s) = array(i)
+ array(i) = temp(1)
+ s = s + 1
+ end if
+ end do
+
+ temp(1) = array(s)
+ array(s) = array(r)
+ array(r) = temp(1)
+
+ 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 mpas_quicksort_1dreal!}}}
+
+ subroutine mpas_quicksort_2dint(nArray, array)!{{{
+
+ implicit none
+
+ integer, intent(in) :: nArray
integer, dimension(2,nArray), intent(inout) :: array
integer :: i, j, top, l, r, pivot, s
@@ -134,11 +257,10 @@
end if
end do
- end subroutine mpas_quicksort_int
+ end subroutine mpas_quicksort_2dint!}}}
+ subroutine mpas_quicksort_2dreal(nArray, array)!{{{
- subroutine mpas_quicksort_real(nArray, array)
-
implicit none
integer, intent(in) :: nArray
@@ -197,11 +319,10 @@
end if
end do
- end subroutine mpas_quicksort_real
+ end subroutine mpas_quicksort_2dreal!}}}
+ integer function mpas_binary_search(array, d1, n1, n2, key)!{{{
- integer function mpas_binary_search(array, d1, n1, n2, key)
-
implicit none
integer, intent(in) :: d1, n1, n2, key
@@ -227,6 +348,6 @@
end if
end do
- end function mpas_binary_search
+ end function mpas_binary_search!}}}
end module mpas_sort
</font>
</pre>