Fortran调用C语言排序函数qsort

module sort
    implicit none

    interface swap !数据交换子程序通用接口 公有
        procedure :: swap_int !交换整数
        procedure :: swap_real !交换实数
    end interface
contains
!NR快速排序法 默认升序
subroutine sort_quick_real( arr, id ) !模块方法
    real, dimension(:), intent(inout) :: arr
    integer, dimension(:), intent(out), allocatable, optional :: id !排序前对应的索引值
    !最小值对应排序前下标为id(1) 最大值对应排序前下标为id(size(id))
    integer, parameter :: NN=15, NSTACK=50 !第3版修改nn=7, nstack=64 第2版为nn=15, nstack=50
    real :: a !临时变量
    integer :: i,j,k,jstack,l,r
    integer, dimension(NSTACK) :: istack
    integer, pointer :: b => null() !临时变量

    if ( present(id) ) then
      allocate(b,source=0) !分配内存
      id = [(i, i=1, size(arr))] !初值即排序前下标
    end if

    jstack=0
    l=1
    r=size(arr)

    do
      if (r-l < NN) then
        do j=l+1,r
          a=arr(j)
          if ( present(id) ) b=id(j) !新增
          do i=j-1,l,-1
            if (arr(i) < a) exit
            arr(i+1)=arr(i)
            if ( present(id) ) id(i+1)=id(i) !新增
          end do
          arr(i+1)=a
          if ( present(id) ) id(i+1)=b !新增
        end do
        if (jstack == 0) return
        r=istack(jstack)
        l=istack(jstack-1)
        jstack=jstack-2
      else
        k=(l+r)/2
        call swap(arr(k),arr(l+1))
        if ( present(id) ) call swap(id(k),id(l+1)) !交换元素索引值

        if(arr(r)<arr(l))then !小于号已重载
          call swap(arr(l),arr(r))
          if ( present(id) ) call swap(id(l),id(r)) !交换元素索引值
        end if

        if(arr(r)<arr(l+1))then !小于号已重载
          call swap(arr(l+1),arr(r))
          if ( present(id) ) call swap(id(l+1),id(r)) !交换元素索引值
        end if

        if(arr(l+1)<arr(l))then !小于号已重载
          call swap(arr(l),arr(l+1))
          if ( present(id) ) call swap(id(l),id(l+1)) !交换元素索引值
        end if
        i=l+1
        j=r
        a=arr(l+1)
        if ( present(id) ) b=id(l+1) !新增
        do
          do
            i=i+1
            if (a < arr(i)) exit !小于号已重载
          end do
          do
            j=j-1
            if (arr(j) < a) exit !小于号已重载
          end do
          if (j < i) exit
          call swap(arr(i),arr(j))
          if ( present(id) ) call swap(id(i),id(j)) !交换元素索引值
        end do
        arr(l+1)=arr(j)
        if ( present(id) ) id(l+1)=id(j) !新增
        arr(j)=a
        if ( present(id) ) id(j)=b !新增
        jstack=jstack+2
        if (jstack > NSTACK) error stop 'sort: NSTACK too small'

        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

    if ( present(id) ) deallocate(b) !释放内存
end subroutine sort_quick_real

subroutine swap_int(a,b) !交换整数
    integer,intent(inout)::a,b
    integer::temp

    temp=b
    b=a
    a=temp
end subroutine swap_int

subroutine swap_real(a,b) !交换浮点数
    real,intent(inout)::a,b
    real::temp

    temp=b
    b=a
    a=temp
end subroutine swap_real
end module sort

program test_sort
    use sort
    use iso_c_binding
    implicit none

    interface
        !void* qsort(void* arr,c_size_t n,c_size_t size,int (*comp)(const void* ,const void*));
        subroutine qsort(arr,n,size,comp) bind(c) !C语言qsort函数
            use iso_c_binding
            type(c_ptr),value::arr !void* arr
            integer(c_size_t),value::n,size !元素个数 元素所占大小
            type(c_funptr),intent(in),value::comp !C语言函数指针 用于比较
        end subroutine

        function comp_real(a,b) bind(c) !浮点数比较
            use iso_c_binding
            real(c_float),intent(in)::a,b
            integer(c_int)::comp_real
        end function

        function comp_int(a,b) bind(c) !整数比较
            use iso_c_binding
            integer(c_int),intent(in)::a,b
            integer(c_int)::comp_int
        end function
    end interface

    integer,parameter::N=5
    real::d(N) !原始数组
    real,target::s(N) !排序后数组 c_loc()必需
    integer,target::iarry(N) !排序后数组 c_loc()必需
    real::t1,t2
    type(c_funptr)::f_ptr !c语言函数指针

    call init_seed(999) !固定种子 每次计算时随机数组相同
    call RANDOM_NUMBER(d) !生成随机数组
    !以下实数类型
    s=d
    call cpu_time(t1)
    call sort_quick_real(s) !NR快速排序法
    call cpu_time(t2)
    write(*,*)"NR quick_sort:",t2-t1
    write(*,*)s

    s=d
    f_ptr=c_funloc(comp_real) !获取Fortran函数的C地址
    call cpu_time(t1)
    call qsort(c_loc(s),int(N,kind=c_size_t),c_sizeof(s(1)),f_ptr) !C语言标准库qsort 混合编程
    call cpu_time(t2)
    write(*,*)"qsort_float:",t2-t1
    write(*,*)s

    iarry=[3,2,5,1,4]
    f_ptr=c_funloc(comp_int) !获取Fortran函数的C地址
    call cpu_time(t1)
    call qsort(c_loc(iarry),int(N,kind=c_size_t),c_sizeof(iarry(1)),f_ptr) !C语言标准库qsort 混合编程
    call cpu_time(t2)
    write(*,*)"qsort_int:",t2-t1
    write(*,*)iarry

    read(*,*)
contains
subroutine init_seed(i) !固定种子
    integer,intent(in)::i
    integer,allocatable::seed(:)
    integer::nn

    call random_seed(size=nn)
    allocate(seed(nn))
    seed = i
    call random_seed(put=seed)
    deallocate(seed)
end subroutine
end program test_sort

function comp_real(a,b) bind(c) !浮点数比较
    use iso_c_binding
    implicit none
    real(c_float),intent(in)::a,b
    integer(c_int)::comp_real

    comp_real=0
    if(a>b)then
        comp_real=1 !升序
    else if(a<b)then
        comp_real=-1
    end if
end function

function comp_int(a,b) bind(c) !整数比较
    use iso_c_binding
    implicit none
    integer(c_int),intent(in)::a,b
    integer(c_int)::comp_int

    comp_int=0
    if(a>b)then
        comp_int=1 !升序
    else if(a<b)then
        comp_int=-1
    end if
end function