调用 C_F_POINTER 时 Fortran 指针和 Fortran 可分配的区别

Question

问题是，'C_F_POINTER' 成功编译（ifort 版本 19.0.5.281），使用 'allocatable arrays' 作为参数，它的工作方式与使用 'pointer' 作为参数的情况完全相同。

program test1
    use mkl_spblas
    use omp_lib
    use iso_c_binding

    implicit none
    integer, parameter    :: DIM_ = 4, DIM_2 = 6
    integer               :: stat, i
    integer               :: irn(DIM_2), jcn(DIM_2)
    real*8                :: val(DIM_2)
    integer(c_int)        :: indexing
    integer               :: DIM_r, DIM_c
    type(c_ptr)           :: rows_start_c, rows_end_c, col_indx_c, values_c
(*1)!integer,allocatable   :: rows_start_f(:), rows_end_f(:), col_indx_f(:)
    !real*8 ,allocatable   :: values_f(:)
(*2)integer  ,pointer     :: rows_start_f(:), rows_end_f(:), col_indx_f(:)
    real*8   ,pointer     :: values_f(:)
    type(SPARSE_MATRIX_T) :: mat1, mat2

    irn = (/ 2, 2, 3, 4, 0, 0 /)
    jcn = (/ 1, 2, 3, 2, 0, 0 /)
    val = (/ 5, 8, 3, 6, 0, 0 /)

    call omp_set_num_threads(1)

    stat = mkl_sparse_d_create_coo (A=mat1, indexing=SPARSE_INDEX_BASE_ONE, &
                                    rows=DIM_, cols=DIM_, nnz=DIM_,&
                                    row_indx=irn, col_indx=jcn, values=val  )
    if (stat /= 0) stop 'Error in mkl_sparse_d_create_coo'

    stat = mkl_sparse_convert_csr (source=mat1,&
                                   operation=SPARSE_OPERATION_NON_TRANSPOSE, &
                                   dest = mat2 )
    if (stat /= 0) stop 'Error in mkl_sparse_convert_csr'

    stat = mkl_sparse_d_export_csr(mat2, indexing, DIM_r, DIM_c,  &
                                   rows_start_c, rows_end_c, col_indx_c, values_c)

(*3)call c_f_pointer(rows_start_c, rows_start_f, [DIM_r])
    call c_f_pointer(rows_end_c  , rows_end_f  , [DIM_c])
    call c_f_pointer(col_indx_c  , col_indx_f  , [rows_end_f(DIM_r)-1])
    call c_f_pointer(values_c    , values_f    , [rows_end_f(DIM_r)-1])

    stat = mkl_sparse_destroy (A=mat1)
    if (stat /= 0) stop 'Error in mkl_sparse_destroy (mat1)'

    stat = mkl_sparse_destroy (A=mat2)
    if (stat /= 0) stop 'Error in mkl_sparse_destroy (mat2)'

    call mkl_free_buffers

(*4)print *, 'rows_start'
    print *, rows_start_f
    print *, 'rows_end'
    print *, rows_end_f
    print *, 'col_indx'
    print *, col_indx_f
    print *, 'values'
    print *, values_f
    print *, 'indexing'
    print *, indexing
    print *, 'size(values_f,1)'
    print *, size(values_f,1)

end program test1

在上面的测试代码中，我在代码的左侧将一些点标记为（*1）、（*2）等。

(*1) & (*2) : 可分配数组版本和指针版本的代码 (*3) : 我称之为 'C_F_POINTER' (*4) : 打印语句以查看输出

在 (*1) 和 (*2) 情况下，结果“完全相同”，并且所有值都正确转换为所需的 CSR 格式。

 rows_start
           1           1           3           4
 rows_end
           1           3           4           5
 col_indx
           1           2           3           2
 values
   5.00000000000000        8.00000000000000        3.00000000000000     
   6.00000000000000     
 indexing
           1
 size(values_f,1)
           4

两年前我在 StackOverflow 中发现了一个类似的问题 (difference between fortran pointers or allocatable arrays for c_f_pointer call)。

这个问题和我现在想的完全一样。

如果我用我的话重新排列问题，

1. 指针和可分配数组的区别？

• 据我所知，在 C 语言中，数组存储在连续内存中，可以用指向其第一个元素的指针表示。在 Fortran90 中，如果我将数组作为“假定大小的数组”传递给子程序，代码的行为就像它从不关心它的分配方式、大小如何，并将数组视为存储在连续站点中的一维数组。
• 在下面的代码中，子程序 'assign_A' 只是将 'tot_array(1,2)' 作为其起点，并在连续的站点上执行它的工作，并且似乎甚至超出了 'tot_array' 的范围！ （tot_array 是 2x2 矩阵，assign_A 的 do 循环从 tot_array(1,2) 开始运行 5 次）我“感觉”指针和可分配数组在这个意义上是类似的东西。但显然，正如difference between fortran pointers or allocatable arrays for c_f_pointer call 中的答案，它们是不同的东西。为什么数组作为“假定大小”传递给子程序时就像指针一样？

