scnc
/
quantum-espresso
mirror of https://gitlab.com/QEF/q-e.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

Ford-modules part 10

This commit is contained in:
fabrizio22 2021-02-24 18:22:39 +01:00
parent 739fd43748
commit 244657bc6e
6 changed files with 159 additions and 134 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

12
Modules/qeh5_module.f90
View File

@ -9,14 +9,12 @@
!------------------------------------------------------ !------------------------------------------------------
MODULE qeh5_base_module MODULE qeh5_base_module
!--------------------------------------------------- !---------------------------------------------------
!! This module contains the basic interface for basic operation for
!! serial I/O in HDF5 format. The parallel interface remains in file
!! hdf5_qe.f90 file.
! !
! this module contains the basic interface for basic operation for !! author N. Varini, P. Delugas.
! serial I/O in HDF5 format. The parallel interface remains in file !! Last revision June 2017
! hdf5_qe.f90 file.
!
! author N. Varini, P. Delugas
! last revision June 2017
!
! !
#if defined(__HDF5) #if defined(__HDF5)
USE KINDS, ONLY: DP USE KINDS, ONLY: DP

29
Modules/recips.f90
View File

@ -10,23 +10,26 @@
subroutine recips (a1, a2, a3, b1, b2, b3) subroutine recips (a1, a2, a3, b1, b2, b3)
!--------------------------------------------------------------------- !---------------------------------------------------------------------
! !! This routine generates the reciprocal lattice vectors \(b_1,b_2,b_3\)
! This routine generates the reciprocal lattice vectors b1,b2,b3 !! given the real space vectors \(a_1,a_2,a_3\). The \(b\)'s are units of
! given the real space vectors a1,a2,a3. The b's are units of 2 pi/a. !! \(2 \pi/a\).
!
! first the input variables
! !
use kinds, ONLY: DP use kinds, ONLY: DP
implicit none implicit none
real(DP) :: a1 (3), a2 (3), a3 (3), b1 (3), b2 (3), b3 (3) real(DP) :: a1 (3)
! input: first direct lattice vector !! input: first direct lattice vector
! input: second direct lattice vector real(DP) :: a2 (3)
! input: third direct lattice vector !! input: second direct lattice vector
! output: first reciprocal lattice vector real(DP) :: a3 (3)
! output: second reciprocal lattice vector !! input: third direct lattice vector
! output: third reciprocal lattice vector real(DP) :: b1 (3)
!! output: first reciprocal lattice vector
real(DP) :: b2 (3)
!! output: second reciprocal lattice vector
real(DP) :: b3 (3)
!! output: third reciprocal lattice vector
! !
! then the local variables ! ... local variables
! !
real(DP) :: den, s real(DP) :: den, s
! the denominator ! the denominator

120
Modules/recvec.f90
View File

