TABLE OF CONTENTS


ABINIT/m_sort [ Modules ]

[ Top ] [ Modules ]

NAME

 m_sort

FUNCTION

 Sorting algorithms.

COPYRIGHT

  Copyright (C) 2008-2018 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

17 #if defined HAVE_CONFIG_H
18 #include "config.h"
19 #endif
20 
21 #include "abi_common.h"
22 
23 module m_sort
24 
25  use defs_basis, only : std_out
26  use m_errors
27 
28  implicit none
29 
30  private 
31 
32  public :: sort_dp       ! Sort double precision array
33  public :: sort_int      ! Sort integer array
34 
35 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        intent(in)    dimension of the list
  tol      intent(in)    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 given the right ascending order

PARENTS

      atomden,cpdrv,critics,denfgr,finddistrproc,invacuum,listkk,m_bz_mesh
      m_chi0,m_cut3d,m_exc_diago,m_gsphere,m_ifc,m_io_screening
      m_paw_pwaves_lmn,m_phonons,m_polynomial_coeff,m_screen,m_skw,m_use_ga
      m_vcoul,mkcore,mkcore_inner,mkcore_wvl,mlwfovlp_qp,outscfcv
      partial_dos_fractions,shellstruct,symkpt,tddft,thmeig,wvl_initro

CHILDREN

SOURCE

 69 subroutine sort_dp(n,list,iperm,tol)
 70 
 71 
 72 !This section has been created automatically by the script Abilint (TD).
 73 !Do not modify the following lines by hand.
 74 #undef ABI_FUNC
 75 #define ABI_FUNC 'sort_dp'
 76 !End of the abilint section
 77 
 78  implicit none
 79 
 80 !Arguments ------------------------------------
 81 !scalars
 82  integer, intent(in) :: n
 83  integer, intent(inout) :: iperm(n)
 84  double precision, intent(inout) :: list(n)
 85  double precision, intent(in) :: tol
 86 
 87 !Local variables-------------------------------
 88 !scalars
 89  integer :: l,ir,iap,i,j
 90  double precision :: ap
 91 
 92  if (n==1) then
 93 
 94 ! Accomodate case of array of length 1: already sorted!
 95   return
 96 
 97  else if (n<1) then
 98 
 99 ! Should not call with n<1
100   write(std_out,1000) n
101   1000  format(/,' sort_dp has been called with array length n=',i12,/, &
102 &  ' having a value less than 1.  This is not allowed.')
103   MSG_ERROR("Aborting now")
104 
105  else ! n>1
106 
107 ! Conduct the usual sort
108 
109   l=n/2+1
110   ir=n
111 
112   do   ! Infinite do-loop
113 
114    if (l>1) then
115 
116     l=l-1
117     ap=list(l)
118     iap=iperm(l)
119 
120    else ! l<=1
121 
122     ap=list(ir)
123     iap=iperm(ir)
124     list(ir)=list(1)
125     iperm(ir)=iperm(1)
126     ir=ir-1
127 
128     if (ir==1) then
129      list(1)=ap
130      iperm(1)=iap
131      exit   ! This is the end of this algorithm
132     end if
133 
134    end if ! l>1
135 
136    i=l
137    j=l+l
138 
139    do while (j<=ir) 
140     if (j<ir) then
141      if ( list(j)<list(j+1)-tol .or.  &
142 &        (list(j)<list(j+1)+tol.and.iperm(j)<iperm(j+1))) j=j+1
143     endif
144     if (ap<list(j)-tol .or. (ap<list(j)+tol.and.iap<iperm(j))) then
145      list(i)=list(j)
146      iperm(i)=iperm(j)
147      i=j
148      j=j+j
149     else
150      j=ir+1
151     end if
152    enddo
153 
154    list(i)=ap
155    iperm(i)=iap
156 
157   enddo ! End infinite do-loop
158 
159  end if ! n>1
160 
161 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        intent(in)    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

PARENTS

      getng,getngrec,initmpi_img,invars2,irrzg,m_dvdb,m_hdr,m_nesting,m_wfk
      mkfskgrid,shellstruct

CHILDREN

SOURCE

191 subroutine sort_int(n,list,iperm)
192 
193 
194 !This section has been created automatically by the script Abilint (TD).
195 !Do not modify the following lines by hand.
196 #undef ABI_FUNC
197 #define ABI_FUNC 'sort_int'
198 !End of the abilint section
199 
200  implicit none
201 
202 !Arguments ------------------------------------
203 !scalars
204  integer,intent(in) :: n 
205  integer,intent(inout) :: list(n),iperm(n)
206 
207 !Local variables-------------------------------
208 !scalars
209  integer :: l,ir,i,j,ip,ipp
210 ! *************************************************************************
211 
212  if (n==1) then
213 
214 ! Accomodate case of array of length 1: already sorted!
215   return
216 
217  else if (n<1) then
218 
219 ! Should not call with n<1
220   write(std_out,1000) n
221   1000  format(/,' sort_int has been called with array length n=',i12,/, &
222 &  ' having a value less than 1.  This is not allowed.')
223   MSG_ERROR("Aborting now")
224 
225  else ! n>1
226 
227 ! Conduct the usual sort
228 
229   l=n/2+1
230   ir=n
231 
232   do   ! Infinite do-loop
233  
234    if (l>1) then
235 
236     l=l-1
237     ip=list(l)
238     ipp=iperm(l)
239 
240    else
241 
242     ip=list(ir)
243     ipp=iperm(ir)
244     list(ir)=list(1)
245     iperm(ir)=iperm(1)
246     ir=ir-1
247 
248     if (ir==1) then
249      list(1)=ip
250      iperm(1)=ipp
251      exit   ! This is the end of this algorithm
252     end if
253 
254    end if ! l>1
255 
256    i=l
257    j=l+l
258 
259    do while (j<=ir)
260     if (j<ir) then
261      if (list(j).lt.list(j+1)) j=j+1
262     end if
263     if (ip.lt.list(j)) then
264      list(i)=list(j)
265      iperm(i)=iperm(j)
266      i=j
267      j=j+j
268     else
269      j=ir+1
270     end if
271    enddo
272 
273    list(i)=ip
274    iperm(i)=ipp
275 
276   enddo ! End infinite do-loop
277 
278  end if ! n>1
279 
280 end subroutine sort_int