Module M_orderpack__uniinv use,intrinsic :: iso_fortran_env, only : int8, int16, int32, int64, real32, real64, real128 implicit none Private integer,parameter :: f_char=selected_char_kind("DEFAULT") public :: uniinv !> !!##NAME !! Rank_Decreasing(3f) - [M_orderpack:RANK:UNIQUE] ranks an array !! in decreasing order, with duplicate entries assigned the same !! rank(Merge-Sort) !! !!##SYNOPSIS !! !! Subroutine Rank_Decreasing (INVALS, IGOEST) !! !! ${TYPE} (kind=${KIND}), Intent (In) :: INVALS(:) !! Integer, Intent (Out) :: IGOEST(:) !! !! Where ${TYPE}(kind=${KIND}) may be !! !! o Real(kind=real32) !! o Real(kind=real64) !! o Integer(kind=int16) !! o Integer(kind=int32) !! o Integer(kind=int64) !! o Character(kind=selected_char_kind("DEFAULT"),len=*) !! !!##DESCRIPTION !! !! RANK_DECREASING(3f) generates an inverse ranking of an array, but !! with duplicate entries assigned the same rank. !! !! Internally, the routine is similar to pure merge-sort ranking, but on !! the last pass, it sets indices in IGOEST to the rank of the original !! value in an ordered set with duplicates removed. For performance !! reasons, the first two passes are taken out of the standard loop, !! and use dedicated coding. !! !!##OPTIONS !! INVALS array to rank !! IGOEST returned rank array !! !!##EXAMPLES !! !! Sample program: !! !! program demo_rank_decreasing !! ! rank input array ranking duplicates the same !! use M_orderpack, only : rank_decreasing !! implicit none !! character(len=*),parameter :: fmt='(a,*(g3.3,1x))' !! integer,allocatable,dimension(:) :: INVALS, igoest, distinct, count !! integer :: imx, i !! ! create an input array !! INVALS=[11, 11, 22, 11, 33, 33, 22, 33, 33] !! ! make an index array of the same size !! if(allocated(igoest))deallocate(igoest) !! allocate(igoest(size(INVALS))) !! print fmt, 'Original: ',INVALS !! print fmt, 'Number of indices to sort:',size(INVALS) !! ! rank input array ranking duplicates the same !! call rank_decreasing(INVALS,igoest) !! print fmt, 'Returned Indices: ',igoest(:) !! ! !! ! interrogate the results !! ! !! imx=maxval(igoest) !! print fmt, 'Number of unique indices :',imx !! ! squeeze it down to just IMX unique values !! count=[(0,i=1,imx)] ! count how many times a value occurs !! distinct=count ! array to set of unique values !! do i=1,size(INVALS) !! distinct(igoest(i))=INVALS(i) !! count(igoest(i))= count(igoest(i))+1 !! enddo !! print fmt, 'Sorted unique values: ',distinct !! print fmt, 'count of occurrences: ',count !! end program demo_rank_decreasing !! !! Results: !! !! Original: 11 11 22 11 33 33 22 33 33 !! Number of indices to sort: 9 !! Returned Indices: 1 1 2 1 3 3 2 3 3 !! Number of unique indices : 3 !! Sorted unique values: 11 22 33 !! count of occurrences: 3 2 4 !! !!##AUTHOR !! Michel Olagnon, 2000-2012 !!##MAINTAINER !! John Urban, 2022.04.16 !!##LICENSE !! CC0-1.0 interface uniinv module procedure real64_uniinv, real32_uniinv, int16_uniinv, int32_uniinv, int64_uniinv, f_char_uniinv end interface uniinv interface nearless module procedure real64_nearless, real32_nearless, int16_nearless, int32_nearless, int64_nearless, f_char_nearless end interface nearless contains Subroutine real64_uniinv (INVALS, IGOEST) Real (kind=real64), Dimension (:), Intent (In) :: INVALS Integer, Dimension (:), Intent (Out) :: IGOEST ! __________________________________________________________ Real (kind=real64) :: XTST, XDONA, XDONB ! ! __________________________________________________________ Integer, Dimension (SIZE(IGOEST)) :: JWRKT, IRNGT Integer :: LMTNA, LMTNC, IRNG, IRNG1, IRNG2, NUNI Integer :: NVAL, IIND, IWRKD, IWRK, IWRKF, JINDA, IINDA, IINDB ! NVAL = Min (SIZE(INVALS), SIZE(IGOEST)) ! Select Case (NVAL) Case (:0) Return Case (1) IGOEST (1) = 1 Return Case Default Continue End Select ! ! Fill-in the index array, creating ordered couples ! Do IIND = 2, NVAL, 2 If (INVALS(IIND-1) < INVALS(IIND)) Then IRNGT (IIND-1) = IIND - 1 IRNGT (IIND) = IIND Else IRNGT (IIND-1) = IIND IRNGT (IIND) = IIND - 1 End If End Do If (Modulo (NVAL, 2) /= 0) Then IRNGT (NVAL) = NVAL End If ! ! We will now have ordered subsets A - B - A - B - ... ! and merge A and B couples into C - C - ... ! LMTNA = 2 LMTNC = 4 ! ! First iteration. The length of the ordered subsets goes from 2 to 4 ! Do If (NVAL <= 4) Exit ! ! Loop on merges of A and B into C ! Do IWRKD = 0, NVAL - 1, 4 If ((IWRKD+4) > NVAL) Then If ((IWRKD+2) >= NVAL) Exit ! ! 1 2 3 ! If (INVALS(IRNGT(IWRKD+2)) <= INVALS(IRNGT(IWRKD+3))) Exit ! ! 1 3 2 ! If (INVALS(IRNGT(IWRKD+1)) <= INVALS(IRNGT(IWRKD+3))) Then IRNG2 = IRNGT (IWRKD+2) IRNGT (IWRKD+2) = IRNGT (IWRKD+3) IRNGT (IWRKD+3) = IRNG2 ! ! 3 1 2 ! Else IRNG1 = IRNGT (IWRKD+1) IRNGT (IWRKD+1) = IRNGT (IWRKD+3) IRNGT (IWRKD+3) = IRNGT (IWRKD+2) IRNGT (IWRKD+2) = IRNG1 End If Exit End If ! ! 1 2 3 4 ! If (INVALS(IRNGT(IWRKD+2)) <= INVALS(IRNGT(IWRKD+3))) Cycle ! ! 1 3 x x ! If (INVALS(IRNGT(IWRKD+1)) <= INVALS(IRNGT(IWRKD+3))) Then IRNG2 = IRNGT (IWRKD+2) IRNGT (IWRKD+2) = IRNGT (IWRKD+3) If (INVALS(IRNG2) <= INVALS(IRNGT(IWRKD+4))) Then ! 1 3 2 4 IRNGT (IWRKD+3) = IRNG2 Else ! 1 3 4 2 IRNGT (IWRKD+3) = IRNGT (IWRKD+4) IRNGT (IWRKD+4) = IRNG2 End If ! ! 3 x x x ! Else IRNG1 = IRNGT (IWRKD+1) IRNG2 = IRNGT (IWRKD+2) IRNGT (IWRKD+1) = IRNGT (IWRKD+3) If (INVALS(IRNG1) <= INVALS(IRNGT(IWRKD+4))) Then IRNGT (IWRKD+2) = IRNG1 If (INVALS(IRNG2) <= INVALS(IRNGT(IWRKD+4))) Then ! 3 1 2 4 IRNGT (IWRKD+3) = IRNG2 Else ! 3 1 4 2 IRNGT (IWRKD+3) = IRNGT (IWRKD+4) IRNGT (IWRKD+4) = IRNG2 End If Else ! 3 4 1 2 IRNGT (IWRKD+2) = IRNGT (IWRKD+4) IRNGT (IWRKD+3) = IRNG1 IRNGT (IWRKD+4) = IRNG2 End If End If End Do ! ! The Cs become As and Bs ! LMTNA = 4 Exit End Do ! ! Iteration loop. Each time, the length of the ordered subsets ! is doubled. ! Do If (2*LMTNA >= NVAL) Exit IWRKF = 0 LMTNC = 2 * LMTNC ! ! Loop on merges of A and B into C ! Do IWRK = IWRKF IWRKD = IWRKF + 1 JINDA = IWRKF + LMTNA IWRKF = IWRKF + LMTNC If (IWRKF >= NVAL) Then If (JINDA >= NVAL) Exit IWRKF = NVAL End If IINDA = 1 IINDB = JINDA + 1 ! ! One steps in the C subset, that we create in the final rank array ! ! Make a copy of the rank array for the iteration ! JWRKT (1:LMTNA) = IRNGT (IWRKD:JINDA) XDONA = INVALS (JWRKT(IINDA)) XDONB = INVALS (IRNGT(IINDB)) ! Do IWRK = IWRK + 1 ! ! We still have unprocessed values in both A and B ! If (XDONA > XDONB) Then IRNGT (IWRK) = IRNGT (IINDB) IINDB = IINDB + 1 If (IINDB > IWRKF) Then ! Only A still with unprocessed values IRNGT (IWRK+1:IWRKF) = JWRKT (IINDA:LMTNA) Exit End If XDONB = INVALS (IRNGT(IINDB)) Else IRNGT (IWRK) = JWRKT (IINDA) IINDA = IINDA + 1 If (IINDA > LMTNA) Exit! Only B still with unprocessed values XDONA = INVALS (JWRKT(IINDA)) End If ! End Do End Do ! ! The Cs become As and Bs ! LMTNA = 2 * LMTNA End Do ! ! Last merge of A and B into C, with removal of duplicates. ! IINDA = 1 IINDB = LMTNA + 1 NUNI = 0 ! ! One steps in the C subset, that we create in the final rank array ! JWRKT (1:LMTNA) = IRNGT (1:LMTNA) If (IINDB <= NVAL) Then If (INVALS(JWRKT(1)) < INVALS(IRNGT(IINDB)) ) Then XTST = NEARLESS (INVALS(JWRKT(1))) Else XTST = NEARLESS (INVALS(IRNGT(IINDB))) Endif Else XTST = NEARLESS (INVALS(JWRKT(1))) Endif Do IWRK = 1, NVAL ! ! We still have unprocessed values in both A and B ! If (IINDA <= LMTNA) Then If (IINDB <= NVAL) Then If (INVALS(JWRKT(IINDA)) > INVALS(IRNGT(IINDB))) Then IRNG = IRNGT (IINDB) IINDB = IINDB + 1 Else IRNG = JWRKT (IINDA) IINDA = IINDA + 1 End If Else ! ! Only A still with unprocessed values ! IRNG = JWRKT (IINDA) IINDA = IINDA + 1 End If Else ! ! Only B still with unprocessed values ! IRNG = IRNGT (IWRK) End If If (INVALS(IRNG) > XTST) Then XTST = INVALS (IRNG) NUNI = NUNI + 1 End If IGOEST (IRNG) = NUNI ! End Do ! End Subroutine real64_uniinv Function real64_nearless (XVAL) result (real64_nl) !! Nearest value less than given value ! __________________________________________________________ Real (kind=real64), Intent (In) :: XVAL Real (kind=real64) :: real64_nl ! __________________________________________________________ real64_nl = nearest (XVAL, -1.0_real64) ! End Function real64_nearless Subroutine real32_uniinv (INVALS, IGOEST) Real (kind=real32), Dimension (:), Intent (In) :: INVALS Integer, Dimension (:), Intent (Out) :: IGOEST ! __________________________________________________________ Real (kind=real32) :: XTST, XDONA, XDONB ! ! __________________________________________________________ Integer, Dimension (SIZE(IGOEST)) :: JWRKT, IRNGT Integer :: LMTNA, LMTNC, IRNG, IRNG1, IRNG2, NUNI Integer :: NVAL, IIND, IWRKD, IWRK, IWRKF, JINDA, IINDA, IINDB ! NVAL = Min (SIZE(INVALS), SIZE(IGOEST)) ! Select Case (NVAL) Case (:0) Return Case (1) IGOEST (1) = 1 Return Case Default Continue End Select ! ! Fill-in the index array, creating ordered couples ! Do IIND = 2, NVAL, 2 If (INVALS(IIND-1) < INVALS(IIND)) Then IRNGT (IIND-1) = IIND - 1 IRNGT (IIND) = IIND Else IRNGT (IIND-1) = IIND IRNGT (IIND) = IIND - 1 End If End Do If (Modulo (NVAL, 2) /= 0) Then IRNGT (NVAL) = NVAL End If ! ! We will now have ordered subsets A - B - A - B - ... ! and merge A and B couples into C - C - ... ! LMTNA = 2 LMTNC = 4 ! ! First iteration. The length of the ordered subsets goes from 2 to 4 ! Do If (NVAL <= 4) Exit ! ! Loop on merges of A and B into C ! Do IWRKD = 0, NVAL - 1, 4 If ((IWRKD+4) > NVAL) Then If ((IWRKD+2) >= NVAL) Exit ! ! 1 2 3 ! If (INVALS(IRNGT(IWRKD+2)) <= INVALS(IRNGT(IWRKD+3))) Exit ! ! 1 3 2 ! If (INVALS(IRNGT(IWRKD+1)) <= INVALS(IRNGT(IWRKD+3))) Then IRNG2 = IRNGT (IWRKD+2) IRNGT (IWRKD+2) = IRNGT (IWRKD+3) IRNGT (IWRKD+3) = IRNG2 ! ! 3 1 2 ! Else IRNG1 = IRNGT (IWRKD+1) IRNGT (IWRKD+1) = IRNGT (IWRKD+3) IRNGT (IWRKD+3) = IRNGT (IWRKD+2) IRNGT (IWRKD+2) = IRNG1 End If Exit End If ! ! 1 2 3 4 ! If (INVALS(IRNGT(IWRKD+2)) <= INVALS(IRNGT(IWRKD+3))) Cycle ! ! 1 3 x x ! If (INVALS(IRNGT(IWRKD+1)) <= INVALS(IRNGT(IWRKD+3))) Then IRNG2 = IRNGT (IWRKD+2) IRNGT (IWRKD+2) = IRNGT (IWRKD+3) If (INVALS(IRNG2) <= INVALS(IRNGT(IWRKD+4))) Then ! 1 3 2 4 IRNGT (IWRKD+3) = IRNG2 Else ! 1 3 4 2 IRNGT (IWRKD+3) = IRNGT (IWRKD+4) IRNGT (IWRKD+4) = IRNG2 End If ! ! 3 x x x ! Else IRNG1 = IRNGT (IWRKD+1) IRNG2 = IRNGT (IWRKD+2) IRNGT (IWRKD+1) = IRNGT (IWRKD+3) If (INVALS(IRNG1) <= INVALS(IRNGT(IWRKD+4))) Then IRNGT (IWRKD+2) = IRNG1 If (INVALS(IRNG2) <= INVALS(IRNGT(IWRKD+4))) Then ! 3 1 2 4 IRNGT (IWRKD+3) = IRNG2 Else ! 3 1 4 2 IRNGT (IWRKD+3) = IRNGT (IWRKD+4) IRNGT (IWRKD+4) = IRNG2 End If Else ! 3 4 1 2 IRNGT (IWRKD+2) = IRNGT (IWRKD+4) IRNGT (IWRKD+3) = IRNG1 IRNGT (IWRKD+4) = IRNG2 End If End If End Do ! ! The Cs become As and Bs ! LMTNA = 4 Exit End Do ! ! Iteration loop. Each time, the length of the ordered subsets ! is doubled. ! Do If (2*LMTNA >= NVAL) Exit IWRKF = 0 LMTNC = 2 * LMTNC ! ! Loop on merges of A and B into C ! Do IWRK = IWRKF IWRKD = IWRKF + 1 JINDA = IWRKF + LMTNA IWRKF = IWRKF + LMTNC If (IWRKF >= NVAL) Then If (JINDA >= NVAL) Exit IWRKF = NVAL End If IINDA = 1 IINDB = JINDA + 1 ! ! One steps in the C subset, that we create in the final rank array ! ! Make a copy of the rank array for the iteration ! JWRKT (1:LMTNA) = IRNGT (IWRKD:JINDA) XDONA = INVALS (JWRKT(IINDA)) XDONB = INVALS (IRNGT(IINDB)) ! Do IWRK = IWRK + 1 ! ! We still have unprocessed values in both A and B ! If (XDONA > XDONB) Then IRNGT (IWRK) = IRNGT (IINDB) IINDB = IINDB + 1 If (IINDB > IWRKF) Then ! Only A still with unprocessed values IRNGT (IWRK+1:IWRKF) = JWRKT (IINDA:LMTNA) Exit End If XDONB = INVALS (IRNGT(IINDB)) Else IRNGT (IWRK) = JWRKT (IINDA) IINDA = IINDA + 1 If (IINDA > LMTNA) Exit! Only B still with unprocessed values XDONA = INVALS (JWRKT(IINDA)) End If ! End Do End Do ! ! The Cs become As and Bs ! LMTNA = 2 * LMTNA End Do ! ! Last merge of A and B into C, with removal of duplicates. ! IINDA = 1 IINDB = LMTNA + 1 NUNI = 0 ! ! One steps in the C subset, that we create in the final rank array ! JWRKT (1:LMTNA) = IRNGT (1:LMTNA) If (IINDB <= NVAL) Then If (INVALS(JWRKT(1)) < INVALS(IRNGT(IINDB)) ) Then XTST = NEARLESS (INVALS(JWRKT(1))) Else XTST = NEARLESS (INVALS(IRNGT(IINDB))) Endif Else XTST = NEARLESS (INVALS(JWRKT(1))) Endif Do IWRK = 1, NVAL ! ! We still have unprocessed values in both A and B ! If (IINDA <= LMTNA) Then If (IINDB <= NVAL) Then If (INVALS(JWRKT(IINDA)) > INVALS(IRNGT(IINDB))) Then IRNG = IRNGT (IINDB) IINDB = IINDB + 1 Else IRNG = JWRKT (IINDA) IINDA = IINDA + 1 End If Else ! ! Only A still with unprocessed values ! IRNG = JWRKT (IINDA) IINDA = IINDA + 1 End If Else ! ! Only B still with unprocessed values ! IRNG = IRNGT (IWRK) End If If (INVALS(IRNG) > XTST) Then XTST = INVALS (IRNG) NUNI = NUNI + 1 End If IGOEST (IRNG) = NUNI ! End Do ! End Subroutine real32_uniinv Function real32_nearless (XVAL) result (real32_nl) !! Nearest value less than given value ! __________________________________________________________ Real (kind=real32), Intent (In) :: XVAL Real (kind=real32) :: real32_nl ! __________________________________________________________ real32_nl = nearest (XVAL, -1.0_real32) ! End Function real32_nearless Subroutine int16_uniinv (INVALS, IGOEST) Integer (kind=int16), Dimension (:), Intent (In) :: INVALS Integer, Dimension (:), Intent (Out) :: IGOEST ! __________________________________________________________ Integer (kind=int16) :: XTST, XDONA, XDONB ! ! __________________________________________________________ Integer, Dimension (SIZE(IGOEST)) :: JWRKT, IRNGT Integer :: LMTNA, LMTNC, IRNG, IRNG1, IRNG2, NUNI Integer :: NVAL, IIND, IWRKD, IWRK, IWRKF, JINDA, IINDA, IINDB ! NVAL = Min (SIZE(INVALS), SIZE(IGOEST)) ! Select Case (NVAL) Case (:0) Return Case (1) IGOEST (1) = 1 Return Case Default Continue End Select ! ! Fill-in the index array, creating ordered couples ! Do IIND = 2, NVAL, 2 If (INVALS(IIND-1) < INVALS(IIND)) Then IRNGT (IIND-1) = IIND - 1 IRNGT (IIND) = IIND Else IRNGT (IIND-1) = IIND IRNGT (IIND) = IIND - 1 End If End Do If (Modulo (NVAL, 2) /= 0) Then IRNGT (NVAL) = NVAL End If ! ! We will now have ordered subsets A - B - A - B - ... ! and merge A and B couples into C - C - ... ! LMTNA = 2 LMTNC = 4 ! ! First iteration. The length of the ordered subsets goes from 2 to 4 ! Do If (NVAL <= 4) Exit ! ! Loop on merges of A and B into C ! Do IWRKD = 0, NVAL - 1, 4 If ((IWRKD+4) > NVAL) Then If ((IWRKD+2) >= NVAL) Exit ! ! 1 2 3 ! If (INVALS(IRNGT(IWRKD+2)) <= INVALS(IRNGT(IWRKD+3))) Exit ! ! 1 3 2 ! If (INVALS(IRNGT(IWRKD+1)) <= INVALS(IRNGT(IWRKD+3))) Then IRNG2 = IRNGT (IWRKD+2) IRNGT (IWRKD+2) = IRNGT (IWRKD+3) IRNGT (IWRKD+3) = IRNG2 ! ! 3 1 2 ! Else IRNG1 = IRNGT (IWRKD+1) IRNGT (IWRKD+1) = IRNGT (IWRKD+3) IRNGT (IWRKD+3) = IRNGT (IWRKD+2) IRNGT (IWRKD+2) = IRNG1 End If Exit End If ! ! 1 2 3 4 ! If (INVALS(IRNGT(IWRKD+2)) <= INVALS(IRNGT(IWRKD+3))) Cycle ! ! 1 3 x x ! If (INVALS(IRNGT(IWRKD+1)) <= INVALS(IRNGT(IWRKD+3))) Then IRNG2 = IRNGT (IWRKD+2) IRNGT (IWRKD+2) = IRNGT (IWRKD+3) If (INVALS(IRNG2) <= INVALS(IRNGT(IWRKD+4))) Then ! 1 3 2 4 IRNGT (IWRKD+3) = IRNG2 Else ! 1 3 4 2 IRNGT (IWRKD+3) = IRNGT (IWRKD+4) IRNGT (IWRKD+4) = IRNG2 End If ! ! 3 x x x ! Else IRNG1 = IRNGT (IWRKD+1) IRNG2 = IRNGT (IWRKD+2) IRNGT (IWRKD+1) = IRNGT (IWRKD+3) If (INVALS(IRNG1) <= INVALS(IRNGT(IWRKD+4))) Then IRNGT (IWRKD+2) = IRNG1 If (INVALS(IRNG2) <= INVALS(IRNGT(IWRKD+4))) Then ! 3 1 2 4 IRNGT (IWRKD+3) = IRNG2 Else ! 3 1 4 2 IRNGT (IWRKD+3) = IRNGT (IWRKD+4) IRNGT (IWRKD+4) = IRNG2 End If Else ! 3 4 1 2 IRNGT (IWRKD+2) = IRNGT (IWRKD+4) IRNGT (IWRKD+3) = IRNG1 IRNGT (IWRKD+4) = IRNG2 End If End If End Do ! ! The Cs become As and Bs ! LMTNA = 4 Exit End Do ! ! Iteration loop. Each time, the length of the ordered subsets ! is doubled. ! Do If (2*LMTNA >= NVAL) Exit IWRKF = 0 LMTNC = 2 * LMTNC ! ! Loop on merges of A and B into C ! Do IWRK = IWRKF IWRKD = IWRKF + 1 JINDA = IWRKF + LMTNA IWRKF = IWRKF + LMTNC If (IWRKF >= NVAL) Then If (JINDA >= NVAL) Exit IWRKF = NVAL End If IINDA = 1 IINDB = JINDA + 1 ! ! One steps in the C subset, that we create in the final rank array ! ! Make a copy of the rank array for the iteration ! JWRKT (1:LMTNA) = IRNGT (IWRKD:JINDA) XDONA = INVALS (JWRKT(IINDA)) XDONB = INVALS (IRNGT(IINDB)) ! Do IWRK = IWRK + 1 ! ! We still have unprocessed values in both A and B ! If (XDONA > XDONB) Then IRNGT (IWRK) = IRNGT (IINDB) IINDB = IINDB + 1 If (IINDB > IWRKF) Then ! Only A still with unprocessed values IRNGT (IWRK+1:IWRKF) = JWRKT (IINDA:LMTNA) Exit End If XDONB = INVALS (IRNGT(IINDB)) Else IRNGT (IWRK) = JWRKT (IINDA) IINDA = IINDA + 1 If (IINDA > LMTNA) Exit! Only B still with unprocessed values XDONA = INVALS (JWRKT(IINDA)) End If ! End Do End Do ! ! The Cs become As and Bs ! LMTNA = 2 * LMTNA End Do ! ! Last merge of A and B into C, with removal of duplicates. ! IINDA = 1 IINDB = LMTNA + 1 NUNI = 0 ! ! One steps in the C subset, that we create in the final rank array ! JWRKT (1:LMTNA) = IRNGT (1:LMTNA) If (IINDB <= NVAL) Then If (INVALS(JWRKT(1)) < INVALS(IRNGT(IINDB)) ) Then XTST = NEARLESS (INVALS(JWRKT(1))) Else XTST = NEARLESS (INVALS(IRNGT(IINDB))) Endif Else XTST = NEARLESS (INVALS(JWRKT(1))) Endif Do IWRK = 1, NVAL ! ! We still have unprocessed values in both A and B ! If (IINDA <= LMTNA) Then If (IINDB <= NVAL) Then If (INVALS(JWRKT(IINDA)) > INVALS(IRNGT(IINDB))) Then IRNG = IRNGT (IINDB) IINDB = IINDB + 1 Else IRNG = JWRKT (IINDA) IINDA = IINDA + 1 End If Else ! ! Only A still with unprocessed values ! IRNG = JWRKT (IINDA) IINDA = IINDA + 1 End If Else ! ! Only B still with unprocessed values ! IRNG = IRNGT (IWRK) End If If (INVALS(IRNG) > XTST) Then XTST = INVALS (IRNG) NUNI = NUNI + 1 End If IGOEST (IRNG) = NUNI ! End Do ! End Subroutine int16_uniinv Function int16_nearless (XVAL) result (int16_nl) !! Nearest value less than given value ! __________________________________________________________ Integer (kind=int16), Intent (In) :: XVAL Integer (kind=int16) :: int16_nl ! __________________________________________________________ int16_nl = XVAL -1_int16 ! End Function int16_nearless Subroutine int32_uniinv (INVALS, IGOEST) Integer (kind=int32), Dimension (:), Intent (In) :: INVALS Integer, Dimension (:), Intent (Out) :: IGOEST ! __________________________________________________________ Integer (kind=int32) :: XTST, XDONA, XDONB ! ! __________________________________________________________ Integer, Dimension (SIZE(IGOEST)) :: JWRKT, IRNGT Integer :: LMTNA, LMTNC, IRNG, IRNG1, IRNG2, NUNI Integer :: NVAL, IIND, IWRKD, IWRK, IWRKF, JINDA, IINDA, IINDB ! NVAL = Min (SIZE(INVALS), SIZE(IGOEST)) ! Select Case (NVAL) Case (:0) Return Case (1) IGOEST (1) = 1 Return Case Default Continue End Select ! ! Fill-in the index array, creating ordered couples ! Do IIND = 2, NVAL, 2 If (INVALS(IIND-1) < INVALS(IIND)) Then IRNGT (IIND-1) = IIND - 1 IRNGT (IIND) = IIND Else IRNGT (IIND-1) = IIND IRNGT (IIND) = IIND - 1 End If End Do If (Modulo (NVAL, 2) /= 0) Then IRNGT (NVAL) = NVAL End If ! ! We will now have ordered subsets A - B - A - B - ... ! and merge A and B couples into C - C - ... ! LMTNA = 2 LMTNC = 4 ! ! First iteration. The length of the ordered subsets goes from 2 to 4 ! Do If (NVAL <= 4) Exit ! ! Loop on merges of A and B into C ! Do IWRKD = 0, NVAL - 1, 4 If ((IWRKD+4) > NVAL) Then If ((IWRKD+2) >= NVAL) Exit ! ! 1 2 3 ! If (INVALS(IRNGT(IWRKD+2)) <= INVALS(IRNGT(IWRKD+3))) Exit ! ! 1 3 2 ! If (INVALS(IRNGT(IWRKD+1)) <= INVALS(IRNGT(IWRKD+3))) Then IRNG2 = IRNGT (IWRKD+2) IRNGT (IWRKD+2) = IRNGT (IWRKD+3) IRNGT (IWRKD+3) = IRNG2 ! ! 3 1 2 ! Else IRNG1 = IRNGT (IWRKD+1) IRNGT (IWRKD+1) = IRNGT (IWRKD+3) IRNGT (IWRKD+3) = IRNGT (IWRKD+2) IRNGT (IWRKD+2) = IRNG1 End If Exit End If ! ! 1 2 3 4 ! If (INVALS(IRNGT(IWRKD+2)) <= INVALS(IRNGT(IWRKD+3))) Cycle ! ! 1 3 x x ! If (INVALS(IRNGT(IWRKD+1)) <= INVALS(IRNGT(IWRKD+3))) Then IRNG2 = IRNGT (IWRKD+2) IRNGT (IWRKD+2) = IRNGT (IWRKD+3) If (INVALS(IRNG2) <= INVALS(IRNGT(IWRKD+4))) Then ! 1 3 2 4 IRNGT (IWRKD+3) = IRNG2 Else ! 1 3 4 2 IRNGT (IWRKD+3) = IRNGT (IWRKD+4) IRNGT (IWRKD+4) = IRNG2 End If ! ! 3 x x x ! Else IRNG1 = IRNGT (IWRKD+1) IRNG2 = IRNGT (IWRKD+2) IRNGT (IWRKD+1) = IRNGT (IWRKD+3) If (INVALS(IRNG1) <= INVALS(IRNGT(IWRKD+4))) Then IRNGT (IWRKD+2) = IRNG1 If (INVALS(IRNG2) <= INVALS(IRNGT(IWRKD+4))) Then ! 3 1 2 4 IRNGT (IWRKD+3) = IRNG2 Else ! 3 1 4 2 IRNGT (IWRKD+3) = IRNGT (IWRKD+4) IRNGT (IWRKD+4) = IRNG2 End If Else ! 3 4 1 2 IRNGT (IWRKD+2) = IRNGT (IWRKD+4) IRNGT (IWRKD+3) = IRNG1 IRNGT (IWRKD+4) = IRNG2 End If End If End Do ! ! The Cs become As and Bs ! LMTNA = 4 Exit End Do ! ! Iteration loop. Each time, the length of the ordered subsets ! is doubled. ! Do If (2*LMTNA >= NVAL) Exit IWRKF = 0 LMTNC = 2 * LMTNC ! ! Loop on merges of A and B into C ! Do IWRK = IWRKF IWRKD = IWRKF + 1 JINDA = IWRKF + LMTNA IWRKF = IWRKF + LMTNC If (IWRKF >= NVAL) Then If (JINDA >= NVAL) Exit IWRKF = NVAL End If IINDA = 1 IINDB = JINDA + 1 ! ! One steps in the C subset, that we create in the final rank array ! ! Make a copy of the rank array for the iteration ! JWRKT (1:LMTNA) = IRNGT (IWRKD:JINDA) XDONA = INVALS (JWRKT(IINDA)) XDONB = INVALS (IRNGT(IINDB)) ! Do IWRK = IWRK + 1 ! ! We still have unprocessed values in both A and B ! If (XDONA > XDONB) Then IRNGT (IWRK) = IRNGT (IINDB) IINDB = IINDB + 1 If (IINDB > IWRKF) Then ! Only A still with unprocessed values IRNGT (IWRK+1:IWRKF) = JWRKT (IINDA:LMTNA) Exit End If XDONB = INVALS (IRNGT(IINDB)) Else IRNGT (IWRK) = JWRKT (IINDA) IINDA = IINDA + 1 If (IINDA > LMTNA) Exit! Only B still with unprocessed values XDONA = INVALS (JWRKT(IINDA)) End If ! End Do End Do ! ! The Cs become As and Bs ! LMTNA = 2 * LMTNA End Do ! ! Last merge of A and B into C, with removal of duplicates. ! IINDA = 1 IINDB = LMTNA + 1 NUNI = 0 ! ! One steps in the C subset, that we create in the final rank array ! JWRKT (1:LMTNA) = IRNGT (1:LMTNA) If (IINDB <= NVAL) Then If (INVALS(JWRKT(1)) < INVALS(IRNGT(IINDB)) ) Then XTST = NEARLESS (INVALS(JWRKT(1))) Else XTST = NEARLESS (INVALS(IRNGT(IINDB))) Endif Else XTST = NEARLESS (INVALS(JWRKT(1))) Endif Do IWRK = 1, NVAL ! ! We still have unprocessed values in both A and B ! If (IINDA <= LMTNA) Then If (IINDB <= NVAL) Then If (INVALS(JWRKT(IINDA)) > INVALS(IRNGT(IINDB))) Then IRNG = IRNGT (IINDB) IINDB = IINDB + 1 Else IRNG = JWRKT (IINDA) IINDA = IINDA + 1 End If Else ! ! Only A still with unprocessed values ! IRNG = JWRKT (IINDA) IINDA = IINDA + 1 End If Else ! ! Only B still with unprocessed values ! IRNG = IRNGT (IWRK) End If If (INVALS(IRNG) > XTST) Then XTST = INVALS (IRNG) NUNI = NUNI + 1 End If IGOEST (IRNG) = NUNI ! End Do ! End Subroutine int32_uniinv Function int32_nearless (XVAL) result (int32_nl) !! Nearest value less than given value ! __________________________________________________________ Integer (kind=int32), Intent (In) :: XVAL Integer (kind=int32) :: int32_nl ! __________________________________________________________ int32_nl = XVAL -1_int32 ! End Function int32_nearless Subroutine int64_uniinv (INVALS, IGOEST) Integer (kind=int64), Dimension (:), Intent (In) :: INVALS Integer, Dimension (:), Intent (Out) :: IGOEST ! __________________________________________________________ Integer (kind=int64) :: XTST, XDONA, XDONB ! ! __________________________________________________________ Integer, Dimension (SIZE(IGOEST)) :: JWRKT, IRNGT Integer :: LMTNA, LMTNC, IRNG, IRNG1, IRNG2, NUNI Integer :: NVAL, IIND, IWRKD, IWRK, IWRKF, JINDA, IINDA, IINDB ! NVAL = Min (SIZE(INVALS), SIZE(IGOEST)) ! Select Case (NVAL) Case (:0) Return Case (1) IGOEST (1) = 1 Return Case Default Continue End Select ! ! Fill-in the index array, creating ordered couples ! Do IIND = 2, NVAL, 2 If (INVALS(IIND-1) < INVALS(IIND)) Then IRNGT (IIND-1) = IIND - 1 IRNGT (IIND) = IIND Else IRNGT (IIND-1) = IIND IRNGT (IIND) = IIND - 1 End If End Do If (Modulo (NVAL, 2) /= 0) Then IRNGT (NVAL) = NVAL End If ! ! We will now have ordered subsets A - B - A - B - ... ! and merge A and B couples into C - C - ... ! LMTNA = 2 LMTNC = 4 ! ! First iteration. The length of the ordered subsets goes from 2 to 4 ! Do If (NVAL <= 4) Exit ! ! Loop on merges of A and B into C ! Do IWRKD = 0, NVAL - 1, 4 If ((IWRKD+4) > NVAL) Then If ((IWRKD+2) >= NVAL) Exit ! ! 1 2 3 ! If (INVALS(IRNGT(IWRKD+2)) <= INVALS(IRNGT(IWRKD+3))) Exit ! ! 1 3 2 ! If (INVALS(IRNGT(IWRKD+1)) <= INVALS(IRNGT(IWRKD+3))) Then IRNG2 = IRNGT (IWRKD+2) IRNGT (IWRKD+2) = IRNGT (IWRKD+3) IRNGT (IWRKD+3) = IRNG2 ! ! 3 1 2 ! Else IRNG1 = IRNGT (IWRKD+1) IRNGT (IWRKD+1) = IRNGT (IWRKD+3) IRNGT (IWRKD+3) = IRNGT (IWRKD+2) IRNGT (IWRKD+2) = IRNG1 End If Exit End If ! ! 1 2 3 4 ! If (INVALS(IRNGT(IWRKD+2)) <= INVALS(IRNGT(IWRKD+3))) Cycle ! ! 1 3 x x ! If (INVALS(IRNGT(IWRKD+1)) <= INVALS(IRNGT(IWRKD+3))) Then IRNG2 = IRNGT (IWRKD+2) IRNGT (IWRKD+2) = IRNGT (IWRKD+3) If (INVALS(IRNG2) <= INVALS(IRNGT(IWRKD+4))) Then ! 1 3 2 4 IRNGT (IWRKD+3) = IRNG2 Else ! 1 3 4 2 IRNGT (IWRKD+3) = IRNGT (IWRKD+4) IRNGT (IWRKD+4) = IRNG2 End If ! ! 3 x x x ! Else IRNG1 = IRNGT (IWRKD+1) IRNG2 = IRNGT (IWRKD+2) IRNGT (IWRKD+1) = IRNGT (IWRKD+3) If (INVALS(IRNG1) <= INVALS(IRNGT(IWRKD+4))) Then IRNGT (IWRKD+2) = IRNG1 If (INVALS(IRNG2) <= INVALS(IRNGT(IWRKD+4))) Then ! 3 1 2 4 IRNGT (IWRKD+3) = IRNG2 Else ! 3 1 4 2 IRNGT (IWRKD+3) = IRNGT (IWRKD+4) IRNGT (IWRKD+4) = IRNG2 End If Else ! 3 4 1 2 IRNGT (IWRKD+2) = IRNGT (IWRKD+4) IRNGT (IWRKD+3) = IRNG1 IRNGT (IWRKD+4) = IRNG2 End If End If End Do ! ! The Cs become As and Bs ! LMTNA = 4 Exit End Do ! ! Iteration loop. Each time, the length of the ordered subsets ! is doubled. ! Do If (2*LMTNA >= NVAL) Exit IWRKF = 0 LMTNC = 2 * LMTNC ! ! Loop on merges of A and B into C ! Do IWRK = IWRKF IWRKD = IWRKF + 1 JINDA = IWRKF + LMTNA IWRKF = IWRKF + LMTNC If (IWRKF >= NVAL) Then If (JINDA >= NVAL) Exit IWRKF = NVAL End If IINDA = 1 IINDB = JINDA + 1 ! ! One steps in the C subset, that we create in the final rank array ! ! Make a copy of the rank array for the iteration ! JWRKT (1:LMTNA) = IRNGT (IWRKD:JINDA) XDONA = INVALS (JWRKT(IINDA)) XDONB = INVALS (IRNGT(IINDB)) ! Do IWRK = IWRK + 1 ! ! We still have unprocessed values in both A and B ! If (XDONA > XDONB) Then IRNGT (IWRK) = IRNGT (IINDB) IINDB = IINDB + 1 If (IINDB > IWRKF) Then ! Only A still with unprocessed values IRNGT (IWRK+1:IWRKF) = JWRKT (IINDA:LMTNA) Exit End If XDONB = INVALS (IRNGT(IINDB)) Else IRNGT (IWRK) = JWRKT (IINDA) IINDA = IINDA + 1 If (IINDA > LMTNA) Exit! Only B still with unprocessed values XDONA = INVALS (JWRKT(IINDA)) End If ! End Do End Do ! ! The Cs become As and Bs ! LMTNA = 2 * LMTNA End Do ! ! Last merge of A and B into C, with removal of duplicates. ! IINDA = 1 IINDB = LMTNA + 1 NUNI = 0 ! ! One steps in the C subset, that we create in the final rank array ! JWRKT (1:LMTNA) = IRNGT (1:LMTNA) If (IINDB <= NVAL) Then If (INVALS(JWRKT(1)) < INVALS(IRNGT(IINDB)) ) Then XTST = NEARLESS (INVALS(JWRKT(1))) Else XTST = NEARLESS (INVALS(IRNGT(IINDB))) Endif Else XTST = NEARLESS (INVALS(JWRKT(1))) Endif Do IWRK = 1, NVAL ! ! We still have unprocessed values in both A and B ! If (IINDA <= LMTNA) Then If (IINDB <= NVAL) Then If (INVALS(JWRKT(IINDA)) > INVALS(IRNGT(IINDB))) Then IRNG = IRNGT (IINDB) IINDB = IINDB + 1 Else IRNG = JWRKT (IINDA) IINDA = IINDA + 1 End If Else ! ! Only A still with unprocessed values ! IRNG = JWRKT (IINDA) IINDA = IINDA + 1 End If Else ! ! Only B still with unprocessed values ! IRNG = IRNGT (IWRK) End If If (INVALS(IRNG) > XTST) Then XTST = INVALS (IRNG) NUNI = NUNI + 1 End If IGOEST (IRNG) = NUNI ! End Do ! End Subroutine int64_uniinv Function int64_nearless (XVAL) result (int64_nl) !! Nearest value less than given value ! __________________________________________________________ Integer (kind=int64), Intent (In) :: XVAL Integer (kind=int64) :: int64_nl ! __________________________________________________________ int64_nl = XVAL -1_int64 ! End Function int64_nearless Subroutine f_char_uniinv (INVALS, IGOEST) character (kind=f_char,len=*), Dimension (:), Intent (In) :: INVALS Integer, Dimension (:), Intent (Out) :: IGOEST ! __________________________________________________________ character (kind=f_char,len=len(INVALS)) :: XTST, XDONA, XDONB ! ! __________________________________________________________ Integer, Dimension (SIZE(IGOEST)) :: JWRKT, IRNGT Integer :: LMTNA, LMTNC, IRNG, IRNG1, IRNG2, NUNI Integer :: NVAL, IIND, IWRKD, IWRK, IWRKF, JINDA, IINDA, IINDB ! NVAL = Min (SIZE(INVALS), SIZE(IGOEST)) ! Select Case (NVAL) Case (:0) Return Case (1) IGOEST (1) = 1 Return Case Default Continue End Select ! ! Fill-in the index array, creating ordered couples ! Do IIND = 2, NVAL, 2 If (INVALS(IIND-1) < INVALS(IIND)) Then IRNGT (IIND-1) = IIND - 1 IRNGT (IIND) = IIND Else IRNGT (IIND-1) = IIND IRNGT (IIND) = IIND - 1 End If End Do If (Modulo (NVAL, 2) /= 0) Then IRNGT (NVAL) = NVAL End If ! ! We will now have ordered subsets A - B - A - B - ... ! and merge A and B couples into C - C - ... ! LMTNA = 2 LMTNC = 4 ! ! First iteration. The length of the ordered subsets goes from 2 to 4 ! Do If (NVAL <= 4) Exit ! ! Loop on merges of A and B into C ! Do IWRKD = 0, NVAL - 1, 4 If ((IWRKD+4) > NVAL) Then If ((IWRKD+2) >= NVAL) Exit ! ! 1 2 3 ! If (INVALS(IRNGT(IWRKD+2)) <= INVALS(IRNGT(IWRKD+3))) Exit ! ! 1 3 2 ! If (INVALS(IRNGT(IWRKD+1)) <= INVALS(IRNGT(IWRKD+3))) Then IRNG2 = IRNGT (IWRKD+2) IRNGT (IWRKD+2) = IRNGT (IWRKD+3) IRNGT (IWRKD+3) = IRNG2 ! ! 3 1 2 ! Else IRNG1 = IRNGT (IWRKD+1) IRNGT (IWRKD+1) = IRNGT (IWRKD+3) IRNGT (IWRKD+3) = IRNGT (IWRKD+2) IRNGT (IWRKD+2) = IRNG1 End If Exit End If ! ! 1 2 3 4 ! If (INVALS(IRNGT(IWRKD+2)) <= INVALS(IRNGT(IWRKD+3))) Cycle ! ! 1 3 x x ! If (INVALS(IRNGT(IWRKD+1)) <= INVALS(IRNGT(IWRKD+3))) Then IRNG2 = IRNGT (IWRKD+2) IRNGT (IWRKD+2) = IRNGT (IWRKD+3) If (INVALS(IRNG2) <= INVALS(IRNGT(IWRKD+4))) Then ! 1 3 2 4 IRNGT (IWRKD+3) = IRNG2 Else ! 1 3 4 2 IRNGT (IWRKD+3) = IRNGT (IWRKD+4) IRNGT (IWRKD+4) = IRNG2 End If ! ! 3 x x x ! Else IRNG1 = IRNGT (IWRKD+1) IRNG2 = IRNGT (IWRKD+2) IRNGT (IWRKD+1) = IRNGT (IWRKD+3) If (INVALS(IRNG1) <= INVALS(IRNGT(IWRKD+4))) Then IRNGT (IWRKD+2) = IRNG1 If (INVALS(IRNG2) <= INVALS(IRNGT(IWRKD+4))) Then ! 3 1 2 4 IRNGT (IWRKD+3) = IRNG2 Else ! 3 1 4 2 IRNGT (IWRKD+3) = IRNGT (IWRKD+4) IRNGT (IWRKD+4) = IRNG2 End If Else ! 3 4 1 2 IRNGT (IWRKD+2) = IRNGT (IWRKD+4) IRNGT (IWRKD+3) = IRNG1 IRNGT (IWRKD+4) = IRNG2 End If End If End Do ! ! The Cs become As and Bs ! LMTNA = 4 Exit End Do ! ! Iteration loop. Each time, the length of the ordered subsets ! is doubled. ! Do If (2*LMTNA >= NVAL) Exit IWRKF = 0 LMTNC = 2 * LMTNC ! ! Loop on merges of A and B into C ! Do IWRK = IWRKF IWRKD = IWRKF + 1 JINDA = IWRKF + LMTNA IWRKF = IWRKF + LMTNC If (IWRKF >= NVAL) Then If (JINDA >= NVAL) Exit IWRKF = NVAL End If IINDA = 1 IINDB = JINDA + 1 ! ! One steps in the C subset, that we create in the final rank array ! ! Make a copy of the rank array for the iteration ! JWRKT (1:LMTNA) = IRNGT (IWRKD:JINDA) XDONA = INVALS (JWRKT(IINDA)) XDONB = INVALS (IRNGT(IINDB)) ! Do IWRK = IWRK + 1 ! ! We still have unprocessed values in both A and B ! If (XDONA > XDONB) Then IRNGT (IWRK) = IRNGT (IINDB) IINDB = IINDB + 1 If (IINDB > IWRKF) Then ! Only A still with unprocessed values IRNGT (IWRK+1:IWRKF) = JWRKT (IINDA:LMTNA) Exit End If XDONB = INVALS (IRNGT(IINDB)) Else IRNGT (IWRK) = JWRKT (IINDA) IINDA = IINDA + 1 If (IINDA > LMTNA) Exit! Only B still with unprocessed values XDONA = INVALS (JWRKT(IINDA)) End If ! End Do End Do ! ! The Cs become As and Bs ! LMTNA = 2 * LMTNA End Do ! ! Last merge of A and B into C, with removal of duplicates. ! IINDA = 1 IINDB = LMTNA + 1 NUNI = 0 ! ! One steps in the C subset, that we create in the final rank array ! JWRKT (1:LMTNA) = IRNGT (1:LMTNA) If (IINDB <= NVAL) Then If (INVALS(JWRKT(1)) < INVALS(IRNGT(IINDB)) ) Then XTST = NEARLESS (INVALS(JWRKT(1))) Else XTST = NEARLESS (INVALS(IRNGT(IINDB))) Endif Else XTST = NEARLESS (INVALS(JWRKT(1))) Endif Do IWRK = 1, NVAL ! ! We still have unprocessed values in both A and B ! If (IINDA <= LMTNA) Then If (IINDB <= NVAL) Then If (INVALS(JWRKT(IINDA)) > INVALS(IRNGT(IINDB))) Then IRNG = IRNGT (IINDB) IINDB = IINDB + 1 Else IRNG = JWRKT (IINDA) IINDA = IINDA + 1 End If Else ! ! Only A still with unprocessed values ! IRNG = JWRKT (IINDA) IINDA = IINDA + 1 End If Else ! ! Only B still with unprocessed values ! IRNG = IRNGT (IWRK) End If If (INVALS(IRNG) > XTST) Then XTST = INVALS (IRNG) NUNI = NUNI + 1 End If IGOEST (IRNG) = NUNI ! End Do ! IGOEST=IGOEST+1 ! End Subroutine f_char_uniinv Function f_char_nearless (XVAL) result (f_char_nl) !! Nearest value less than given value ! __________________________________________________________ character (kind=f_char,len=*), Intent (In) :: XVAL character (kind=f_char,len=len(XVAL)) :: f_char_nl ! __________________________________________________________ f_char_nl = XVAL ! End Function f_char_nearless end module M_orderpack__uniinv