M_sort Module

sort_heap(3f) - [M_sort:sort:heapsort] indexed sort of an array
(LICENSE:PD)

  subroutine sort_heap(dat,indx)

   TYPE,intent(in) :: dat
   integer,intent(out) :: indx(size(dat))

An indexed sort of an array. The data is not moved. An integer array
is generated instead with values that are indices to the sorted
order of the data. This requires a second array the size of the input
array, which for large arrays could require a significant amount of
memory. One major advantage of this method is that any element of a
user-defined type that is a scalar intrinsic can be used to provide the
sort data and subsequently the indices can be used to access the entire
user-defined type in sorted order. This makes this seemingly simple
sort procedure usuable with the vast majority of user-defined types.

 DAT    an array of type REAL, INTEGER, or CHARACTER(KIND=kind('A')
        to be sorted

 INDX   an INTEGER array of default kind that contains the sorted
        indices.

program demo_sort_heap
use M_sort, only : sort_heap
implicit none
integer,parameter            :: isz=10000
real                         :: rr(isz)
integer                      :: ii(isz)
character(len=63)            :: cc(isz)
integer                      :: indx(isz)
integer                      :: i
write(*,*)'initializing array with ',isz,' random numbers'
CALL RANDOM_NUMBER(RR)
rr=rr*450000.0
ii=rr
do i=1,size(cc)
   cc(i)=random_string(&
   & 'abcdefghijklmnopqrstuvwxyzABCDEFGHIJKLMNOPQRSTUVWXYZ0123456789 ', &
   & len(cc))
enddo

write(*,*)'checking if real values are sorted(3f)'
call sort_heap(rr,indx)
! use the index array to actually move the input array into a sorted order
rr=rr(indx)
do i=1,isz-1
   if(rr(i).gt.rr(i+1))then
      write(*,*)'Error in sorting reals small to large ',i,rr(i),rr(i+1)
   endif
enddo
write(*,*)'test of real sort_heap(3f) complete'

write(*,*)'checking if integer values are sorted(3f)'
call sort_heap(ii,indx)
! use the index array to actually move the input array into a sorted order
ii=ii(indx)
do i=1,isz-1
   if(ii(i).gt.ii(i+1))then
      write(*,*)'Error sorting integers small to large ',i,ii(i),ii(i+1)
   endif
enddo
write(*,*)'test of integer sort_heap(3f) complete'

write(*,*)'checking if character values are sorted(3f)'
call sort_heap(cc,indx)
! use the index array to actually move the input array into a sorted order
cc=cc(indx)
do i=1,isz-1
   if(cc(i).gt.cc(i+1))then
      write(*,*)'Error sorting characters small to large ',i,cc(i),cc(i+1)
   endif
enddo
write(*,*)'test of character sort_heap(3f) complete'

contains

function random_string(chars,length) result(out)

! create random string from provided chars

character(len=*),intent(in)     :: chars
integer,intent(in)              :: length
character(len=:),allocatable    :: out
   real                         :: x
   integer                      :: ilen   ! length of list of characters
   integer                      :: which
   integer                      :: i
   ilen=len(chars)
   out=''
   if(ilen.gt.0)then
      do i=1,length
         call random_number(x)
         which=nint(real(ilen-1)*x)+1
         out=out//chars(which:which)
      enddo
   endif
end function random_string

end program demo_sort_heap

 initializing array with        10000  random numbers
 checking if real values are sorted(3f)
 test of real sort_heap(3f) complete
 checking if integer values are sorted(3f)
 test of integer sort_heap(3f) complete
 checking if character values are sorted(3f)
 test of character sort_heap(3f) complete

sort_indexed(3f) - [M_sort] indexed sort of an array
(LICENSE:PD)

  function sort_indexed(data) result(indx)

   TYPE,intent(in) :: data
   integer :: indx(size(data))

This routine is very slow on large arrays but I liked writing a sort
routine with one executable line!

An indexed sort of an array. The data is not moved. An integer array is
generated instead with values that are indices to the sorted order of
the data. This requires a second array the size of the input array,
which for large arrays could require a significant amount of memory. One
major advantage of this method is that any element of a user-defined type
that is a scalar intrinsic can be used to provide the sort data and
subsequently the indices can be used to access the entire user-defined
type in sorted order. This makes this seemingly simple sort procedure
usuable with the vast majority of user-defined types.

 DATA   an array of type REAL, INTEGER, or CHARACTER to be sorted

 INDEX  an INTEGER array of default kind that contains the sorted
        indices.

program demo_sort_indexed
use M_sort, only : sort_indexed
implicit none
integer,parameter            :: isz=10000
real                         :: rr(isz)
integer                      :: i
write(*,*)'initializing array with ',isz,' random numbers'
CALL RANDOM_NUMBER(RR)
rr=rr*450000.0
! use the index array to actually move the input array into a sorted order
rr=rr(sort_indexed(rr))
! or
!rr(sort_indexed(rr))=rr
write(*,*)'checking if values are sorted(3f)'
do i=1,isz-1
   if(rr(i).gt.rr(i+1))then
      write(*,*)'Error in sorting reals small to large ',i,rr(i),rr(i+1)
   endif
enddo
write(*,*)'test of sort_indexed(3f) complete'
end program demo_sort_indexed

sort_quick_compact(3f) - [M_sort:sort:quicksort] recursive quicksort of an array
(LICENSE: CC BY 3.0)

  function sort_quick_compact(data) result(sorted)

    type(TYPE(KIND=**),intent(in) :: data(*)
    type(TYPE(KIND=**)            :: sorted(size(data))

where TYPE may be real, doubleprecision, integer, character
or complex and of any standard kind except the character type
must be the default.

A quicksort from high to low (descending order) using a
compact recursive algorithm.

  program demo_sort_quick_compact
  use M_sort, only : sort_quick_compact
  implicit none
  integer,parameter            :: isz=10000
  real                         :: rrin(isz)
  real                         :: rrout(isz)
  integer                      :: i
  write(*,*)'initializing array with ',isz,' random numbers'
  CALL RANDOM_NUMBER(rrin)
  rrin=rrin*450000.0
  write(*,*)'sort real array with sort_quick_compact(3f)'
  rrout=sort_quick_compact(rrin)
  write(*,*)'checking '
  do i=1,isz-1
     if(rrout(i).lt.rrout(i+1))then
        write(*,*)'Error in sorting reals', &
        & i,rrout(i),rrout(i+1)
     endif
  enddo
  write(*,*)'test of sort_quick_compact(3f) complete'
  end program demo_sort_quick_compact

 initializing array with        10000  random numbers
 sort real array with sort_quick_compact(3f)
 checking index of sort_quick_compact(3f)
 test of sort_quick_compact(3f) complete

sort_quick_rx(3f) - [M_sort:sort:quicksort] indexed hybrid quicksort of an array
(LICENSE:PD)

  subroutine sort_quick_rx(data,index)

   ! one of
      real,intent(in)            :: data(:)
      doubleprecision,intent(in) :: data(:)
      integer,intent(in)         :: data(:)
      character,intent(in)       :: data(:)
      complex,intent(in)         :: data(:)

   integer,intent(out)           :: indx(size(data))

A rank hybrid quicksort. The data is not moved. An integer array is
generated instead with values that are indices to the sorted order
of the data. This requires a second array the size of the input
array, which for large arrays could require a significant amount of
memory. One major advantage of this method is that any element of a
user-defined type that is a scalar intrinsic can be used to provide the
sort data and subsequently the indices can be used to access the entire
user-defined type in sorted order. This makes this seemingly simple
sort procedure usuable with the vast majority of user-defined types.

From Leonard J. Moss of SLAC:

Here's a hybrid QuickSort I wrote a number of years ago. It's
based on suggestions in Knuth, Volume 3, and performs much better
than a pure QuickSort on short or partially ordered input arrays.

This routine performs an in-memory sort of the first N elements of
array DATA, returning into array INDEX the indices of elements of
DATA arranged in ascending order. Thus,

   DATA(INDX(1)) will be the smallest number in array DATA;
   DATA(INDX(N)) will be the largest number in DATA.

The original data is not physically rearranged. The original order
of equal input values is not necessarily preserved.

sort_quick_rx(3f) uses a hybrid QuickSort algorithm, based on several
suggestions in Knuth, Volume 3, Section 5.2.2. In particular, the
"pivot key" [my term] for dividing each subsequence is chosen to be
the median of the first, last, and middle values of the subsequence;
and the QuickSort is cut off when a subsequence has 9 or fewer
elements, and a straight insertion sort of the entire array is done
at the end. The result is comparable to a pure insertion sort for
very short arrays, and very fast for very large arrays (of order 12
micro-sec/element on the 3081K for arrays of 10K elements). It is
also not subject to the poor performance of the pure QuickSort on
partially ordered data.

Complex values are sorted by the magnitude of sqrt(r**2+i**2).

o Created: sortrx(3f): 15 Jul 1986, Len Moss
o saved from url=(0044)http://www.fortran.com/fortran/quick_sort2.f
o changed to update syntax from F77 style; John S. Urban 20161021
o generalized from only real values to include other intrinsic types;
  John S. Urban 20210110

program demo_sort_quick_rx
use M_sort, only : sort_quick_rx
implicit none
integer,parameter            :: isz=10000
real                         :: rr(isz)
integer                      :: ii(isz)
integer                      :: i
write(*,*)'initializing array with ',isz,' random numbers'
CALL RANDOM_NUMBER(RR)
rr=rr*450000.0
write(*,*)'sort real array with sort_quick_rx(3f)'
call sort_quick_rx(rr,ii)
write(*,*)'checking index of sort_quick_rx(3f)'
do i=1,isz-1
   if(rr(ii(i)).gt.rr(ii(i+1)))then
      write(*,*)'Error in sorting reals small to large ', &
      & i,rr(ii(i)),rr(ii(i+1))
   endif
enddo
write(*,*)'test of sort_quick_rx(3f) complete'
! use the index array to actually move the input array into a sorted
! order
rr=rr(ii)
do i=1,isz-1
   if(rr(i).gt.rr(i+1))then
      write(*,*)'Error in sorting reals small to large ', &
      & i,rr(i),rr(i+1)
   endif
enddo
write(*,*)'test of sort_quick_rx(3f) complete'
end program demo_sort_quick_rx

 initializing array with        10000  random numbers
 sort real array with sort_quick_rx(3f)
 checking index of sort_quick_rx(3f)
 test of sort_quick_rx(3f) complete
 test of sort_quick_rx(3f) complete

M_sort(3fm) - [M_sort::INTRO] Fortran module containing sorting
              algorithms for arrays of standard scalar types
              (LICENSE:PD)

 use M_sort, only : sort_quick_rx, sort_quick_compact
 use M_sort, only : sort_shell, sort_heap
 use M_sort, only : unique

Under development. Currently only provides a few common routines,
but it is intended that other procedures will provide a variety of
sort methods, including ...

Exchange sorts      Bubble sort, Cocktail shaker sort, Odd-even sort,
                    Comb sort, Gnome sort, Quicksort, Stooge sort,
                    Bogosort
Selection sorts     Selection sort, Heapsort, Smoothsort, Cartesian
                    tree sort, Tournament sort, Cycle sort
Insertion sorts     Insertion sort, Shellsort, Splaysort, Tree sort,
                    Library sort, Patience sorting
Merge sorts         Merge sort, Cascade merge sort, Oscillating merge
                    sort, Polyphase merge sort
Distribution sorts  American flag sort, Bead sort, Bucket sort,
                    Burstsort, Counting sort, Pigeonhole sort,
                    Proxmap sort, Radix sort, Flashsort
Concurrent sorts    Bitonic sorter, Batcher odd-even mergesort,
                    Pairwise sorting network
Hybrid sorts        Block merge sortTimsort, Introsort, Spreadsort
Other               Topological sorting, Pancake sorting, Spaghetti
                    sort

and an overview of topics concerning sorting

Theory              Computational complexity theory, Big O notation,
                    Total orderLists, InplacementStabilityComparison
                    sort, Adaptive sort, Sorting network, Integer
                    sorting, X + Y sorting, Transdichotomous model,
                    Quantum sort

In the mean time those keywords can be useful in locating materials
on the WWW, especially in Wikipedia.

Sorting generally refers to rearranging data in a specific order.
Ranking consists in finding, for each element of a set, its rank in
the sorted set, without changing the initial order of the set. In many
instances ranking is more flexible in that the order may be based
on one value of a user-defined type, and the indices can be used to
reorder virtually any type or related set; but it requires creating
an array of indices as well as the data being sorted.

Ranking is also useful when the sizes of the elements being sorted
are large, and therefore moving them around is resource-consuming.

Quicksort, also known as partition-exchange sort, uses these steps

 o Choose any element of the array to be the pivot.
 o Divide all other elements (except the pivot) into two partitions.
 o All elements less than the pivot must be in the first partition.
 o All elements greater than the pivot must be in the second partition.
 o Use recursion to sort both partitions.
 o Join the first sorted partition, the pivot, and the second sorted
   partition.

The best pivot creates partitions of equal length (or lengths differing
by 1).

The worst pivot creates an empty partition (for example, if the pivot
is the first or last element of a sorted array).

The run-time of Quicksort ranges from O(n log n) with the best pivots,
to O(n2) with the worst pivots, where n is the number of elements in
the array.

Quicksort has a reputation as the fastest sort. Optimized variants
of quicksort are common features of many languages and libraries.

sort_shell(3f) - [M_sort:sort:shellsort] Generic subroutine sorts the array X using
                 Shell's Method
(LICENSE:PD)

subroutine sort_shell(lines,order[,startcol,endcol])

 character(len=*),intent(inout) :: lines(:)
 character(len=*),intent(in)    :: order
 integer,optional,intent(in)    :: startcol, endcol

subroutine sort_shell(ints,order)

 integer,intent(inout)          :: ints(:)
 character(len=*),intent(in)    :: order

subroutine sort_shell(reals,order)

 real,intent(inout)             :: reals(:)
 character(len=*),intent(in)    :: order

subroutine sort_shell(complexs,order,type)

 character(len=*),intent(inout) :: lines(:)
 character(len=*),intent(in)    :: order
 character(len=*),intent(in)    :: type

   subroutine sort_shell(3f) sorts an array over a specified field
   in numeric or alphanumeric order.

   From Wikipedia, the free encyclopedia:

   The step-by-step process of replacing pairs of items during the shell
   sorting algorithm. Shellsort, also known as Shell sort or Shell's
   method, is an in-place comparison sort. It can be seen as either a
   generalization of sorting by exchange (bubble sort) or sorting by
   insertion (insertion sort).[3] The method starts by sorting pairs of
   elements far apart from each other, then progressively reducing the gap
   between elements to be compared. Starting with far apart elements, it
   can move some out-of-place elements into position faster than a simple
   nearest neighbor exchange. Donald Shell published the first version
   of this sort in 1959.[4][5] The running time of Shellsort is heavily
   dependent on the gap sequence it uses. For many practical variants,
   determining their time complexity remains an open problem.

   Shellsort is a generalization of insertion sort that allows the
   exchange of items that are far apart. The idea is to arrange the list
   of elements so that, starting anywhere, considering every hth element
   gives a sorted list. Such a list is said to be h-sorted. Equivalently,
   it can be thought of as h interleaved lists, each individually sorted.[6]
   Beginning with large values of h, this rearrangement allows elements
   to move long distances in the original list, reducing large amounts
   of disorder quickly, and leaving less work for smaller h-sort steps to
   do. If the file is then k-sorted for some smaller integer k, then the
   file remains h-sorted. Following this idea for a decreasing sequence of
   h values ending in 1 is guaranteed to leave a sorted list in the end.

 F90 NOTES:

  o  procedure names are declared private in this module so they are
     not accessible except by their generic name
  o  procedures must include a "use M_sort" to access the generic
     name SORT_SHELL
  o  if these routines are recompiled, routines with the use statement
     should then be recompiled and reloaded.

Usage:

 X          input/output array to sort of type CHARACTER, INTEGER,
            REAL, DOUBLEPRECISION, COMPLEX, or DOUBLEPRECISION COMPLEX.
 order      sort order
            o A for Ascending  (a-z for strings, small to large for values)
            o D for Descending (z-a for strings, large to small for
              values, default)

FOR CHARACTER DATA:

 startcol   character position in strings which starts sort field.
            Only applies to character values. Defaults to 1. Optional.
 endcol     character position in strings which ends sort field
            Only applies to character values. Defaults to end of string.
            Optional.

FOR COMPLEX AND COMPLEX(KIND=KIND(0.0D0)) DATA:

 type       Sort by

              R  for Real component,
              I  for Imaginary component,
              S  for the magnitude Sqrt(R**2+I**2)

  program demo_sort_shell
  use M_sort, only : sort_shell
  implicit none
  character(len=:),allocatable :: array(:)
  integer :: i

  array = [                                                     &
  & 'red    ','green  ','blue   ','yellow ','orange ','black  ',&
  & 'white  ','brown  ','gray   ','cyan   ','magenta',          &
  & 'purple ']

  write(*,'(a,*(a:,","))')'BEFORE ',(trim(array(i)),i=1,size(array))
  call sort_shell(array,order='a')
  write(*,'(a,*(a:,","))')'A-Z    ',(trim(array(i)),i=1,size(array))
  do i=1,size(array)-1
     if(array(i).gt.array(i+1))then
        write(*,*)'Error in sorting strings a-z'
     endif
  enddo

  array= [                                                      &
  & 'RED    ','GREEN  ','BLUE   ','YELLOW ','ORANGE ','BLACK  ',&
  & 'WHITE  ','BROWN  ','GRAY   ','CYAN   ','MAGENTA',          &
  & 'PURPLE ']

  write(*,'(a,*(a:,","))')'Before ',(trim(array(i)),i=1,size(array))
  call sort_shell(array,order='d')
  write(*,'(a,*(a:,","))')'Z-A    ',(trim(array(i)),i=1,size(array))
  do i=1,size(array)-1
     if(array(i).lt.array(i+1))then
        write(*,*)'Error in sorting strings z-a'
     endif
  enddo

  end program demo_sort_shell

   Before
   red,green,blue,yellow,orange,black,white,brown,gray,cyan,magenta,purple
   A-Z
   black,blue,brown,cyan,gray,green,magenta,orange,purple,red,white,yellow
   Before
   RED,GREEN,BLUE,YELLOW,ORANGE,BLACK,WHITE,BROWN,GRAY,CYAN,MAGENTA,PURPLE
   Z-A
   YELLOW,WHITE,RED,PURPLE,ORANGE,MAGENTA,GREEN,GRAY,CYAN,BROWN,BLUE,BLACK

1. Algorithm 201, SHELLSORT, J. Boothroyd, CACM Vol. 6, No. 8, P 445, (1963)
2. D. L. Shell, CACM, Vol. 2, P. 30, (1959)

  John S. Urban, 19970201

swap(3f) - [M_sort] elemental subroutine swaps two standard type
           variables of like type
(LICENSE:PD)

subroutine swap(X,Y)

Generic subroutine SWAP(GEN1,GEN2) swaps two variables of like type
(real, integer, complex, character, double, logical).

On output, the values of X and Y have been interchanged. Swapping is
commonly required in procedures that sort data.

SWAP(3f) is elemental, so it can operate on vectors and arrays as
well as scalar values.

program demo_swap
use M_sort, only : swap
integer             :: iarray(2)=[10,20]
real                :: rarray(2)=[11.11,22.22]
doubleprecision     :: darray(2)=[1234.56789d0,9876.54321d0]
complex             :: carray(2)=[(1234,56789),(9876,54321)]
logical             :: larray(2)=[.true.,.false.]
character(len=16)   :: string(2)=["First string    ","The other string"]

integer             :: one(13)=1
integer             :: two(13)=2

integer             :: one2(3,3)=1
integer             :: two2(3,3)=2

   print *, "integers before swap", iarray
   call swap (iarray(1), iarray(2))
   print *, "integers after swap ", iarray

   print *, "reals before swap", rarray
   call swap (rarray(1), rarray(2))
   print *, "reals after swap ", rarray

   print *, "doubles before swap", darray
   call swap (darray(1), darray(2))
   print *, "doubles after swap ", darray

   print *, "complexes before swap", carray
   call swap (carray(1), carray(2))
   print *, "complexes after swap ", carray

   print *, "logicals before swap", larray
   call swap (larray(1), larray(2))
   print *, "logicals after swap ", larray

   print *, "strings before swap", string
   call swap (string(1), string(2))
   print *, "strings after swap ", string

   write(*,*)'swap two vectors'
   write(*,'("one before: ",*(i0,:","))') one
   write(*,'("two before: ",*(i0,:","))') two
   call swap(one,two)
   write(*,'("one after: ",*(i0,:","))') one
   write(*,'("two after: ",*(i0,:","))') two

   write(*,*)'given these arrays initially each time '
   one2=1
   two2=2
   call printarrays()

   write(*,*)'swap two rows'
   one2=1
   two2=2
   call swap(one2(2,:),two2(3,:))
   call printarrays()

   write(*,*)'swap two columns'
   one2=1
   two2=2
   call swap(one2(:,2),two2(:,2))
   call printarrays()

   write(*,*)'swap two arrays with same number of elements'
   one2=1
   two2=2
   call swap(one2,two2)
   call printarrays()

   contains
   subroutine printarrays()
   integer :: i
   do i=1,size(one2(1,:))
      write(*,'(*(i0,:","))') one2(i,:)
   enddo
   write(*,*)
   do i=1,size(two2(1,:))
      write(*,'(*(i0,:","))') two2(i,:)
   enddo
   end subroutine printarrays

end program demo_swap

> integers before swap          10          20
> integers after swap           20          10
> reals before swap   11.1099997       22.2199993
> reals after swap    22.2199993       11.1099997
> doubles before swap   1234.5678900000000        9876.5432099999998
> doubles after swap    9876.5432099999998        1234.5678900000000
> complexes before swap (1234.00000,56789.0000) (9876.00000,54321.0000)
> complexes after swap  (9876.00000,54321.0000) (1234.00000,56789.0000)
> logicals before swap T F
> logicals after swap  F T
> strings before swap First string    The other string
> strings after swap  The other stringFirst string
> swap two vectors
>one before: 1,1,1,1,1,1,1,1,1,1,1,1,1
>two before: 2,2,2,2,2,2,2,2,2,2,2,2,2
>one after: 2,2,2,2,2,2,2,2,2,2,2,2,2
>two after: 1,1,1,1,1,1,1,1,1,1,1,1,1
> given these arrays initially each time
>1,1,1
>1,1,1
>1,1,1
>
>2,2,2
>2,2,2
>2,2,2
> swap two rows
>1,1,1
>2,2,2
>1,1,1
>
>2,2,2
>2,2,2
>1,1,1
> swap two columns
>1,2,1
>1,2,1
>1,2,1
>
>2,1,2
>2,1,2
>2,1,2
> swap two arrays with same number of elements
>2,2,2
>2,2,2
>2,2,2
>
>1,1,1
>1,1,1
>1,1,1