@ -8,71 +8,72 @@
!=----------------------------------------------------------------------------=! !=----------------------------------------------------------------------------=!
MODULE gvect MODULE gvect
!=----------------------------------------------------------------------------=! !=----------------------------------------------------------------------------=!
!! Variables describing the reciprocal lattice vectors
! ... variables describing the reciprocal lattice vectors !! G vectors with |G|^2 < ecutrho, cut-off for charge density
! ... G vectors with |G|^2 < ecutrho, cut-off for charge density !! With gamma tricks, G-vectors are divided into two half-spheres,
! ... With gamma tricks, G-vectors are divided into two half-spheres, !! G> and G<, containing G and -G (G=0 is in G>)
! ... G> and G<, containing G and -G (G=0 is in G>) !! This is referred to as the "dense" (or "hard", or "thick") grid
! ... This is referred to as the "dense" (or "hard", or "thick") grid
USE kinds, ONLY: DP USE kinds, ONLY: DP
IMPLICIT NONE IMPLICIT NONE
SAVE SAVE
INTEGER :: ngm = 0 ! local number of G vectors (on this processor) INTEGER :: ngm = 0
! with gamma tricks, only vectors in G> !! local number of G vectors (on this processor) with gamma tricks, only
INTEGER :: ngm_g= 0 ! global number of G vectors (summed on all procs) !! vectors in G>
! in serial execution, ngm_g = ngm INTEGER :: ngm_g = 0
INTEGER :: ngl = 0 ! number of G-vector shells !! global number of G vectors (summed on all procs) in serial execution,
INTEGER :: ngmx = 0 ! local number of G vectors, maximum across all procs !! \(\text{ngm}_g = \text{ngm}\).
INTEGER :: ngl = 0
!! number of G-vector shells
INTEGER :: ngmx = 0
!! local number of G vectors, maximum across all procs
REAL(DP) :: ecutrho = 0.0_DP ! energy cut-off for charge density REAL(DP) :: ecutrho = 0.0_DP
REAL(DP) :: gcutm = 0.0_DP ! ecutrho/(2 pi/a)^2, cut-off for |G|^2 !! energy cut-off for charge density
REAL(DP) :: gcutm = 0.0_DP
!! \(\text{ecutrho}/(2 \pi/a)^2\), cut-off for \(\|G\|^2\)
INTEGER :: gstart = 2 ! index of the first G vector whose module is > 0 INTEGER :: gstart = 2
! Needed in parallel execution: gstart=2 for the !! index of the first G vector whose module is > 0. Needed in parallel
! proc that holds G=0, gstart=1 for all others !execution: gstart=2 for the proc that holds G=0, gstart=1 for all others.
! G^2 in increasing order (in units of tpiba2=(2pi/a)^2)
!
REAL(DP), ALLOCATABLE, TARGET :: gg(:) REAL(DP), ALLOCATABLE, TARGET :: gg(:)
REAL(DP), ALLOCATABLE :: gg_d(:) ! device !! \(G^2\) in increasing order (in units of \(\text{tpiba2}=(2\pi/a)^2\) )
REAL(DP), ALLOCATABLE :: gg_d(:)
!! device
! gl(i) = i-th shell of G^2 (in units of tpiba2)
! igtongl(n) = shell index for n-th G-vector
!
REAL(DP), POINTER, PROTECTED :: gl(:) REAL(DP), POINTER, PROTECTED :: gl(:)
!! gl(i) = i-th shell of G^2 (in units of tpiba2)
INTEGER, ALLOCATABLE, TARGET, PROTECTED :: igtongl(:) INTEGER, ALLOCATABLE, TARGET, PROTECTED :: igtongl(:)
!! igtongl(n) = shell index for n-th G-vector
! Duplicate of the above variables (new style duplication). ! Duplicate of the above variables (new style duplication).
REAL(DP), ALLOCATABLE :: gl_d(:) ! device REAL(DP), ALLOCATABLE :: gl_d(:) ! device
INTEGER, ALLOCATABLE, TARGET :: igtongl_d(:) ! device INTEGER, ALLOCATABLE, TARGET :: igtongl_d(:) ! device
! !
! G-vectors cartesian components ( in units tpiba =(2pi/a) )
!
REAL(DP), ALLOCATABLE, TARGET :: g(:,:) REAL(DP), ALLOCATABLE, TARGET :: g(:,:)
!! G-vectors cartesian components ( in units \(\text{tpiba} =(2\pi/a)\) )
REAL(DP), ALLOCATABLE :: g_d(:,:) ! device REAL(DP), ALLOCATABLE :: g_d(:,:) ! device
! mill = miller index of G vectors (local to each processor)
! G(:) = mill(1)*bg(:,1)+mill(2)*bg(:,2)+mill(3)*bg(:,3)
! where bg are the reciprocal lattice basis vectors
! !
INTEGER, ALLOCATABLE, TARGET :: mill(:,:) INTEGER, ALLOCATABLE, TARGET :: mill(:,:)
!! miller index of G vectors (local to each processor)
!! G(:) = mill(1)*bg(:,1)+mill(2)*bg(:,2)+mill(3)*bg(:,3)
!! where bg are the reciprocal lattice basis vectors.
INTEGER, ALLOCATABLE :: mill_d(:,:) ! device INTEGER, ALLOCATABLE :: mill_d(:,:) ! device
! ig_l2g = converts a local G-vector index into the global index
! ("l2g" means local to global): ig_l2g(i) = index of i-th
! local G-vector in the global array of G-vectors
! !
INTEGER, ALLOCATABLE, TARGET :: ig_l2g(:) INTEGER, ALLOCATABLE, TARGET :: ig_l2g(:)
! !! converts a local G-vector index into the global index
! mill_g = miller index of all G vectors !! ("l2g" means local to global): ig\_l2g(i) = index of i-th
!! local G-vector in the global array of G-vectors
! !
INTEGER, ALLOCATABLE, TARGET :: mill_g(:,:) INTEGER, ALLOCATABLE, TARGET :: mill_g(:,:)
!! Miller index of all G vectors
! !
! the phases e^{-iG*tau_s} used to calculate structure factors COMPLEX(DP), ALLOCATABLE :: eigts1(:,:)
! !! the phases \(e^{-iG\text{tau}_s}\) used to calculate structure factors.
COMPLEX(DP), ALLOCATABLE :: eigts1(:,:), eigts2(:,:), eigts3(:,:) COMPLEX(DP), ALLOCATABLE :: eigts2(:,:), eigts3(:,:)
COMPLEX(DP), ALLOCATABLE :: eigts1_d(:,:) ! device COMPLEX(DP), ALLOCATABLE :: eigts1_d(:,:) ! device
COMPLEX(DP), ALLOCATABLE :: eigts2_d(:,:) COMPLEX(DP), ALLOCATABLE :: eigts2_d(:,:)
COMPLEX(DP), ALLOCATABLE :: eigts3_d(:,:) COMPLEX(DP), ALLOCATABLE :: eigts3_d(:,:)
@ -85,13 +86,14 @@
SUBROUTINE gvect_init( ngm_ , comm ) SUBROUTINE gvect_init( ngm_ , comm )
! !
! Set local and global dimensions, allocate arrays !! Set local and global dimensions, allocate arrays
! !
USE mp, ONLY: mp_max, mp_sum USE mp, ONLY: mp_max, mp_sum
USE control_flags, ONLY : use_gpu USE control_flags, ONLY : use_gpu
IMPLICIT NONE IMPLICIT NONE
INTEGER, INTENT(IN) :: ngm_ INTEGER, INTENT(IN) :: ngm_
INTEGER, INTENT(IN) :: comm ! communicator of the group on which g-vecs are distributed INTEGER, INTENT(IN) :: comm
!! communicator of the group on which g-vecs are distributed
! !
ngm = ngm_ ngm = ngm_
! !
@ -123,6 +125,8 @@
END SUBROUTINE gvect_init END SUBROUTINE gvect_init
SUBROUTINE deallocate_gvect(vc) SUBROUTINE deallocate_gvect(vc)
!! Deallocate G-vector related arrays.
!
USE control_flags, ONLY : use_gpu USE control_flags, ONLY : use_gpu
IMPLICIT NONE IMPLICIT NONE
! !
@ -170,9 +174,8 @@
!----------------------------------------------------------------------- !-----------------------------------------------------------------------
SUBROUTINE gshells ( vc ) SUBROUTINE gshells ( vc )
!---------------------------------------------------------------------- !----------------------------------------------------------------------
! !! Calculate number of G shells: ngl, and the index ng = igtongl(ig)
! calculate number of G shells: ngl, and the index ng = igtongl(ig) !! that gives the shell index ng for (local) G-vector of index ig.
! that gives the shell index ng for (local) G-vector of index ig
! !
USE kinds, ONLY : DP USE kinds, ONLY : DP
USE constants, ONLY : eps8 USE constants, ONLY : eps8
@ -229,6 +232,9 @@
!=----------------------------------------------------------------------------=! !=----------------------------------------------------------------------------=!
MODULE gvecs MODULE gvecs
!=----------------------------------------------------------------------------=! !=----------------------------------------------------------------------------=!
!! G vectors with \(\|G\|^2 < 4*\text{ecutwfc}\), cut-off for wavefunctions
!! ("smooth" grid). Gamma tricks and units as for the "dense" grid.
USE kinds, ONLY: DP USE kinds, ONLY: DP
IMPLICIT NONE IMPLICIT NONE
@ -237,25 +243,33 @@
! ... G vectors with |G|^2 < 4*ecutwfc, cut-off for wavefunctions ! ... G vectors with |G|^2 < 4*ecutwfc, cut-off for wavefunctions
! ... ("smooth" grid). Gamma tricks and units as for the "dense" grid ! ... ("smooth" grid). Gamma tricks and units as for the "dense" grid
! !
INTEGER :: ngms = 0 ! local number of smooth vectors (on this processor) INTEGER :: ngms = 0
INTEGER :: ngms_g=0 ! global number of smooth vectors (summed on procs) !! local number of smooth vectors (on this processor)
! in serial execution this is equal to ngms INTEGER :: ngms_g=0
INTEGER :: ngsx = 0 ! local number of smooth vectors, max across procs !! global number of smooth vectors (summed on procs)
!! in serial execution this is equal to ngms
INTEGER :: ngsx = 0
!! local number of smooth vectors, max across procs
REAL(DP) :: ecuts = 0.0_DP ! energy cut-off = 4*ecutwfc REAL(DP) :: ecuts = 0.0_DP
REAL(DP) :: gcutms= 0.0_DP ! ecuts/(2 pi/a)^2, cut-off for |G|^2 !! energy cut-off = 4*ecutwfc
REAL(DP) :: gcutms= 0.0_DP
!! ecuts/(2 pi/a)^2, cut-off for |G|^2
REAL(DP) :: dual = 0.0_DP ! ecutrho=dual*ecutwfc REAL(DP) :: dual = 0.0_DP
LOGICAL :: doublegrid = .FALSE. ! true if smooth and dense grid differ !! ecutrho=dual*ecutwfc
! doublegrid = (dual > 4) LOGICAL :: doublegrid = .FALSE.
!! true if smooth and dense grid differ. doublegrid = (dual > 4)
CONTAINS CONTAINS
SUBROUTINE gvecs_init( ngs_ , comm ) SUBROUTINE gvecs_init( ngs_ , comm )
!! G-vector initialization
USE mp, ONLY: mp_max, mp_sum USE mp, ONLY: mp_max, mp_sum
IMPLICIT NONE IMPLICIT NONE
INTEGER, INTENT(IN) :: ngs_ INTEGER, INTENT(IN) :: ngs_
INTEGER, INTENT(IN) :: comm ! communicator of the group on which g-vecs are distributed INTEGER, INTENT(IN) :: comm
!! communicator of the group on which g-vecs are distributed
! !
ngms = ngs_ ngms = ngs_
! !

