!!! Downloaded from http://www.fortran-2000.com/rank/index.html
!!! on October 13th, 2010

!PFL 2010 Oct 13->
! I deleted the 'module' part to keep only
! the Double subroutine.

! Module m_refsor
! ! PFL 2010 Oct 13 ->
! !Integer, Parameter :: kdp = selected_real_kind(15)
! Integer, Parameter :: kdp = KIND(1.0D0)
! ! <- PFL 2010 Oct 13
! public :: refsor
! private :: kdp
! private :: R_refsor, I_refsor, D_refsor
! private :: R_inssor, I_inssor, D_inssor
! private :: R_subsor, I_subsor, D_subsor
! interface refsor
!   module procedure d_refsor, r_refsor, i_refsor
! end interface refsor
! contains

! Subroutine D_refsor (XDONT)
 Subroutine Drefsor (XDONT)
!  Sorts XDONT into ascending order - Quicksort
! __________________________________________________________
!  Quicksort chooses a "pivot" in the set, and explores the
!  array from both ends, looking for a value > pivot with the
!  increasing index, for a value <= pivot with the decreasing
!  index, and swapping them when it has found one of each.
!  The array is then subdivided in 2 ([3]) subsets:
!  { values <= pivot} {pivot} {values > pivot}
!  One then call recursively the program to sort each subset.
!  When the size of the subarray is small enough, one uses an
!  insertion sort that is faster for very small sets.
!  Michel Olagnon - Apr. 2000
! __________________________________________________________
! __________________________________________________________

   use VarTypes

      Real (KREAL), Dimension (:), Intent (InOut) :: XDONT
! __________________________________________________________
!
!
      Call D_subsor (XDONT, 1, Size (XDONT))
      Call D_inssor (XDONT)
      Return
    End Subroutine Drefsor

Recursive Subroutine D_subsor (XDONT, IDEB1, IFIN1)
!  Sorts XDONT from IDEB1 to IFIN1
! __________________________________________________________

  use VarTypes 

  IMPLICIT NONE

      Real(KREAL), dimension (:), Intent (InOut) :: XDONT
      Integer(KINT), Intent (In) :: IDEB1, IFIN1
! __________________________________________________________
      Integer(KINT), Parameter :: NINS = 16 ! Max for insertion sort
      Integer(KINT) :: ICRS, IDEB, IDCR, IFIN, IMIL
      Real(KREAL) :: XPIV, XWRK
!
      IDEB = IDEB1
      IFIN = IFIN1
!
!  If we don't have enough values to make it worth while, we leave
!  them unsorted, and the final insertion sort will take care of them
!
      If ((IFIN - IDEB) > NINS) Then
         IMIL = (IDEB+IFIN) / 2
!
!  One chooses a pivot, median of 1st, last, and middle values
!
         If (XDONT(IMIL) < XDONT(IDEB)) Then
            XWRK = XDONT (IDEB)
            XDONT (IDEB) = XDONT (IMIL)
            XDONT (IMIL) = XWRK
         End If
         If (XDONT(IMIL) > XDONT(IFIN)) Then
            XWRK = XDONT (IFIN)
            XDONT (IFIN) = XDONT (IMIL)
            XDONT (IMIL) = XWRK
            If (XDONT(IMIL) < XDONT(IDEB)) Then
               XWRK = XDONT (IDEB)
               XDONT (IDEB) = XDONT (IMIL)
               XDONT (IMIL) = XWRK
            End If
         End If
         XPIV = XDONT (IMIL)
!
!  One exchanges values to put those > pivot in the end and
!  those <= pivot at the beginning
!
         ICRS = IDEB
         IDCR = IFIN
         ECH2: Do
            Do
               ICRS = ICRS + 1
               If (ICRS >= IDCR) Then
!
!  the first  >  pivot is IDCR
!  the last   <= pivot is ICRS-1
!  Note: If one arrives here on the first iteration, then
!        the pivot is the maximum of the set, the last value is equal
!        to it, and one can reduce by one the size of the set to process,
!        as if XDONT (IFIN) > XPIV
!
                  Exit ECH2
