calling fortran subroutine from matlab: keep getting errors

Hello everyone,

my fortran subroutine, called "opt_r_finder(b,u,rall,cr,opt_r)", takes four input and outputs one. Here are the four inputs and one output, its type and its size:

• input 1) real*8 b(14536)
• input 2) real*8 u(1817,8)
• input 3) integer rall(40320,8)
• input 4) integer cr(8)
• output 1) opt_r(1817,8)

By modifying the example that mathworks provided, I have written a mex function to call this subroutine from matlab that looks like follows:

        subroutine mexFunction(nlhs,plhs,nrhs,prhs)

        implicit none

        mwPointer plhs(*),prhs(*)
        integer nlhs, nrhs, j

        mwPointer mxGetPr
        mwPointer mxGetM, mxGetN

        mwPointer x_ptr,y_ptr

        mwPointer mrows,ncols
        mwSize size

        integer,parameter::nst=1817,sN=8,numr=40320
        real*8 b(nst*sN),u(nst,sN)
        integer rall(numr,sN), cr(sN), opt_r(nst,sN)
        mwPointer mxCreateNumericArray,mxClassIDFromClassName

        mrows=mxGetM(prhs(1))
        ncols=mxGetN(prhs(1))
        size=mrows*ncols

        x_ptr=mxGetPr(prhs(1))
        call mxCopyPtrToReal8(x_ptr,b,size)

        mrows=mxGetM(prhs(2))
        ncols=mxGetN(prhs(2))
        size=mrows*ncols
!
        x_ptr=mxGetPr(prhs(2))
        call mxCopyPtrToReal8(x_ptr,u,size)
!
        mrows=mxGetM(prhs(3))
        ncols=mxGetN(prhs(3))
        size=mrows*ncols
!
        x_ptr=mxGetPr(prhs(3))
        call mxCopyPtrToInteger4(x_ptr,rall,size)
!
        mrows=mxGetM(prhs(4))
        ncols=mxGetN(prhs(4))
        size=mrows*ncols
!
        x_ptr=mxGetPr(prhs(4))
        call mxCopyPtrToInteger4(x_ptr,cr,size)

        j=mxClassIDFromClassName('int32')
        plhs(1)=mxCreateNumericArray(2,[nst,sN],j,0)
        y_ptr=mxGetPr(plhs(1))

        call opt_r_finder(b,u,rall,cr,opt_r)

        call mxCopyInteger4ToPtr(opt_r,y_ptr,nst*sN)

        return
        end

However, when I call the subroutine, I get "segmentation violation" error. The subroutine works fine in fortran, so the only reason for the error I can think of is my mex function is written incorrectly. Especially, after checking where in the subroutine is causing the error, I suspect that there is something wrong with the input "rall", which is an integer matrix. Can someone please let me know what is wrong with my mex function? Thank you in advance.

Answer 1

[James Tursa]

For Fortran mex routines, one should always use the exact subroutine signatures that are listed in the doc to make sure there is not a mismatch with the API library functions. With C/C++ you can rely on automatic type promotion for some of the arguments, but you don't get that with Fortran. So my first advice is to rewrite your calling syntax to use variables, not explicit values. And then, if anything is amiss in the prhs it will often cause a crash if your mex routine doesn't handle it. So my second advice is to put in checking code to make sure your routine is robust against unexpected inputs. E.g., something like this (CAUTION: UNTESTED)

          subroutine mexFunction(nlhs,plhs,nrhs,prhs)

          implicit none

          mwPointer plhs(*),prhs(*)
          integer nlhs, nrhs

          mwPointer mxGetPr, mxGetData
          mwPointer mxGetM, mxGetN

          mwPointer x_ptr,y_ptr

          mwPointer mrows, ncols
          mwSize size

          integer,parameter::nst=1817,sN=8,numr=40320
          real*8 b(nst*sN),u(nst,sN)
          integer*4 rall(numr,sN), cr(sN), opt_r(nst,sN)
          mwPointer mxCreateNumericArray, mxClassIDFromClassName
          integer*4, external :: mxIsDouble, mxIsComplex, mxIsSparse
          integer*4, external :: mxIsInt32
          mwPointer, external :: mxGetNumberOfElements
          mwSize ndim
          mwSize dims(2)
          integer*4 classid, ComplexFlag

          if( nrhs /= 4 ) then
              call mexErrMsgTxt('Expected exactly 4 inputs')
          endif
          if( nlhs > 1 ) then
              call mexErrMsgTxt('Too many outputs')
          endif

          mrows=mxGetM(prhs(1))
          ncols=mxGetN(prhs(1))
          size=mrows*ncols
          if( mxGetNumberOfElements(prhs(1)) /= nst*sN .OR.
       #      mxIsDouble(prhs(1)) /= 1 .OR.
       #      mxIsComplex(prhs(1)) == 1 .OR.
       #      mxIsSparse(prhs(1)) == 1 ) then
              call mexErrMsgTxt('1st argument is not correct')
          endif

          x_ptr=mxGetPr(prhs(1))
          call mxCopyPtrToReal8(x_ptr,b,size)

          mrows=mxGetM(prhs(2))
          ncols=mxGetN(prhs(2))
          size=mrows*ncols
          if( mrows /= nst .OR. ncols /= sN .OR.
       #      mxIsDouble(prhs(2)) /= 1 .OR.
       #      mxIsComplex(prhs(2)) == 1 .OR.
       #      mxIsSparse(prhs(2)) == 1 ) then
              call mexErrMsgTxt('2nd argument is not correct')
          endif
  !
          x_ptr=mxGetPr(prhs(2))
          call mxCopyPtrToReal8(x_ptr,u,size)
  !
          mrows=mxGetM(prhs(3))
          ncols=mxGetN(prhs(3))
          size=mrows*ncols
          if( mrows /= numr .OR. ncols /= sN .OR.
       #      mxIsInt32(prhs(3)) /= 1 .OR.
       #      mxIsComplex(prhs(3)) == 1 ) then
              call mexErrMsgTxt('3rd argument is not correct')
          endif
  !
          x_ptr=mxGetData(prhs(3))
          call mxCopyPtrToInteger4(x_ptr,rall,size)
  !
          mrows=mxGetM(prhs(4))
          ncols=mxGetN(prhs(4))
          size=mrows*ncols
          if( mxGetNumberOfElements(prhs(4)) /= sN .OR.
       #      mxIsInt32(prhs(4)) /= 1 .OR.
       #      mxIsComplex(prhs(4)) == 1 ) then
              call mexErrMsgTxt('4th argument is not correct')
          endif
  !
          x_ptr=mxGetData(prhs(4))
          call mxCopyPtrToInteger4(x_ptr,cr,size)

          ndim = 2
          dims(1) = nst
          dims(2) = sN
          classid = mxClassIDFromClassName('int32')
          ComplexFlag = 0
          plhs(1)=mxCreateNumericArray(ndim,dims,classid,ComplexFlag)
          y_ptr=mxGetData(plhs(1))

          call opt_r_finder(b,u,rall,cr,opt_r)

          size = nst*sN
          call mxCopyInteger4ToPtr(opt_r,y_ptr,size)

          return
          end