如何在C中访问(动态分配)Fortran数组
我的主要问题是为什么数组做了如此奇怪的事情以及是否有任何方式以“干净”的方式执行以下操作。
我目前有一个C程序foo.c
通过dlopen/dlsym
连接Fortran程序bar.f90
,大致类似于下面的代码:
foo.c的:
#include #include int main() { int i, k = 4; double arr[k]; char * e; void * bar = dlopen("Code/Test/bar.so", RTLD_NOW | RTLD_LOCAL); void (*allocArray)(int*); *(void **)(&allocArray) = dlsym(bar, "__bar_MOD_allocarray"); void (*fillArray)(double*); *(void **)(&fillArray) = dlsym(bar, "__bar_MOD_fillarray"); void (*printArray)(void); *(void **)(&printArray) = dlsym(bar, "__bar_MOD_printarray"); double *a = (double*)dlsym(bar, "__bar_MOD_a"); for(i = 0; i < k; i++) arr[i] = i * 3.14; (*allocArray)(&k); (*fillArray)(arr); (*printArray)(); for(i = 0; i < 4; i++) printf("%f ", a[i]); printf("\n"); return 0; }
bar.f90:
module bar integer, parameter :: pa = selected_real_kind(15, 307) real(pa), dimension(:), allocatable :: a integer :: as contains subroutine allocArray(asize) integer, intent(in) :: asize as = asize allocate(a(asize)) return end subroutine subroutine fillArray(values) real(pa), dimension(as), intent(in) :: values a = values return end subroutine subroutine printArray() write(*,*) a return end subroutine end module
运行主要收益率
0.0000000000000000 3.1400000000000001 6.2800000000000002 9.4199999999999999 0.000000 -nan 0.000000 0.000000
这表明Fortran正确地分配了数组,甚至正确地存储了给定的值,但它们不再可以通过dlsym访问(处理该数据导致段错误)。 我也试过这个固定大小的数组 – 结果保持不变。
有谁知道这种行为的原因? 我个人认为,事情要么是双向的,要么根本不是 – 这个“Fortran接受Carrays而不是反之亦然”让我想知道我是否以这种方式从C访问arrays时遇到了一些基本错误。
另一个(甚至更重要的)问题是,如何像这些“正确的方式”那样进行数组访问。 目前我甚至不确定是否坚持使用“Fortran as .so”界面是一个好方法 – 我认为在这种情况下也可以尝试混合编程。 尽管如此,arrays问题仍然存在 – 我读到这可以通过使用ISO C绑定以某种方式解决,但我无法弄清楚如何,(我还没有与Fortran一起工作,但是,特别是没有使用Binding) ,所以对此问题的帮助将不胜感激。
编辑:
好的,所以我读了一下ISO C Binding,并在这里找到了一个非常有用的方法。 使用C_LOC
我可以获得指向我的Fortran结构的C指针。 不幸的是,指向数组的指针似乎是指针的指针,需要在它们作为C数组处理之前在C代码中取消引用 – 或类似的东西。
编辑:
让我的程序现在使用C绑定,就像Vladimir F指出的那样,至少在大多数情况下。 C文件和Fortran文件现在链接在一起,所以我可以避免使用libdl接口,至少对于Fortran部分 – 我仍然需要加载动态C库,获取指向其中一个符号的函数指针并传递作为Fortran的函数指针,后者将该函数作为其计算的一部分调用。 由于所述函数需要double * s [数组],我无法使用C_LOC设法传递我的Fortran数组,奇怪的是 – C_LOC(array)
和C_LOC(array(1))
将正确的指针传递回C函数。 array(1)
虽然做了这个技巧。 可悲的是,这不是“最干净”的做法。 如果有人告诉我如何使用C_LOC
函数执行此操作,那将是很好的。 尽管如此,我接受弗拉基米尔F的答案,因为我认为这是更安全的解决方案。
在我看来,尝试访问Fortran库中的全局数据并不是一个好习惯。 它可以使用COMMON块来完成,但它们是邪恶的并且需要静态大小的数组。 通常存储关联是一件坏事。
永远不要将模块符号作为“__bar_MOD_a”访问,它们是编译器特定的,不应直接使用。 使用函数和子例程传递poiters。
将数组作为子例程参数传递。 您还可以在C中分配数组并将其传递给Fortran。 还可以做的是获取指向数组的第一个元素的指针。 它将用于指向数组的C指针。
我的解决方案,为了简单而没有.so,添加它是微不足道的:
bar.f90
module bar use iso_C_binding implicit none integer, parameter :: pa = selected_real_kind(15, 307) real(pa), dimension(:), allocatable,target :: a integer :: as contains subroutine allocArray(asize,ptr) bind(C,name="allocArray") integer, intent(in) :: asize type(c_ptr),intent(out) :: ptr as = asize allocate(a(asize)) ptr = c_loc(a(1)) end subroutine subroutine fillArray(values) bind(C,name="fillArray") real(pa), dimension(as), intent(in) :: values a = values end subroutine subroutine printArray() bind(C,name="printArray") write(*,*) a end subroutine end module
main.c中
#include #include int main() { int i, k = 4; double arr[k]; char * e; double *a; void allocArray(int*,double**); void fillArray(double*); void allocArray(); for(i = 0; i < k; i++) arr[i] = i * 3.14; allocArray(&k,&a); fillArray(arr); printArray(); for(i = 0; i < 4; i++) printf("%f ", a[i]); printf("\n"); return 0; }
编译并运行:
gcc -c -g main.c gfortran -c -g -fcheck=all bar.f90 gfortran main.o bar.o ./a.out 0.0000000000000000 3.1400000000000001 6.2800000000000002 9.4199999999999999 0.000000 3.140000 6.280000 9.420000
注意:Fortran子例程中没有理由返回,它们只是模糊了代码。
许多Fortran编译器在内部使用称为数组描述符的东西 – 包含数组形状的结构(即每个维度的大小和范围以及指向实际数据的指针)。 它允许实现假设形状数组参数,数组指针和可分配数组之类的工作。 您通过__bar_MOD_a
符号访问的是可分配数组的描述符,而不是其数据。
数组描述符是特定于编译器的,依赖于特定描述符格式的代码是不可移植的。 示例描述符:
- GNU Fortran
- 英特尔Fortran
请注意,即使是那些特定于那些编译器的某些版本。 例如,英特尔声称其当前的描述符格式与英特尔Fortran 7.0中使用的格式不兼容。
如果你看两个描述符,你会发现它们非常相似,第一个元素是指向数组数据的指针。 因此,您可以使用double **
而不是double *
轻松读取数据:
double **a_descr = (double**)dlsym(bar, "__bar_MOD_a"); ... for(i = 0; i < 4; i++) printf("%f ", (*a_descr)[i]);
再一次,这是不可移植的,因为这些描述符的格式将来可能会改变(尽管我怀疑数据指针会移动到除描述符开头之外的其他地方)。 有一个规范草案试图统一所有描述符格式,但不清楚它们将如何以及何时被不同的编译器供应商采用。
编辑:以下是如何使用一个访问器函数,该函数使用ISO_C_BINDING
模块中的C_LOC()
来可移植地获取指向可分配数组的指针:
Fortran代码:
module bar use iso_c_binding ... ! Note that the array should be a pointer target real(pa), dimension(:), allocatable, target :: a ... contains ... function getArrayPtr() result(cptr) type(c_ptr) :: cptr cptr = c_loc(a) end function end module
C代码:
... void * (*getArrayPtr)(void); *(void **)(&getArrayPtr) = dlsym(bar, "__bar_MOD_getarrayptr"); ... double *a = (*getArrayPtr)(); for(i = 0; i < 4; i++) printf("%f ", a[i]); ...
结果:
$ ./prog.x 0.0000000000000000 3.1400000000000001 6.2800000000000002 9.4199999999999999 0.000000 3.140000 6.280000 9.420000