Matlab与Fortran混合编程：m…

原文地址：Matlab与Fortran混合编程：mexfunction 作者：fishpjy

一、编程环境

Matlab R2012b + Visual Studio 2010 + Intel Visual Fortran Composer XE 2011

二、混合编程示例

通过下面的简单例子来简单说明Matlab与Fortran的混合编程。

例：计算一个向量的各个元素的乘积。

1、Fortran源程序：

===============================

      program main

      implicit none

      integer n, i

      parameter(n=10)

      real*8 x(n), y

C     给向量赋值 x(i)=2*i

      do i = 1, n

          x(i) = 2.0*i 

      enddo

      call myprod(x,y,n)

 

      write(*,*) "y=", y 

      stop

      end

C-----------------------------------------------------------------------

      subroutine myprod(x, y, n)

      integer i, n

      real*8 x(*), y

c

      y = 1.0

      do i=1, n

          y = y * x(i) 

      end do

    

      return

      end

==========================================================================

2、改写源程序，使得 可供Matlab调用，x由Matlab指定，然后将结果返回给Matlab。程序如下：

#include "fintrf.h"  

C Gateway routine ! 入口子程序

C 子程序必须取名mexFunction，类似于C语言的main函数

C 在Matlab中用mex编译后会生成一个与文件名同名的.mexw32或.mexw64文件

C 如：mytest.mex64，调用时将其看成是Matlab的函数即可，如

C

C y = mytest(x) % 这里 x 是一个向量

C

C-----------------------------------------------------------------------

      subroutine mexFunction(nlhs, plhs, nrhs, prhs)

C

C 注意：子程序名必须为mexFunction

C nlhs：用于记录输出参数（left hand side）的个数

C plhs：用于记录输出的指针数组，指向各个输出参数的首地址

C nrhs：用于记录输入参数（right hand side）的个数

C prhs：用于记录输入的指针数组，指向各个输入参数的首地址

C

C-----------------------------------------------------------------------

C     Declarations ! 变量声明

      implicit none

C     mexFunction arguments:

      integer nlhs, nrhs

      mwPointer plhs(*), prhs(*)

C     Function declarations:

      mwPointer mxGetPr, mxGetM, mxGetN

      mwPointer mxCreateDoubleMatrix

C mxGetPr：获取数据的实数部分

C mxGetPi：获取数据的虚数部分，这里没有用到

C mxGetM：获取数据的行数

C mxGetN：获取数据的列数

C mxCreateDoubleMatrix：创建一个双精度Matlab矩阵

C

C     Array information:

      mwPointer x_pr, y_pr

      mwPointer m, n

      mwSize size

 

      real*8 x(100), y

C 由于不知道 x 的长度，这里规定最长不超过100

C-----------------------------------------------------------------------

C 获取输入数据的行数与列数，并判断数据个数是否超过100

      m = mxGetM(prhs(1))

      n = mxGetN(prhs(1))

      size = m * n

      if(size .gt. 100) then

          call mexErrMsgIdAndTxt ('MATLAB:InputTooBig',

      +               'Input #1: number of elements exceeds buffer')

      endif

C

C 获取从Matlab中输入的数据

      x_pr = mxGetPr(prhs(1))

C 将数据转化为双精度，并复制给Fortran中的 x，size表示数据个数

      call mxCopyPtrToReal8(x_pr,x,size)

C 调用计算子程序

      call myprod(x,y,size)

C 将计算结果传给 Matlab

C 先创建一个双精度实矩阵，这里是1x1的

      plhs(1) = mxCreateDoubleMatrix(1,1,0) !  如果第三个参数是1，则表示带虚部

      y_pr = mxGetPr(plhs(1))

C

C 将y的值（实数部分）复制给 y_pr，也就是传给 Matlab

      call mxCopyReal8ToPtr(y,y_pr,1)      

      return

      end

C-----------------------------------------------------------------------

C 计算子程序，不用做任何改动

      subroutine myprod(x, y, n)

      integer i, n

      real*8 x(*), y

C

      y = 1.0

      do i=1, n

          y = y * x(i) 

      end do      

      return

      end

3、 编译：将上面的源程序保存，取名为 mytest.F， 打开Matlab，输入命令: 

mex mytest.F

回车后Matlab会编译生成 mytest.mexw32 或 mytest.mexw64

这里假定Matlab中的Fortran编译器已经设置好了。设置方法：输入 mex -setup，

根据提示选择相应的Fortran编译器即可，比如

Intel Visual Fortran 12.0 (with Microsoft Visual C++ 2010 linker) 

4、测试：在Matlab命令窗口输入下面的命令：

x=rand(10,1);

y=mytest(x)