C++中的Fortran多维数组

Question

C++中的Fortran多维数组

c++arraysmultidimensional-arrayfortranfortran90

3

我正在尝试在C++ Fortran互操作程序中传递多维Fortran数组给C++程序。我基本知道如何从Fortran传递数组到C++，即从Fortran向C++传递数组的位置。然后，C++接收扁平化的数组，并需要进行一些代数计算以找到给定多维数组中的元素。

我已经成功地在标量数组上测试了这个想法。在C++中找到元素的索引并不难，因为它是从Fortran索引线性映射到C++索引，偏移量为-1。Fortran和C++的示例代码如下：

! Fortran main program
program fprogram

integer :: i
real*8 :: array(2)

array(1) = 1.0
array(2) = 2.0

! call cpp function
call cppfuncarray(array, 2)
write(*,*) array

end program

//cpp file
extern "C" { 
void cppfuncarray_(double *array, int *imax);}

void cppfuncarray_(double *array, int *imax) {
    int iMax = *imax;
    for ( int i = 0; i < iMax; i++ ) {
        array[i] = array + 1.0*i;
    }
};

输出将为array(1) = 1和array(2) =3。现在我正在尝试将多维数组(如A(2,2)或A(2,3,2))从Fortran传递到C++。我知道像A(2,2)这样的二维数组在C++中会很容易地找出压平的数组。但我认为在C++中查找三维或四维数组的元素可能会更具挑战性。

正确的方法是什么，可以在C++中构建多维数组，以便我可以像在Fortran中那样通过A[k][j][i]引用数组A(k,j,i)?

提前致谢！

- J. Brad Maeng

我实际上在考虑更加优雅的解决方案，而不仅仅是最简单的方法。我认为使用节点主序来查找多维数组中的元素可能很简单，但它并不总是最令人满意的结构。如果我无法找到更令人满意的结构，我可能会坚持使用最有效的方法。 - J. Brad Maeng

据我所知，Fortran 数组无论有多少维，都是通过指向第一个元素的指针来引用的。这意味着当我将其传递给 C++ 时，只需要指针和数组大小即可。bg2b，您是否暗示我理解有误？数组是否比我刚才描述的更复杂？ - J. Brad Maeng

@J.BradMaeng 在Fortran 77中曾经是正确的，但自从Fortran 90以后（已经是2015年！甚至90现在也很老了！）这不再是正确的了。然而，当您使用编译器不知道的外部过程时，它将制作一个数组的副本，如果需要，可以用单个指针引用该副本。 - Vladimir F Героям слава

@J.BradMaeng 我猜Fortran >= 90中的可分配数组和数组指针类似于多维数组类。据我所记，一些旧编译器直接将这种“数组类”对象的地址传递给外部例程，导致结果不正确（因为后者期望接收数组的第一个元素的地址）。但是最近的编译器似乎不是这种情况，因此最好检查使用的特定编译器。 - roygvib

@roygvib 再次感谢，今天在这里学到的东西非常有用。 - J. Brad Maeng

显示剩余7条评论

1个回答

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- roygvib · Accepted Answer

将指针强制转换为标量 (在下面的示例中为 int*)，然后将其转换为指向 (N-1) 维数组的指针可能很有用，尽管编写一个数组视图类应该更加灵活。

fortsub.f90:

subroutine fortsub()
    implicit none
    integer a(4,3,2)   !! this line will be changed in the EDIT below
    integer ndims(3), i

    ndims(:) = [ ( size( a, i ), i = 1, 3 ) ]
    call cppfunc( a, ndims )

    print *, "a(1,1,1) = ", a(1,1,1)   !! gives 10101
    print *, "a(2,2,1) = ", a(2,2,1)   !! gives 20201
    print *, "a(4,3,2) = ", a(4,3,2)   !! gives 40302
end subroutine

cppfunc.cpp:

extern "C" {
void fortsub_( void );

void cppfunc_( int *A, int *n )
{
    typedef int (*A3d_t)[ n[1] ][ n[0] ];
    A3d_t A3d = (A3d_t) A;  // get a pointer to 2-dim subarray

    // 3-dim access                                                                  
    for( int k = 0; k < n[2]; k++ )
    for( int j = 0; j < n[1]; j++ )
    for( int i = 0; i < n[0]; i++ ) {
        A3d[ k ][ j ][ i ] = (i+1)*10000 + (j+1)*100 + (k+1); // set test data    
    }
}
}; // extern "C"                                                                     

int main()
{
    fortsub_();
    return 0;
}

编译：

$ g++ fortsub.f90 cppfunc.cpp -lgfortran  # tested with gcc >=4.4.7

编辑：虽然这可能与主题无关，但我也尝试将可分配数组或数组指针传递给同一个cppfunc()例程。具体来说，我将上面的a(4,3,2)声明更改为以下之一：

!! case 1
integer, allocatable :: a(:,:,:)                                               
allocate( a(4,3,2) ) 

!! case 2
integer, pointer :: a(:,:,:)                                                   
allocate( a(4,3,2) )

!! case 3 : an array view for contiguous memory
integer, target :: b(4,3,2,5)
integer, pointer :: a(:,:,:)
a => b( :, :, :, 5 )

!! case 4 : an array view for non-contiguous memory
integer, target :: c(5,4,3,2)                                                  
integer, pointer :: a(:,:,:)                                                   
a => c( 5, :, :, : )

当使用编译器进行编译时

$ g++ fortsub.f90 cppfunc.cpp -lgfortran -fcheck-array-temporaries

所有案例都会给出正确的结果。只有在第4个案例中，编译器会创建一个数组临时变量，将其第一个元素的地址传递给cppfunc()，并将获取的数据复制回实际参数。否则，编译器会像Fortran77一样直接将实际数组的第一个元素的地址传递给cppfunc()。

编辑2：有些程序可能希望将N维数组作为指针数组接收。在这种情况下，更传统的方法将是类似于以下内容：

getptr3d.hpp:

template <typename T>
T*** get_ptr3d( T* A, int* n )
{
    typedef T (*A3d_t)[ n[1] ][ n[0] ];
    A3d_t A3d = (A3d_t) A;

    T*** p = new T** [ n[2] ];

    for( int k = 0; k < n[2]; k++ ) {
        p[ k ] = new T* [ n[1] ];

        for ( int j = 0; j < n[1]; j++ ) {
            p[ k ][ j ] = (T*) &( A3d[ k ][ j ][ 0 ] );
        }
    }
    return p;
}

template <typename T>
void free_ptr3d( T*** p, int*n )
{
    for( int k = 0; k < n[2]; k++ ) { delete[] p[ k ]; }
    delete[] p;
}

修改后的cppfunc.cpp：

#include "getptr3d.hpp"
...
void cppfunc_( int* A, int* n )
{
    int*** A3d = get_ptr3d( A, n );  // can also be used for double***
    ... // use A3d[ k ][ j ][ i ] 
        // or pass A3d to other functions like func( int*** B, ... )
    free_ptr3d( A3d, n ); // when calculation finished
}