33
Modules/recvec_subs.f90
View File

@ -9,9 +9,8 @@
!=----------------------------------------------------------------------= !=----------------------------------------------------------------------=
MODULE recvec_subs MODULE recvec_subs
!=----------------------------------------------------------------------= !=----------------------------------------------------------------------=
!! Subroutines generating G-vectors and variables nl* needed to map
! ... subroutines generating G-vectors and variables nl* needed to map !! G-vector components onto the FFT grid(s) in reciprocal space.
! ... G-vector components onto the FFT grid(s) in reciprocal space
! !
USE kinds, ONLY : dp USE kinds, ONLY : dp
USE fft_types, ONLY: fft_stick_index, fft_type_descriptor USE fft_types, ONLY: fft_stick_index, fft_type_descriptor
@ -30,11 +29,10 @@ CONTAINS
SUBROUTINE ggen ( dfftp, gamma_only, at, bg, gcutm, ngm_g, ngm, & SUBROUTINE ggen ( dfftp, gamma_only, at, bg, gcutm, ngm_g, ngm, &
g, gg, mill, ig_l2g, gstart, no_global_sort ) g, gg, mill, ig_l2g, gstart, no_global_sort )
!---------------------------------------------------------------------- !----------------------------------------------------------------------
! !! This routine generates all the reciprocal lattice vectors
! This routine generates all the reciprocal lattice vectors !! contained in the sphere of radius gcutm. Furthermore it
! contained in the sphere of radius gcutm. Furthermore it !! computes the indices nl which give the correspondence
! computes the indices nl which give the correspondence !! between the fft mesh points and the array of g vectors.
! between the fft mesh points and the array of g vectors.
! !
USE mp, ONLY: mp_rank, mp_size, mp_sum USE mp, ONLY: mp_rank, mp_size, mp_sum
USE constants, ONLY : eps8 USE constants, ONLY : eps8
@ -271,27 +269,24 @@ CONTAINS
SUBROUTINE ggens( dffts, gamma_only, at, g, gg, mill, gcutms, ngms, & SUBROUTINE ggens( dffts, gamma_only, at, g, gg, mill, gcutms, ngms, &
gs, ggs ) gs, ggs )
!----------------------------------------------------------------------- !-----------------------------------------------------------------------
! !! Initialize number and indices of g-vectors for a subgrid,
! Initialize number and indices of g-vectors for a subgrid, !! for exactly the same ordering as for the original FFT grid
! for exactly the same ordering as for the original FFT grid
!
!--------------------------------------------------------------------
! !
IMPLICIT NONE IMPLICIT NONE
! !
LOGICAL, INTENT(IN) :: gamma_only LOGICAL, INTENT(IN) :: gamma_only
TYPE (fft_type_descriptor), INTENT(INOUT) :: dffts TYPE (fft_type_descriptor), INTENT(INOUT) :: dffts
! primitive lattice vectors !! primitive lattice vectors
REAL(dp), INTENT(IN) :: at(3,3) REAL(dp), INTENT(IN) :: at(3,3)
! G-vectors in FFT grid !! G-vectors in FFT grid
REAL(dp), INTENT(IN) :: g(:,:), gg(:) REAL(dp), INTENT(IN) :: g(:,:), gg(:)
! Miller indices for G-vectors of FFT grid !! Miller indices for G-vectors of FFT grid
INTEGER, INTENT(IN) :: mill(:,:) INTEGER, INTENT(IN) :: mill(:,:)
! cutoff for subgrid !! cutoff for subgrid
REAL(DP), INTENT(IN):: gcutms REAL(DP), INTENT(IN):: gcutms
! Local number of G-vectors in subgrid !! Local number of G-vectors in subgrid
INTEGER, INTENT(OUT):: ngms INTEGER, INTENT(OUT):: ngms
! Optionally: G-vectors and modules !! Optionally: G-vectors and modules
REAL(DP), INTENT(INOUT), POINTER, OPTIONAL:: gs(:,:), ggs(:) REAL(DP), INTENT(INOUT), POINTER, OPTIONAL:: gs(:,:), ggs(:)
! !
INTEGER :: i, ng, ngm INTEGER :: i, ng, ngm