!
               End If
               If (XDONT(ICRS) > XPIV) Exit
            End Do
            Do
               If (XDONT(IDCR) <= XPIV) Exit
               IDCR = IDCR - 1
               If (ICRS >= IDCR) Then
!
!  The last value < pivot is always ICRS-1
!
                  Exit ECH2
               End If
            End Do
!
            XWRK = XDONT (IDCR)
            XDONT (IDCR) = XDONT (ICRS)
            XDONT (ICRS) = XWRK
         End Do ECH2
!
!  One now sorts each of the two sub-intervals
!
         Call D_subsor (XDONT, IDEB1, ICRS-1)
         Call D_subsor (XDONT, IDCR, IFIN1)
      End If
      Return
    End Subroutine D_subsor

   Subroutine D_inssor (XDONT)
!  Sorts XDONT into increasing order (Insertion sort)
! __________________________________________________________

  use VarTypes 
  
  IMPLICIT NONE
      Real(KREAL), dimension (:), Intent (InOut) :: XDONT
! __________________________________________________________
      Integer(KINT) :: ICRS, IDCR
      Real(KREAL) :: XWRK
!
      Do ICRS = 2, Size (XDONT)
         XWRK = XDONT (ICRS)
         If (XWRK >= XDONT(ICRS-1)) Cycle
         XDONT (ICRS) = XDONT (ICRS-1)
         Do IDCR = ICRS - 2, 1, - 1
            If (XWRK >= XDONT(IDCR)) Exit
            XDONT (IDCR+1) = XDONT (IDCR)
         End Do
         XDONT (IDCR+1) = XWRK
      End Do
!
      Return
!
End Subroutine D_inssor


! !
! Subroutine R_refsor (XDONT)
! !  Sorts XDONT into ascending order - Quicksort
! ! __________________________________________________________
! !  Quicksort chooses a "pivot" in the set, and explores the
! !  array from both ends, looking for a value > pivot with the
! !  increasing index, for a value <= pivot with the decreasing
! !  index, and swapping them when it has found one of each.
! !  The array is then subdivided in 2 ([3]) subsets:
! !  { values <= pivot} {pivot} {values > pivot}
! !  One then call recursively the program to sort each subset.
! !  When the size of the subarray is small enough, one uses an
! !  insertion sort that is faster for very small sets.
! !  Michel Olagnon - Apr. 2000
! ! __________________________________________________________
! ! _________________________________________________________
!       Real, Dimension (:), Intent (InOut) :: XDONT
! ! __________________________________________________________
! !
! !
!       Call R_subsor (XDONT, 1, Size (XDONT))
!       Call R_inssor (XDONT)
!       Return
! End Subroutine R_refsor