program assumed_size_array_test
  implicit none
  external assign_A
  real*8 :: tot_array(2,2)
  integer:: i
  
  ! Initially 'tot_array' set to be 1.d0
  tot_array = 1.d0
  
  write(*,*) 'Before'
  write(*,'(5f5.2)') tot_array
  
  call assign_A(tot_array(1,2))
  
  write(*,*) 'After'
  write(*,'(5f5.2)') tot_array

end program

subroutine assign_A(A)
  implicit none
  real*8, intent(inout) :: A(*)
  integer :: i
  
  do i = 1,5
    A(i) = 2.d0
  enddo

end subroutine

 Before
 1.00 1.00 1.00 1.00
 After
 1.00 1.00 2.00 2.00

1. 在 Fortran90 中调用“C_F_POINTER”时使用“可分配数组”和“指针”有什么区别吗？
   • 我使用了 ifort 版本 19.0.5.281，据我检查，这个编译器似乎给了我完全相同的结果。如果可以的话，我更喜欢使用可分配数组而不是指针。将“可分配数组”和“指针”与“C_F_POINTER”一起使用有什么区别吗？这样做有什么需要注意的吗？
   • difference between fortran pointers or allocatable arrays for c_f_pointer call 中的答案说我应该使用指针，而不是使用带有 C_F_POINTER 的可分配数组，但它似乎是一些持续存在的问题，当时尚未完全确定。为什么'C_F_POINTER'是为fortran指针设计的，可以很好地与可分配数组一起工作并且结果是一样的，有什么结论吗？

感谢您阅读此问题。

Answer 1

显然，Fortran POINTER 变量和ALLOCATABLE 变量在它们的内部实现上有很多共同点。其中大部分都在后台，不应直接访问。两者都分配一些内存，并且可能使用相同的操作系统或 C 运行时库的分配器。例如，malloc()。

在两者中都有一些内存被分配或指向并由简单地址（对于标量）或数组描述符（对于数组）描述。

指针和可分配变量的主要区别在于您可以使用它们做什么以及编译器将为您做什么。您可以将可分配对象视为一种“智能指针”，与 C++ 中的 std::unique_ptr 非常相似。回想一下 C++ 中发生的事情，你有 new 和 delete，它们依次调用 malloc 和 free，但你不能混合使用它们。当然也不允许你手动修改存储在 C++ 智能指针中的地址。

当您将可分配变量发送到需要指针的过程时，任何事情都可能发生，这是一种未定义的行为。但是，如果内部隐藏结构具有类似的布局，则可能会发生您实际上将可分配内部设置为指向某些未通过可分配分配的内存。然后，您可能会认为一切正常，并且您有一个新功能。但是，当释放的时间到来时，可分配对象通常会自动释放，它很容易以非常不可预测的方式失败。它可能会在代码的非常奇怪的地方崩溃，结果可能是错误的等等。任何事情都有可能发生。

例如，这个非常丑陋的程序也适用于我（在 gfortran 中）：

subroutine point(ptr, x)
  pointer :: ptr
  target :: x
  ptr => x
end subroutine

  interface
    subroutine point(ptr, x)
      allocatable :: ptr
      target :: x
    end subroutine    
  end interface

  allocatable z

  y = 1.0

  call point(z, y)

  print *, z
end

但是你应该永远不要做这样的事情。这真的是非常非常错误的事情。如果您将z 设为局部变量，使其被释放，或者如果您尝试释放它，它将崩溃。那是因为编译器拥有的唯一信息是地址。在内部，可分配对象实际上看起来与指针相同。它只是一个地址（对于标量）。唯一的区别是您可以使用它做什么以及编译器会自动为您做什么。

这甚至不会崩溃，因为我提到了内部实现的相似之处。但同样是错误的。

subroutine point(ptr, x)
  pointer :: ptr
  target :: x
  ptr => x
end subroutine     

  interface
    subroutine point(ptr, x)
      allocatable :: ptr
      target :: x
    end subroutine    
  end interface

  allocatable z
  pointer y

  allocate(y)
  y = 1.0

  call point(z, y)

  print *, z

  deallocate(z)
end

它之所以存在，是因为可分配和指针在 gfortran 中都使用相同的内部分配器 (malloc)，并且它们都实现为简单地址。