prank_special(3f) - [M_orderpack:RANK:PARTIAL] partially ranks an array in ASCENDING order (Insertion Sort)
Synopsis
Description
Options
Returns
Examples
Author
Maintainer
License
Subroutine Prank_Special (INVALS, IRNGT, NORD)
${TYPE} (kind=${KIND}), Intent (In) :: INVALS(:) Integer, Intent (Out) :: IRNGT(:) Integer, Intent (In) :: NORDWhere ${TYPE}(kind=${KIND}) may be
o |
Real(kind=real32)
| ||||||
Partially ranks INVALS(). Returns IRNGT(1:NORD) filled with the indices of the lowest values in the array INVALS. More technically, it does a partial ranking of the array INVALS of order NORD.NORD is restricted to the range 1 to size(IRNGT).
Internally, this version is not optimized for performance, and is thus not as difficult to read as some other ones.
The subroutine uses an insertion sort, limiting insertion to the first NORD values. It does not use any work array and is fastest when NORD is very small (2-5). but worst case behavior can easily happen (ie. if INVALS initially is inverse sorted). Therefore, In many cases, the refined Quick-sort method is faster.
INVALS array to partially sort NORD number of indices to return, restricted to 1 to size(IRNGT)
IRNGT indices of requested number (NORD) of lowest values in INVALS
Sample program:
program demo_prank_special ! partially rank N lowest values in an array use M_orderpack, only : prank_special implicit none character(len=*),parameter :: g=(*(g0,1x)) integer,allocatable :: INVALS(:) integer,allocatable :: irngt(:) integer :: nord INVALS=[10,5,7,1,4,5,6,8,9,10,1] nord=5 allocate(irngt(nord)) write(*,g)ORIGINAL:,INVALS call prank_special(INVALS,irngt,nord) write(*,g)NUMBER OF INDICES TO RETURN:,nord write(*,g)RETURNED INDICES:,irngt write(*,g)nord,SMALLEST VALUES:,INVALS(irngt(:nord)) end program demo_prank_specialResults:
ORIGINAL: 10 5 7 1 4 5 6 8 9 10 1 NUMBER OF INDICES TO RETURN: 5 RETURNED INDICES: 4 11 5 2 6 5 SMALLEST VALUES: 1 1 4 5 5
Michel Olagnon - Feb. 2000
John Urban, 2022.04.16
CC0-1.0
Nemo Release 3.1 | prank_special (3) | February 23, 2025 |