Cpython 和 fortran 链表之间的互操作性

Question

我有一个大致类似的fortran链表

type :: node
    type(node), pointer :: next => null()
    integer :: value
end type node

理想情况下，我希望使用 Cpython 与之交互。我使用 f2py 程序在 python 中使用了 fortran 子例程来创建共享对象。但是，f2py 不能与派生类型一起使用。

我的问题很简单，是否可以使用 cpython 在 Fortran 中访问链表之类的东西。我想我需要遵循 fortran 到 c 到 cpython 的路线。但是，我已经读到，对于与 c 互操作的 fortran 派生类型，“每个组件都必须具有可互操作的类型和类型参数，不能是指针，也不能是可分配的。”同样，c-fortran interoperability - derived types with pointers 的帖子似乎也证实了这一点。

我想知道是否有人知道是否绝对不可能从 cpython 访问 fortran 中的链表。如果可能的话，即使是间接的或迂回的，我也会很感激听到更多。

谢谢你， 标记

Answer 1

正如 Bálint Aradi 在 cmets 中已经提到的，该节点不能与当前形式的 C 互操作。为此，您需要将 fortran 指针更改为 C 指针，但这使得在 fortran 本身内部使用非常痛苦。我能想到的最优雅的解决方案是将 C 可互操作类型放入 fortran 类型中，并保存不同版本的 C 和 fortran 指针。

实现如下所示，其中我还定义了用于在 fortran 中分配、解除分配和初始化节点的便利函数。

module node_mod

    use, intrinsic :: iso_c_binding
    implicit none

    ! the C interoperable type
    type, bind(c) :: cnode
        type(c_ptr) :: self = c_null_ptr
        type(c_ptr) :: next = c_null_ptr
        integer(c_int) :: value
    end type cnode

    ! the type used for work in fortran
    type :: fnode
        type(cnode) :: c
        type(fnode), pointer :: next => null()
    end type fnode

contains

recursive function allocate_nodes(n, v) result(node)

    integer, intent(in) :: n
    integer, optional, intent(in) :: v
    type(fnode), pointer :: node
    integer :: val

    allocate(node)
    if (present(v)) then
        val = v
    else
        val = 1
    end if
    node%c%value = val
    if (n > 1) then
        node%next => allocate_nodes(n-1, val+1)
    end if

end function allocate_nodes

recursive subroutine deallocate_nodes(node)

    type(fnode), pointer, intent(inout) :: node
    if (associated(node%next)) then
        call deallocate_nodes(node%next)
    end if
    deallocate(node)

end subroutine deallocate_nodes

end module node_mod

如您所见，访问“value”元素需要额外的“%c”，这有点麻烦。要在 python 中使用之前定义的例程来检索链表，必须定义 C 互操作包装器并且必须链接 C 指针。

module node_mod_cinter

    use, intrinsic :: iso_c_binding
    use, non_intrinsic :: node_mod

    implicit none

contains

recursive subroutine set_cptr(node)

    type(fnode), pointer, intent(in) :: node

    node%c%self = c_loc(node)
    if (associated(node%next)) then
        node%c%next = c_loc(node%next%c)
        call set_cptr(node%next)
    end if

end subroutine set_cptr

function allocate_nodes_citer(n) bind(c, name="allocate_nodes") result(cptr)

    integer(c_int), value, intent(in) :: n
    type(c_ptr) :: cptr
    type(fnode), pointer :: node

    node => allocate_nodes(n)
    call set_cptr(node)
    cptr = c_loc(node%c)

end function allocate_nodes_citer

subroutine deallocate_nodes_citer(cptr) bind(c, name="deallocate_nodes")

    type(c_ptr), value, intent(in) :: cptr
    type(cnode), pointer :: subnode
    type(fnode), pointer :: node

    call c_f_pointer(cptr, subnode)
    call c_f_pointer(subnode%self, node)
    call deallocate_nodes(node)

end subroutine deallocate_nodes_citer

end module node_mod_cinter

“*_nodes_citer”例程只处理不同的指针类型，set_cptr 子例程根据 fortran 指针链接 C 互操作类型内部的 C 指针。我添加了 node%c%self 元素，以便可以恢复 fortran 指针并将其用于正确的释放，但如果您不太关心这一点，那么它并不是绝对需要的。

此代码需要编译为共享库以供其他程序使用。我在我的 linux 机器上使用了以下带有 gfortran 的命令。

gfortran -fPIC -shared -o libnode.so node.f90

最后，分配 10 个节点的列表的 python 代码，打印出每个节点的 node%c%value，然后再次释放所有内容。此外，还显示了 fortran 和 C 节点的内存位置。

#!/usr/bin/env python
import ctypes
from ctypes import POINTER, c_int, c_void_p
class Node(ctypes.Structure):
    pass
Node._fields_ = (
        ("self", c_void_p),
        ("next", POINTER(Node)),
        ("value", c_int),
        )

def define_function(res, args, paramflags, name, lib):

    prot = ctypes.CFUNCTYPE(res, *args)
    return prot((name, lib), paramflags)

def main():

    import os.path

    libpath = os.path.abspath("libnode.so")
    lib = ctypes.cdll.LoadLibrary(libpath)

    allocate_nodes = define_function(
            res=POINTER(Node),
            args=(
                c_int,
                ),
            paramflags=(
                (1, "n"),
                ),
            name="allocate_nodes",
            lib=lib,
            )

    deallocate_nodes = define_function(
            res=None,
            args=(
                POINTER(Node),
                ),
            paramflags=(
                (1, "cptr"),
                ),
            name="deallocate_nodes",
            lib=lib,
            )

    node_ptr = allocate_nodes(10)

    n = node_ptr[0]
    print "value", "f_ptr", "c_ptr"
    while True:
        print n.value, n.self, ctypes.addressof(n)
        if n.next:
            n = n.next[0]
        else:
            break

    deallocate_nodes(node_ptr)

if __name__ == "__main__":
    main()

执行此操作会得到以下输出：

value f_ptr c_ptr
1 15356144 15356144
2 15220144 15220144
3 15320384 15320384
4 14700384 14700384
5 15661152 15661152
6 15661200 15661200
7 15661248 15661248
8 14886672 14886672
9 14886720 14886720
10 14886768 14886768

有趣的是，两种节点类型都从相同的内存位置开始，所以 node%c%self 并不是真正需要的，但这只是因为我对类型定义很小心，这真的不应该指望。

你有它。即使不必处理链表也很麻烦，但是 ctypes 比 f2py 更强大和健壮得多。希望这会带来一些好处。