list - Fortran 中的排序函数
问题描述
是否有一个 Fortran 库具有排序函数的实现,即一个函数ordering(list)
(如MathematicaOrdering[]
中的函数),它给出排序列表的每个连续元素出现的位置?list
我可以实现它,但我不想重新发明轮子(我的轮子可能远非完美......)。由于它是如此基本,我正在搜索包含此类列表操作的库,但未能找到。
你有什么建议吗?
解决方案
由于我很久以前就已经实现了这个(它依赖并大量借鉴了 Bill Press 等人的数字食谱书),这里有一个在 Fortran 中的独立实现:
module index_mod
use, intrinsic :: iso_fortran_env, only: IK=>int32, RK=>real64
implicit none
contains
subroutine indexArrayReal(n,Array,Index)
implicit none
integer(IK), intent(in) :: n
real(RK) , intent(in) :: Array(n)
integer(IK), intent(out) :: Index(n)
integer(IK), parameter :: nn=15, nstack=50
integer(IK) :: k,i,j,indext,jstack,l,r
integer(IK) :: istack(nstack)
real(RK) :: a
do j = 1,n
Index(j) = j
end do
jstack=0
l=1
r=n
do
if (r-l < nn) then
do j=l+1,r
indext=Index(j)
a=Array(indext)
do i=j-1,l,-1
if (Array(Index(i)) <= a) exit
Index(i+1)=Index(i)
end do
Index(i+1)=indext
end do
if (jstack == 0) return
r=istack(jstack)
l=istack(jstack-1)
jstack=jstack-2
else
k=(l+r)/2
call swap(Index(k),Index(l+1))
call exchangeIndex(Index(l),Index(r))
call exchangeIndex(Index(l+1),Index(r))
call exchangeIndex(Index(l),Index(l+1))
i=l+1
j=r
indext=Index(l+1)
a=Array(indext)
do
do
i=i+1
if (Array(Index(i)) >= a) exit
end do
do
j=j-1
if (Array(Index(j)) <= a) exit
end do
if (j < i) exit
call swap(Index(i),Index(j))
end do
Index(l+1)=Index(j)
Index(j)=indext
jstack=jstack+2
if (jstack > nstack) then
write(*,*) 'NSTACK too small in indexArrayReal()' ! xxx
error stop
end if
if (r-i+1 >= j-l) then
istack(jstack)=r
istack(jstack-1)=i
r=j-1
else
istack(jstack)=j-1
istack(jstack-1)=l
l=i
end if
end if
end do
contains
subroutine exchangeIndex(i,j)
integer(IK), intent(inout) :: i,j
integer(IK) :: swp
if (Array(j) < Array(i)) then
swp=i
i=j
j=swp
end if
end subroutine exchangeIndex
pure elemental subroutine swap(a,b)
implicit none
integer(IK), intent(inout) :: a,b
integer(IK) :: dum
dum=a
a=b
b=dum
end subroutine swap
end subroutine indexArrayReal
end module Index_mod
program Index_prog
use Index_mod, only: IK, RK, indexArrayReal
implicit none
integer(IK), parameter :: n = 5
integer(IK) :: Index(n)
real(RK) :: Array(n) = [ 1.,3.,4.,2.,-1. ]
call indexArrayReal(n,Array,Index)
write(*,*) "Index: ", Index
write(*,*) "Array(Index): ", Array(Index)
end program Index_prog
使用 GFortran 2008 编译,输出如下:
$gfortran -std=f2008 *.f95 -o main
$main
Index: 5 1 4 2 3
Array(Index): -1.0000000000000000 1.0000000000000000 2.0000000000000000 3.0000000000000000 4.0000000000000000
上述例程用于对实值数组进行排序。要对整数数组进行排序,只需更改toreal(RK) :: Array(n)
的接口subroutine indexArrayReal()
即可integer(IK) :: Array(n)
。