24
Modules/remove_tot_torque.f90
View File

@ -8,12 +8,20 @@
!---------------------------------------------------------------------------- !----------------------------------------------------------------------------
SUBROUTINE remove_tot_torque( nat, tau, mass, force ) SUBROUTINE remove_tot_torque( nat, tau, mass, force )
!---------------------------------------------------------------------------- !----------------------------------------------------------------------------
!! This routine sets to zero the total torque associated to the internal
!! forces acting on the atoms by correcting the force vector.
! !
! ... This routine sets to zero the total torque associated to the internal !! The algorithm is based on the following expressions ( F' is the
! ... forces acting on the atoms by correcting the force vector. !! torqueless force ) :
!
!! $$ {\bf m} = \frac{1}{N} \sum_i d {\bf R}_i\times {\bf F}_i $$
!
!! $$ {\bf F}_i' = {\bf F}_i - \frac{1}{\|d {\bf R}_i\|^2} {\bf m} \times d {\bf R}_i $$
!
!! with \( d {\bf R}_i = {\bf R}_i - {\bf R}_\text{cm} \)
!
!! Written by Carlo Sbraccia (2006).
! !
! ... The algorithm is based on the following expressions ( F' is the
! ... torqueless force ) :
! _ ! _
! _ 1 \ __ _ __ _ _ ! _ 1 \ __ _ __ _ _
! ... m = --- /_ dR_i /\ F_i , dR_i = ( R_i - R_cm ) , ! ... m = --- /_ dR_i /\ F_i , dR_i = ( R_i - R_cm ) ,
@ -24,16 +32,20 @@ SUBROUTINE remove_tot_torque( nat, tau, mass, force )
! |dR_i|^2 ! |dR_i|^2
! !
! !
! ... written by carlo sbraccia (2006)
!
USE kinds, ONLY : DP USE kinds, ONLY : DP
! !
IMPLICIT NONE IMPLICIT NONE
! !
INTEGER, INTENT(IN) :: nat INTEGER, INTENT(IN) :: nat
!! number of atoms
REAL(DP), INTENT(IN) :: tau(3,nat) REAL(DP), INTENT(IN) :: tau(3,nat)
!! atomic positions
REAL(DP), INTENT(IN) :: mass(nat) REAL(DP), INTENT(IN) :: mass(nat)
!! atomic mass
REAL(DP), INTENT(INOUT) :: force(3,nat) REAL(DP), INTENT(INOUT) :: force(3,nat)
!! force on atoms
!
! ... local variables
! !
INTEGER :: ia INTEGER :: ia
REAL(DP) :: m(3), mo(3), tauref(3), delta(3), sumf(3) REAL(DP) :: m(3), mo(3), tauref(3), delta(3), sumf(3)

