TABLE OF CONTENTS
ABINIT/m_sort [ Modules ]
NAME
m_sort
FUNCTION
Sorting algorithms.
COPYRIGHT
Copyright (C) 2008-2022 ABINIT group (XG) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 module m_sort 23 24 use defs_basis 25 use m_errors 26 use m_profiling_abi 27 28 implicit none 29 30 private 31 32 ! Low-level routines 33 public :: sort_dp ! Sort double precision array 34 public :: sort_int ! Sort integer array 35 36 ! Helper functions to perform common operations. 37 public :: sort_rpts ! Sort list of real points by |r| 38 public :: sort_rvals ! Out-of-place sort of real values 39 40 CONTAINS !====================================================================================================
m_sort/sort_dp [ Functions ]
[ Top ] [ m_sort ] [ Functions ]
NAME
sort_dp
FUNCTION
Sort double precision array list(n) into ascending numerical order using Heapsort algorithm, while making corresponding rearrangement of the integer array iperm. Consider that two double precision numbers within tolerance tol are equal.
INPUTS
n: dimension of the list tol: numbers within tolerance are equal list(n) intent(inout) list of double precision numbers to be sorted iperm(n) intent(inout) iperm(i)=i (very important)
OUTPUT
list(n) sorted list iperm(n) index of permutation giving the right ascending order: the i-th element of the ouput ordered list had index iperm(i) in the input list.
SOURCE
65 subroutine sort_dp(n, list, iperm, tol) 66 67 !Arguments ------------------------------------ 68 !scalars 69 integer, intent(in) :: n 70 integer, intent(inout) :: iperm(n) 71 real(dp), intent(inout) :: list(n) 72 real(dp), intent(in) :: tol 73 74 !Local variables------------------------------- 75 !scalars 76 integer :: l,ir,iap,i,j 77 real(dp) :: ap 78 character(len=500) :: msg 79 80 if (n==1) then 81 82 ! Accomodate case of array of length 1: already sorted! 83 return 84 85 else if (n<1) then 86 87 ! Should not call with n<1 88 write(msg, "(a,i0,2a)")& 89 "sort_dp has been called with array length n= ",n, ch10, & 90 "having a value less than 1. This is not allowed." 91 ABI_ERROR(msg) 92 93 else ! n>1 94 95 ! Conduct the usual sort 96 97 l=n/2+1 98 ir=n 99 100 do ! Infinite do-loop 101 102 if (l>1) then 103 104 l=l-1 105 ap=list(l) 106 iap=iperm(l) 107 108 else ! l<=1 109 110 ap=list(ir) 111 iap=iperm(ir) 112 list(ir)=list(1) 113 iperm(ir)=iperm(1) 114 ir=ir-1 115 116 if (ir==1) then 117 list(1)=ap 118 iperm(1)=iap 119 exit ! This is the end of this algorithm 120 end if 121 122 end if ! l>1 123 124 i=l 125 j=l+l 126 127 do while (j<=ir) 128 if (j<ir) then 129 if ( list(j)<list(j+1)-tol .or. & 130 & (list(j)<list(j+1)+tol.and.iperm(j)<iperm(j+1))) j=j+1 131 endif 132 if (ap<list(j)-tol .or. (ap<list(j)+tol.and.iap<iperm(j))) then 133 list(i)=list(j) 134 iperm(i)=iperm(j) 135 i=j 136 j=j+j 137 else 138 j=ir+1 139 end if 140 enddo 141 142 list(i)=ap 143 iperm(i)=iap 144 145 enddo ! End infinite do-loop 146 147 end if ! n>1 148 149 end subroutine sort_dp
m_sort/sort_int [ Functions ]
[ Top ] [ m_sort ] [ Functions ]
NAME
sort_int
FUNCTION
Sort integer array list(n) into ascending numerical order using Heapsort algorithm, while making corresponding rearrangement of the integer array iperm.
INPUTS
n: dimension of the list list(n) intent(inout) list of double precision numbers to be sorted iperm(n) intent(inout) iperm(i)=i (very important)
OUTPUT
list(n): sorted list iperm(n): index of permutation given the right ascending order the i-th element of the ouput ordered list had index iperm(i) in the input list.
SOURCE
172 subroutine sort_int(n, list, iperm) 173 174 !Arguments ------------------------------------ 175 !scalars 176 integer,intent(in) :: n 177 integer,intent(inout) :: list(n),iperm(n) 178 179 !Local variables------------------------------- 180 !scalars 181 integer :: l,ir,i,j,ip,ipp 182 character(len=500) :: msg 183 ! ************************************************************************* 184 185 if (n==1) then 186 187 ! Accomodate case of array of length 1: already sorted! 188 return 189 190 else if (n<1) then 191 192 ! Should not call with n<1 193 write(msg, "(a,i0,2a)")& 194 "sort_int has been called with array length n= ",n, ch10, & 195 "having a value less than 1. This is not allowed." 196 ABI_ERROR(msg) 197 198 else ! n>1 199 200 ! Conduct the usual sort 201 202 l=n/2+1 203 ir=n 204 205 do ! Infinite do-loop 206 207 if (l>1) then 208 209 l=l-1 210 ip=list(l) 211 ipp=iperm(l) 212 213 else 214 215 ip=list(ir) 216 ipp=iperm(ir) 217 list(ir)=list(1) 218 iperm(ir)=iperm(1) 219 ir=ir-1 220 221 if (ir==1) then 222 list(1)=ip 223 iperm(1)=ipp 224 exit ! This is the end of this algorithm 225 end if 226 227 end if ! l>1 228 229 i=l 230 j=l+l 231 232 do while (j<=ir) 233 if (j<ir) then 234 if (list(j).lt.list(j+1)) j=j+1 235 end if 236 if (ip.lt.list(j)) then 237 list(i)=list(j) 238 iperm(i)=iperm(j) 239 i=j 240 j=j+j 241 else 242 j=ir+1 243 end if 244 enddo 245 246 list(i)=ip 247 iperm(i)=ipp 248 249 enddo ! End infinite do-loop 250 251 end if ! n>1 252 253 end subroutine sort_int
m_sort/sort_rpts [ Functions ]
[ Top ] [ m_sort ] [ Functions ]
NAME
sort_rpts
FUNCTION
Sort list of real space 3d-points by norm (ascending order) Input list is not modified.
INPUTS
n: dimension of the list rpts(3, n): points in reduced coordinates metric: Metric used to compute |r|. [tol]: numbers within tolerance are equal.
OUTPUT
iperm(n) index of permutation giving the right ascending order: the i-th element of the ordered list had index iperm(i) in rpts. [rmod(n)]= list of sorted |r| values.
SOURCE
279 subroutine sort_rpts(n, rpts, metric, iperm, tol, rmod) 280 281 !Arguments ------------------------------------ 282 !scalars 283 integer,intent(in) :: n 284 integer,allocatable,intent(out) :: iperm(:) 285 real(dp),optional,allocatable,intent(out) :: rmod(:) 286 real(dp),optional,intent(in) :: tol 287 !arrays 288 real(dp),intent(in) :: rpts(3,n), metric(3,3) 289 290 !Local variables------------------------------- 291 !scalars 292 integer :: ii 293 real(dp) :: my_tol 294 !arrays 295 real(dp),allocatable :: my_rmod(:) 296 297 !************************************************************************ 298 299 my_tol = tol12; if (present(tol)) my_tol = tol 300 301 ABI_MALLOC(my_rmod, (n)) 302 do ii=1,n 303 my_rmod(ii) = sqrt(dot_product(rpts(:, ii), matmul(metric, rpts(:, ii)))) 304 end do 305 ABI_MALLOC(iperm, (n)) 306 iperm = [(ii, ii=1,n)] 307 call sort_dp(n, my_rmod, iperm, my_tol) 308 309 if (present(rmod)) then 310 call move_alloc(my_rmod, rmod) 311 else 312 ABI_FREE(my_rmod) 313 end if 314 315 end subroutine sort_rpts
m_sort/sort_rvals [ Functions ]
[ Top ] [ m_sort ] [ Functions ]
NAME
sort_rvals
FUNCTION
Sort list of real values (ascending order) Input list is not modified.
INPUTS
n: dimension of the list in_vals(n): input weigts. [tol]: tolerance for comparison
OUTPUT
iperm(n) index of permutation giving the right ascending order: the i-th element of the ordered list had index iperm(i) in in_vals. [sorted_in_vals(n)]= list of sorted weigts.
SOURCE
340 subroutine sort_rvals(n, in_vals, iperm, sorted_vals, tol) 341 342 !Arguments ------------------------------------ 343 !scalars 344 integer,intent(in) :: n 345 real(dp),intent(in) :: in_vals(n) 346 integer,allocatable,intent(out) :: iperm(:) 347 real(dp),allocatable,intent(out) :: sorted_vals(:) 348 real(dp),optional,intent(in) :: tol 349 350 !Local variables------------------------------- 351 !scalars 352 integer :: ii 353 real(dp) :: my_tol 354 355 !************************************************************************ 356 357 my_tol = tol12; if (present(tol)) my_tol = tol 358 359 ABI_MALLOC(sorted_vals, (n)) 360 sorted_vals = in_vals 361 ABI_MALLOC(iperm, (n)) 362 iperm = [(ii, ii=1,n)] 363 call sort_dp(n, sorted_vals, iperm, my_tol) 364 365 end subroutine sort_rvals