! Recursive Subroutine R_subsor (XDONT, IDEB1, IFIN1)
! !  Sorts XDONT from IDEB1 to IFIN1
! ! __________________________________________________________
!       Real, dimension (:), Intent (InOut) :: XDONT
!       Integer, Intent (In) :: IDEB1, IFIN1
! ! __________________________________________________________
!       Integer, Parameter :: NINS = 16 ! Max for insertion sort
!       Integer :: ICRS, IDEB, IDCR, IFIN, IMIL
!       Real :: XPIV, XWRK
! !
!       IDEB = IDEB1
!       IFIN = IFIN1
! !
! !  If we don't have enough values to make it worth while, we leave
! !  them unsorted, and the final insertion sort will take care of them
! !
!       If ((IFIN - IDEB) > NINS) Then
!          IMIL = (IDEB+IFIN) / 2
! !
! !  One chooses a pivot, median of 1st, last, and middle values
! !
!          If (XDONT(IMIL) < XDONT(IDEB)) Then
!             XWRK = XDONT (IDEB)
!             XDONT (IDEB) = XDONT (IMIL)
!             XDONT (IMIL) = XWRK
!          End If
!          If (XDONT(IMIL) > XDONT(IFIN)) Then
!             XWRK = XDONT (IFIN)
!             XDONT (IFIN) = XDONT (IMIL)
!             XDONT (IMIL) = XWRK
!             If (XDONT(IMIL) < XDONT(IDEB)) Then
!                XWRK = XDONT (IDEB)
!                XDONT (IDEB) = XDONT (IMIL)
!                XDONT (IMIL) = XWRK
!             End If
!          End If
!          XPIV = XDONT (IMIL)
! !
! !  One exchanges values to put those > pivot in the end and
! !  those <= pivot at the beginning
! !
!          ICRS = IDEB
!          IDCR = IFIN
!          ECH2: Do
!             Do
!                ICRS = ICRS + 1
!                If (ICRS >= IDCR) Then
! !
! !  the first  >  pivot is IDCR
! !  the last   <= pivot is ICRS-1
! !  Note: If one arrives here on the first iteration, then
! !        the pivot is the maximum of the set, the last value is equal
! !        to it, and one can reduce by one the size of the set to process,
! !        as if XDONT (IFIN) > XPIV
! !
!                   Exit ECH2
! !
!                End If
!                If (XDONT(ICRS) > XPIV) Exit
!             End Do
!             Do
!                If (XDONT(IDCR) <= XPIV) Exit
!                IDCR = IDCR - 1
!                If (ICRS >= IDCR) Then
! !
! !  The last value < pivot is always ICRS-1
! !
!                   Exit ECH2
!                End If
!             End Do
! !
!             XWRK = XDONT (IDCR)
!             XDONT (IDCR) = XDONT (ICRS)
!             XDONT (ICRS) = XWRK
!          End Do ECH2
! !
! !  One now sorts each of the two sub-intervals
! !
!          Call R_subsor (XDONT, IDEB1, ICRS-1)
!          Call R_subsor (XDONT, IDCR, IFIN1)
!       End If
!       Return
!    End Subroutine R_subsor
!    Subroutine R_inssor (XDONT)
! !  Sorts XDONT into increasing order (Insertion sort)
! ! __________________________________________________________
!       Real, dimension (:), Intent (InOut) :: XDONT
! ! __________________________________________________________
!       Integer :: ICRS, IDCR
!       Real :: XWRK
! !
!       Do ICRS = 2, Size (XDONT)
!          XWRK = XDONT (ICRS)
!          If (XWRK >= XDONT(ICRS-1)) Cycle
!          XDONT (ICRS) = XDONT (ICRS-1)
!          Do IDCR = ICRS - 2, 1, - 1
!             If (XWRK >= XDONT(IDCR)) Exit
!             XDONT (IDCR+1) = XDONT (IDCR)
!          End Do
!          XDONT (IDCR+1) = XWRK
!       End Do
! !
!       Return
! !
! End Subroutine R_inssor
! !
! Subroutine I_refsor (XDONT)
! !  Sorts XDONT into ascending order - Quicksort
! ! __________________________________________________________
! !  Quicksort chooses a "pivot" in the set, and explores the
! !  array from both ends, looking for a value > pivot with the
! !  increasing index, for a value <= pivot with the decreasing
! !  index, and swapping them when it has found one of each.
! !  The array is then subdivided in 2 ([3]) subsets:
! !  { values <= pivot} {pivot} {values > pivot}
! !  One then call recursively the program to sort each subset.
! !  When the size of the subarray is small enough, one uses an
! !  insertion sort that is faster for very small sets.
! !  Michel Olagnon - Apr. 2000
! ! __________________________________________________________
! ! __________________________________________________________
!       Integer, Dimension (:), Intent (InOut)  :: XDONT
! ! __________________________________________________________
! !
! !
!       Call I_subsor (XDONT, 1, Size (XDONT))
!       Call I_inssor (XDONT)
!       Return
! End Subroutine I_refsor
! Recursive Subroutine I_subsor (XDONT, IDEB1, IFIN1)
! !  Sorts XDONT from IDEB1 to IFIN1
! ! __________________________________________________________
!       Integer, dimension (:), Intent (InOut) :: XDONT
!       Integer, Intent (In) :: IDEB1, IFIN1
! ! __________________________________________________________
!       Integer, Parameter :: NINS = 16 ! Max for insertion sort
!       Integer :: ICRS, IDEB, IDCR, IFIN, IMIL
!       Integer :: XPIV, XWRK
! !
!       IDEB = IDEB1
!       IFIN = IFIN1
! !
! !  If we don't have enough values to make it worth while, we leave
! !  them unsorted, and the final insertion sort will take care of them
! !
!       If ((IFIN - IDEB) > NINS) Then
!          IMIL = (IDEB+IFIN) / 2
! !
! !  One chooses a pivot, median of 1st, last, and middle values
! !
!          If (XDONT(IMIL) < XDONT(IDEB)) Then
!             XWRK = XDONT (IDEB)
!             XDONT (IDEB) = XDONT (IMIL)
!             XDONT (IMIL) = XWRK
!          End If
!          If (XDONT(IMIL) > XDONT(IFIN)) Then
!             XWRK = XDONT (IFIN)
!             XDONT (IFIN) = XDONT (IMIL)
!             XDONT (IMIL) = XWRK
!             If (XDONT(IMIL) < XDONT(IDEB)) Then
!                XWRK = XDONT (IDEB)
!                XDONT (IDEB) = XDONT (IMIL)
!                XDONT (IMIL) = XWRK
!             End If
!          End If
!          XPIV = XDONT (IMIL)
! !
! !  One exchanges values to put those > pivot in the end and
! !  those <= pivot at the beginning
! !
!          ICRS = IDEB
!          IDCR = IFIN
!          ECH2: Do
!             Do
!                ICRS = ICRS + 1
!                If (ICRS >= IDCR) Then
! !
! !  the first  >  pivot is IDCR
! !  the last   <= pivot is ICRS-1
! !  Note: If one arrives here on the first iteration, then
! !        the pivot is the maximum of the set, the last value is equal
! !        to it, and one can reduce by one the size of the set to process,
! !        as if XDONT (IFIN) > XPIV
! !
!                   Exit ECH2
! !
!                End If
!                If (XDONT(ICRS) > XPIV) Exit
!             End Do
!             Do
!                If (XDONT(IDCR) <= XPIV) Exit
!                IDCR = IDCR - 1
!                If (ICRS >= IDCR) Then
! !
! !  The last value < pivot is always ICRS-1
! !
!                   Exit ECH2
!                End If
!             End Do
! !
!             XWRK = XDONT (IDCR)
!             XDONT (IDCR) = XDONT (ICRS)
!             XDONT (ICRS) = XWRK
!          End Do ECH2
! !
! !  One now sorts each of the two sub-intervals
! !
!          Call I_subsor (XDONT, IDEB1, ICRS-1)
!          Call I_subsor (XDONT, IDCR, IFIN1)
!       End If
!       Return
!    End Subroutine I_subsor
!    Subroutine I_inssor (XDONT)
! !  Sorts XDONT into increasing order (Insertion sort)
! ! __________________________________________________________
!       Integer, dimension (:), Intent (InOut) :: XDONT
! ! __________________________________________________________
!       Integer :: ICRS, IDCR
!       Integer :: XWRK
! !
!       Do ICRS = 2, Size (XDONT)
!          XWRK = XDONT (ICRS)
!          If (XWRK >= XDONT(ICRS-1)) Cycle
!          XDONT (ICRS) = XDONT (ICRS-1)
!          Do IDCR = ICRS - 2, 1, - 1
!             If (XWRK >= XDONT(IDCR)) Exit
!             XDONT (IDCR+1) = XDONT (IDCR)
!          End Do
!          XDONT (IDCR+1) = XWRK
!       End Do
! !
!       Return
! !
! End Subroutine I_inssor
! !
! end module m_refsor