67
Modules/sort.f90
View File

@ -8,18 +8,19 @@
!--------------------------------------------------------------------- !---------------------------------------------------------------------
subroutine hpsort_eps (n, ra, ind, eps) subroutine hpsort_eps (n, ra, ind, eps)
!--------------------------------------------------------------------- !---------------------------------------------------------------------
! sort an array ra(1:n) into ascending order using heapsort algorithm, !! Sort an array ra(1:n) into ascending order using heapsort algorithm,
! and considering two elements being equal if their values differ !! and considering two elements being equal if their values differ
! for less than "eps". !! for less than "eps".
! n is input, ra is replaced on output by its sorted rearrangement. !! \(\text{n}\) is input, \(\text{ra}\) is replaced on output by its
! create an index table (ind) by making an exchange in the index array !! sorted rearrangement.
! whenever an exchange is made on the sorted data array (ra). !! Create an index table (ind) by making an exchange in the index array
! in case of equal values in the data array (ra) the values in the !! whenever an exchange is made on the sorted data array (\(\text{ra}\)).
! index array (ind) are used to order the entries. !! In case of equal values in the data array (\(\text{ra}\)) the values
! if on input ind(1) = 0 then indices are initialized in the routine, !! in the index array (ind) are used to order the entries.
! if on input ind(1) != 0 then indices are assumed to have been !! If on input ind(1) = 0 then indices are initialized in the routine,
! initialized before entering the routine and these !! if on input ind(1) != 0 then indices are assumed to have been
! indices are carried around during the sorting process !! initialized before entering the routine and these indices are carried
!! around during the sorting process.
! !
! no work space needed ! ! no work space needed !
! free us from machine-dependent sorting-routines ! ! free us from machine-dependent sorting-routines !
@ -134,16 +135,17 @@ end subroutine hpsort_eps
!--------------------------------------------------------------------- !---------------------------------------------------------------------
subroutine hpsort (n, ra, ind) subroutine hpsort (n, ra, ind)
!--------------------------------------------------------------------- !---------------------------------------------------------------------
! sort an array ra(1:n) into ascending order using heapsort algorithm. !! Sort an array ra(1:n) into ascending order using heapsort algorithm.
! n is input, ra is replaced on output by its sorted rearrangement. !! \(\text{n}\) is input, \(\text{ra}\) is replaced on output by its
! create an index table (ind) by making an exchange in the index array !! sorted rearrangement.
! whenever an exchange is made on the sorted data array (ra). !! Create an index table (ind) by making an exchange in the index array
! in case of equal values in the data array (ra) the values in the !! whenever an exchange is made on the sorted data array (ra).
! index array (ind) are used to order the entries. !! In case of equal values in the data array (ra) the values in the
! if on input ind(1) = 0 then indices are initialized in the routine, !! index array (ind) are used to order the entries.
! if on input ind(1) != 0 then indices are assumed to have been !! If on input ind(1) = 0 then indices are initialized in the routine,
! initialized before entering the routine and these !! If on input ind(1) != 0 then indices are assumed to have been
! indices are carried around during the sorting process !! initialized before entering the routine and these indices are carried
!! around during the sorting process.
! !
! no work space needed ! ! no work space needed !
! free us from machine-dependent sorting-routines ! ! free us from machine-dependent sorting-routines !
@ -251,16 +253,17 @@ end subroutine hpsort
!--------------------------------------------------------------------- !---------------------------------------------------------------------
subroutine ihpsort (n, ia, ind) subroutine ihpsort (n, ia, ind)
!--------------------------------------------------------------------- !---------------------------------------------------------------------
! sort an integer array ia(1:n) into ascending order using heapsort algorithm. !! Sort an integer array ia(1:n) into ascending order using heapsort
! n is input, ia is replaced on output by its sorted rearrangement. !! algorithm. \(\text{n}\) is input, \(\text{ia}\) is replaced on output
! create an index table (ind) by making an exchange in the index array !! by its sorted rearrangement.
! whenever an exchange is made on the sorted data array (ia). !! Create an index table (ind) by making an exchange in the index array
! in case of equal values in the data array (ia) the values in the !! whenever an exchange is made on the sorted data array (ia).
! index array (ind) are used to order the entries. !! in case of equal values in the data array (ia) the values in the
! if on input ind(1) = 0 then indices are initialized in the routine, !! index array (ind) are used to order the entries.
! if on input ind(1) != 0 then indices are assumed to have been !! If on input ind(1) = 0 then indices are initialized in the routine,
! initialized before entering the routine and these !! If on input ind(1) != 0 then indices are assumed to have been
! indices are carried around during the sorting process !! initialized before entering the routine and these indices are carried
!! around during the sorting process.
! !
! no work space needed ! ! no work space needed !
! free us from machine-dependent sorting-routines ! ! free us from machine-dependent sorting-routines !