Add radix_sort

example/sorting/CMakeLists.txt

@@ -1,2 +1,3 @@
 ADD_EXAMPLE(ord_sort)
 ADD_EXAMPLE(sort)
+ADD_EXAMPLE(radix_sort)

example/sorting/example_radix_sort.f90

@@ -0,0 +1,26 @@
+program example_radix_sort
+ use iso_fortran_env, only: int8, int16, dp => real64
+ use stdlib_sorting, only: radix_sort
+ implicit none
+ integer(int8), allocatable :: arri8(:)
+ integer(int16), allocatable :: arri16(:)
+ real(dp) :: x
+ real(dp), allocatable :: arrf64(:)
+
+ arri8 = [-128, 127, 0, -1, 1]
+ call radix_sort(arri8)
+ print *, arri8
+
+ arri16 = [-32767, 32767, 0, 0, -3, 2, -3]
+ call radix_sort(arri16, reverse=.true.)
+ print *, arri16
+
+ allocate (arrf64(10))
+ x = 0.0_dp ! divide zero will arise compile error
+ arrf64 = [1.0_dp/x, 0.0_dp, 0.0_dp/x, -1.0_dp/x, -0.0_dp, 1.0_dp, -1.0_dp, 3.45_dp, -3.14_dp, 3.44_dp]
+ call radix_sort(arrf64)
+ print *, arrf64
+ ! Expect output:
+ ! nan, -inf, -3.14, -1.0, -0.0, 0.0, 1.0, 3.44, 3.45, inf
+ ! Note: the position of nan is undefined
+end program example_radix_sort

src/CMakeLists.txt

@@ -76,6 +76,7 @@ set(SRC
  stdlib_hashmap_chaining.f90
  stdlib_hashmap_open.f90
  stdlib_logger.f90
+ stdlib_sorting_radix_sort.f90
  stdlib_system.F90
  stdlib_specialfunctions.f90
  stdlib_specialfunctions_legendre.f90

src/stdlib_sorting.fypp

@@ -236,6 +236,51 @@ module stdlib_sorting
 !! ! Process the sorted data
 !! call array_search( array, values )
 !! ...
+!!```
+
+ public radix_sort
+!! Version: experimental
+!!
+!! The generic subroutine implementing the LSD radix sort algorithm to return
+!! an input array with its elements sorted in order of (non-)decreasing
+!! value. Its use has the syntax:
+!!
+!! call radix_sort( array[, work, reverse] )
+!!
+!! with the arguments:
+!!
+!! * array: the rank 1 array to be sorted. It is an `intent(inout)`
+!! argument of any of the types `integer(int8)`, `integer(int16)`,
+!! `integer(int32)`, `integer(int64)`, `real(real32)`, `real(real64)`.
+!! If both the type of `array` is real and at least one of the
+!! elements is a `NaN`, then the ordering of the result is undefined.
+!! Otherwise it is defined to be the original elements in
+!! non-decreasing order. Especially, -0.0 is lesser than 0.0.
+!!
+!! * work (optional): shall be a rank 1 array of the same type as
+!! `array`, and shall have at least `size(array)` elements. It is an
+!! `intent(inout)` argument to be used as buffer. Its value on return is
+!! undefined. If it is not present, `radix_sort` will allocate a
+!! buffer for use, and deallocate it before return. If you do several
+!! similar `radix_sort`s, reusing the `work` array is a good parctice.
+!! This argument is not present for `int8_radix_sort` because it use
+!! counting sort, so no buffer is needed.
+!!
+!! * `reverse` (optional): shall be a scalar of type default logical. It
+!! is an `intent(in)` argument. If present with a value of `.true.` then
+!! `array` will be sorted in order of non-increasing values in stable
+!! order. Otherwise index will sort `array` in order of non-decreasing
+!! values in stable order.
+!!
+!!#### Example
+!!
+!!```fortran
+!! ...
+!! ! Read random data from a file
+!! call read_file( 'dummy_file', array )
+!! ! Sort the random data
+!! call radix_sort( array )
+!! ...
 !!```
 
  public sort_index
@@ -379,6 +424,50 @@ module stdlib_sorting
 #:endfor
 
  end interface ord_sort
+ interface radix_sort
+!! Version: experimental
+!!
+!! The generic subroutine interface implementing the LSD radix sort algorithm,
+!! see https://en.wikipedia.org/wiki/Radix_sort for more details.
+!! It is always O(N) in sorting random data, but need a O(N) buffer.
+!! ([Specification](../page/specs/stdlib_sorting.html#radix_sort-sorts-an-input-array))
+!!
+
+ pure module subroutine int8_radix_sort(array, reverse)
+ integer(kind=int8), dimension(:), intent(inout) :: array
+ logical, intent(in), optional :: reverse
+ end subroutine int8_radix_sort
+
+ pure module subroutine int16_radix_sort(array, work, reverse)
+ integer(kind=int16), dimension(:), intent(inout) :: array
+ integer(kind=int16), dimension(:), intent(inout), target, optional :: work
+ logical, intent(in), optional :: reverse
+ end subroutine int16_radix_sort
+
+ pure module subroutine int32_radix_sort(array, work, reverse)
+ integer(kind=int32), dimension(:), intent(inout) :: array
+ integer(kind=int32), dimension(:), intent(inout), target, optional :: work
+ logical, intent(in), optional :: reverse
+ end subroutine int32_radix_sort
+
+ pure module subroutine int64_radix_sort(array, work, reverse)
+ integer(kind=int64), dimension(:), intent(inout) :: array
+ integer(kind=int64), dimension(:), intent(inout), target, optional :: work
+ logical, intent(in), optional :: reverse
+ end subroutine int64_radix_sort
+
+ module subroutine sp_radix_sort(array, work, reverse)
+ real(kind=sp), dimension(:), intent(inout), target :: array
+ real(kind=sp), dimension(:), intent(inout), target, optional :: work
+ logical, intent(in), optional :: reverse
+ end subroutine sp_radix_sort
+
+ module subroutine dp_radix_sort(array, work, reverse)
+ real(kind=dp), dimension(:), intent(inout), target :: array
+ real(kind=dp), dimension(:), intent(inout), target, optional :: work
+ logical, intent(in), optional :: reverse
+ end subroutine dp_radix_sort
+ end interface radix_sort
 
  interface sort
 !! Version: experimental