如何在 fortran mex 文件中使用 mxCalloc

Question

我必须在 mex 文件中使用 mxCalloc 函数而不是常规分配，以避免 matlab 在使用 dgesv 时崩溃。 我尝试了很多方法，但都没有奏效。 这是其中一个示例

 #include "fintrf.h"

C     Gateway subroutine
      subroutine mexfunction(nlhs, plhs, nrhs, prhs)

C     Declarations
      implicit none

C     mexFunction arguments:
      mwPointer plhs(*), prhs(*)
      integer nlhs, nrhs

C     Function declarations:
      mwPointer mxGetPr
      mwPointer mxCreateDoubleMatrix
      mwPointer mxGetM

C     Pointers to input/output mxArrays:
      mwPointer pr_A, pr_B

C     Array information:

      mwPointer sizea,mxCalloc
      real*8 :: A,B
      character*120 :: line

C     Get the size of the input array.
      sizea = mxGetM(prhs(1))

      A=mxCalloc(sizea*sizea,8)
      B=mxCalloc(sizea*sizea,8)

C     Create Fortran array from the input argument.
      pr_A = mxGetPr(prhs(1))

      call mxCopyPtrToReal8(pr_A,A,sizea**2)

C     Create matrix for the return argument.
      plhs(1) = mxCreateDoubleMatrix(sizea, sizea, 0)
      pr_B = mxGetPr(plhs(1))

      write(line,*), sizea
      call mexPrintf(line) 
      B=A

      call mxCopyReal8ToPtr(B,pr_B,sizea*sizea)

      return
      end

当我运行这段代码时，我得到以下结果

A = [0.9575 , 0.1576 ; 0.9649 , 0.9706]

测试(A) = [0.9575 , 0 ; 0.9649 , 0]

但是如果我换行 call mxCopyReal8ToPtr(B,pr_B,sizea*sizea) 到

call mxCopyReal8ToPtr(A,pr_B,sizea*sizea), 结果是正确的

变量 sizea 等于 2，这是正确的，但我无法访问 A 的任何成员，比如 A(1,1)，这是我遇到的错误：

错误 #6410：此名称尚未声明为数组或 功能。 [一]

Answer 1

（免责声明：我没有使用 Mex 的经验，所以请使用以下代码...）

这是 @ftiaronsem 和 Dave 之前的 one 答案的一些延续。通过查看这个tutorial，似乎可以在mexfunction() 中使用Fortran 可分配数组。但是，OP 的目的是使用mxCalloc() 而不是可分配的数组进行内存管理。

根据Matlab的在线手册，mxCalloc()似乎返回分配内存的地址为mwpointer（这只是integer*8在64位机器上的别名，根据头文件@987654323 @)。所以，我们首先得到地址为

mwpointer iA
iA = mxCalloc( sizea * sizea, 8 )

（B 也是如此）。接下来，我们要以二维 Fortran 数组的形式访问从 iA 开始的内存。这可以通过使用c_f_pointer()（在 F2003 中）来完成。但是，由于这个函数的第一个参数接收到type(c_ptr)，我们继续如下：

use iso_c_binding
type(c_ptr) cA
real*8, pointer :: A(:,:)

cA = transfer( iA, cA )                      !! cast integer*8 to c_ptr
call c_f_pointer( cA, A, [ sizea, sizea ] )  !! init an array pointer with c_ptr

在此语句之后，我们可以将A 用作普通的二维 Fortran 数组（例如，我们可以将其传递给 LAPACK 例程）。将这些修改包含在 OP 的代码中会给出以下内容。那你能不能试试看它是否有效...？

---

更新代码（ver1）：

#include "fintrf.h"

C     Gateway subroutine
      subroutine mexfunction(nlhs, plhs, nrhs, prhs)

      use iso_c_binding    !<---

C     Declarations
      implicit none

C     mexFunction arguments:
      mwPointer plhs(*), prhs(*)
      integer nlhs, nrhs

C     Function declarations:
      mwPointer mxGetPr
      mwPointer mxCreateDoubleMatrix
      mwsize    mxGetM, sizea          !<--- changed to mwsize (may not be necessary, though)

C     Pointers to input/output mxArrays:
      mwPointer pr_A, pr_B

C     Array information:

      mwPointer mxCalloc

      mwPointer       :: iA, iB           !<---
      type(c_ptr)     :: cA, cB           !<---
      real*8, pointer :: A(:,:), B(:,:)   !<---

      character*120 :: line

C     Get the size of the input array.
      sizea = mxGetM(prhs(1))

      iA = mxCalloc( sizea**2, 8 )   !<---
      iB = mxCalloc( sizea**2, 8 )   !<---

      cA = transfer( iA, cA )        !<---
      cB = transfer( iB, cB )        !<---

      call c_f_pointer( cA, A, [ nsizea, nsizea ] )   !<---
      call c_f_pointer( cB, B, [ nsizea, nsizea ] )   !<---

C     Create Fortran array from the input argument.
      pr_A = mxGetPr(prhs(1))

      call mxCopyPtrToReal8( pr_A, A, sizea**2 )

C     Create matrix for the return argument.
      plhs(1) = mxCreateDoubleMatrix(sizea, sizea, 0)
      pr_B = mxGetPr(plhs(1))

      write(line,*), sizea
      call mexPrintf(line) 
      B = A

      call mxCopyReal8ToPtr( B, pr_B, sizea**2 )

      !! may need to call mxFree( iA ) and mxFree( iB ) manually?
      !! (or Matlab may automatically do it upon exit)

      return
      end

Answer 2

代码中的 A 被声明为单个变量，但 size(A,1) 试图像访问数组一样访问它。您必须明确告诉 Fortran A 具有数组的形状。快速浏览其余代码后，我假设 A 和 B 应该是二维的

 real*8 :: A(:,:),B(:,:)

现在我不太了解 mxCalloc 的工作原理以及它如何处理内存分配和与 Fortran 的接口，但可能您还需要将 A 和 B 声明为指针。

real*8, pointer :: A(:,:),B(:,:)