How to compile a Fortran subroutine by using the MEX engine

Hi!

I have been trying to use a Matlab function (interp2) in Fortran using the Mex engine.

I have adapted the Fengdemo.F to suit my needs and the test cases work.

However, it is supposed to be part of a larger Fortran project with input coming from the rest of the program.

Therefore, I need the function as a subroutine or such, but mex does not let me compile it as anything other than a main() program, which inherentely does not accept input.

I have tried defining it as a subroutine, as well as a module.

Any ideas on what I might do wrong or how to get around this problem?

I find it hard to imagine the mex engine is only able to run main() programs as this would render this function more or less useless.

I am attaching a minimal example based on Fengdemo.F which is now a subroutine and just supposed to square an input.

Any hints are highly appreciated!

Kind regards

Lars

module InterpMatlab_mod

contains

#include "fintrf.h"
! Interpolation
#if 0
!
!     fengdemo.F
!     .F file need to be preprocessed to generate .for equivalent
!
#endif
!
!     fengdemo.f
!
!     This is a simple program that illustrates how to call the MATLAB
!     Engine functions from a FORTRAN program.
!
! Copyright 1984-2018 The MathWorks, Inc.
!======================================================================
! 

      subroutine InterpolMatlab(x) 

!     Declarations
      implicit none

      mwPointer engOpen, engGetVariable, mxCreateDoubleMatrix
#if MX_HAS_INTERLEAVED_COMPLEX
      mwPointer mxGetDoubles
#else
      mwPointer mxGetPr
#endif
      mwPointer x_point, xOut_point 
      integer engPutVariable, engEvalString, engClose
      


      
      mwSize M
      parameter(M=1) 

      integer, PARAMETER :: sp = KIND(1.0D0)
      !----------------------------------------------------------------
      ! Input
      integer                                   :: x

      ! Output

      integer                                   :: x_square



    !--------------------------------------------------------------------------------------------------!
      ep = engOpen('matlab ')

      if (ep .eq. 0) then
         write(6,*) 'Can''t start MATLAB engine'
         stop
      endif
      write(6,*) 'started MATLAB engine'

!     Place the variables into the MATLAB workspace
      x_point = mxCreateDoubleMatrix(M, M, 0) 
#if MX_HAS_INTERLEAVED_COMPLEX
      call mxCopyReal8ToPtr(x, mxGetDoubles(x_point), M)
#else
      call mxCopyReal8ToPtr(x, mxGetPr(x_point),M)
#endif
      status = engPutVariable(ep, 'x_point', x_point)
      if (status .ne. 0) then 
         write(6,*) 'engPutVariable failed'
         stop
      endif
      
! Calculations in Matlab
      if (engEvalString(ep, 'x_out=x^2') .ne. 0) then
        write(6,*) 'engEvalString failed'
        stop
      endif


      xOut_point = engGetVariable(ep, 'x_out')
      write(6,*) 'Got pointer'

#if MX_HAS_INTERLEAVED_COMPLEX
      call mxCopyPtrToReal8(mxGetDoubles(xOut_point), x_square, 1)
#else
      call mxCopyPtrToReal8(mxGetPr(xOut_point), x_square, 1)
#endif
      write(6,*) 'Result=', x_square
      
     
      print *, 'Type 0 <return> to Exit'
      print *, 'Type 1 <return> to continue'

      read(*,*) temp

      if (temp.eq.0) then
         print *, 'EXIT!'
         status = engClose(ep)

         if (status .ne. 0) then 
            write(6,*) 'engClose failed'
         endif

         stop
      end if
      

       call mxDestroyArray(x_point)
       call mxDestroyArray(xOut_point)

      status = engClose(ep)
      
      if (status .ne. 0) then 
         write(6,*) 'engClose failed'
         stop
      endif

      stop


      end
end module InterpMatlab_mod

[James Tursa]

We need you to post a small sample of code that reproduces the problem. Then we can steer you in the right direction. Do you want to call this Fortran code from MATLAB, i.e. getting inputs from MATLAB and perhaps returning outputs to MATLAB? That is the general outline of a mex function using the subroutine mexFunction( ) interface for the input/output.

[Lars Schmitz]

Hi James,

Thanks for your answer.

No, actually I want it the other way around. I am trying to call a Matlab function (interp2) from Fortran.

I have edited my original question and added a minimal example.

Kind regards

Lars

[Benjamin Thompson]

Probably not worth the effort. Can you find or code up a 2D interpolation algorithm in Fortran directly? Otherwise you can probably call C functions from Fortran. See "Numerical Recipies in C" for possible sample implementations, or I suppose you can use MATLAB Coder to generate the C code, and compile it into a C library that Fortran could use. Seems like a lot of work for a fairly simple algorithm though.

[Lars Schmitz]

Probably right. Was hoping for a bit less effort on this and then I was already in over my head.

Answer 1

[James Tursa]

It looks like you are taking input, calling a MATLAB function, and getting output from that function. I still say it would be easier to use a mex routine for this. However, if you insist on using an Engine application that can be done too. My comments on your posted Engine code:

1)

#if MX_HAS_INTERLEAVED_COMPLEX
    mwPointer mxGetDoubles
#else
    mwPointer mxGetPr
#endif
    :
#if MX_HAS_INTERLEAVED_COMPLEX
      call mxCopyPtrToReal8(mxGetDoubles(xOut_point), x_square, 1)
#else
      call mxCopyPtrToReal8(mxGetPr(xOut_point), x_square, 1)
#endif

I know that MATLAB examples have this type of code now because of the change in how complex numbers are stored. But you can make your code simpler by just using mxGetData( ) instead:

mwPointer, external :: mxGetData
:
call mxCopyPtrToReal8(mxGetData(xOut_point), x_square, 1) ! But don't use literal 1 here, see below

2)

You don't define ep. WIth your implicit none this shouldn't have even compiled. You need something like this:

mwPointer ep

3)

You use an explicit integer argument in an API function call. NEVER do this in Fortran! You don't get automatic type conversion in Fortran like you do in C/C++. Always use typed variables for this to ensure the integer size is the correct size for the interface. E.g., this line

x_point = mxCreateDoubleMatrix(M, M, 0) 

should be something like this instead

integer*4 :: ComplexFlag = 0
:
x_point = mxCreateDoubleMatrix(M, M, ComplexFlag) 

Same comment applies to your mxCopyPtrToReal8(...,1) call. That last input argument needs to be an explicitly typed variable holding the value 1, not a literal integer 1.

4)

You don't use the correct names on the Engine side. E.g.,

status = engPutVariable(ep, 'x_point', x_point) ! Here you put a variable called x_point in the Engine
:
if (engEvalString(ep, 'x_out=x^2') .ne. 0) then ! Here you use a variable called x in the Engine

In the above code, there is no variable called 'x' in the Engine because you haven't put any 'x' there. You put 'x_point'. Either put 'x' in the Engine or change your string to use 'x_point' on the rhs. E.g.,

status = engPutVariable(ep, 'x', x_point) ! Here you put a variable called x in the Engine

5)

You don't use the correct sizes on the data copy, although it does work for scalars. E.g., x_point is an M x M mxArray, but you only copy M elements when you should be copying M*M elements.

6)

You are using floating point copy routines to copy integer variables. Not only could this lead to a crash if the integers are not the correct size, but the bit patterns for the same value are completely different between integer and double precision so the calculations on the MATLAB Engine side will be erroneous. E.g.,

    integer :: x            ! could be 4-byte or 8-byte integer depending on compiler settings
    integer :: x_square     ! could be 4-byte or 8-byte integer depending on compiler settings

But the following code creates a double mxArray and copies the bit patterns of the integers into the double memory locations. Interpreting integer bit patterns as double precision IEEE variables on the MATLAB side is just going to give you garbage.

    x_point = mxCreateDoubleMatrix(M, M, 0) 
    #if MX_HAS_INTERLEAVED_COMPLEX
        call mxCopyReal8ToPtr(x, mxGetDoubles(x_point), M) ! This does a bit pattern memory copy
    #else
        call mxCopyReal8ToPtr(x, mxGetPr(x_point),M)       ! This does a bit pattern memory copy
    #endif

You need x and x_square to be double precision. Using real*8 here simply to match the doc exactly:

    real*8 :: x
    real*8 :: x_square
    

If you really want to use integers, then you should not be creating a double mxArray with mxCreateDoubleMatrix( ). You would need to create an integer mxArray of the correct size with mxCreateNumericMatrix( ) and then use the correct mxCopyEtc( ) routine for that particular integer size.

[Lars Schmitz]

Hi James,

thank you for your input.

Maybe I am mistaken here, but I thought, the Mex routine could only be used the other way. As in compiling fortran code to be called from Matlab. If that is not the case and it is indeed possible with a mex function to call Matlab from Fortran, I am of course open to that.

Thanks for your comments. Many of that originated from me hurriedly trying to produce the minimal example.

However, as stated above, my real code is working like it is supposed to do. My problem is, that it will only compile as a main function, therefore not taking any input. What I need is a subroutine or module or function to call it with input.

Regards

Lars

[James Tursa]

It is true that MATLAB must first call the mex routine. But the mex routine can call back into MATLAB to perform functions via the mexCallMATLAB( ) routine. A mex function doesn't have to take any input from MATLAB. In fact, the mexFunction( ) routine could have only one executable line in it ... a call to your Fortran main program routine.

Where is the input for your Fortran routine coming from?

[Lars Schmitz]

That might be a solution then.

The code is supposed to be part of a large Fortran project. Essentially Fortran is reading in some data from few text files, then performs a large amount of numerical calculations and is supposed to call the Matlab function at the end of each timestep.

How would I get the input from my Fortran main into the mex to call back to Matlab via mexCallMATLAB?

I have used mex functions once or twice but am unfortunately not very familiar with it.

[James Tursa]

Almost all of the Fortran features work the same in an Engine Application vs Mex Routine. In both cases you can read and write data from/to a text file using READ and WRITE statements written exactly the same way. I.E., file I/O code is exactly the same. All of the data manipulation in your Fortran code is the same also. Pretty much all of your existing Fortran code will run the same way, including the text file I/O using normal READ and WRITE statements.

The main differences between the two approaches are console I/O and how you interact with MATLAB. At a high level:

Console I/O:

-----------------

Engine Application:

• CAN use READ(6,etc) and WRITE(6,etc) for console I/O (i.e., user keyboard and display screen)

Mex Routine:

• CAN NOT use READ(6,etc) and WRITE(6,etc) for console I/O (i.e., user keyboard and display screen)
• MUST use call backs into MATLAB via mexCallMATLAB( ) or mexPrintf( ) to get console I/O.

Calling MATLAB Functions:

------------------------------------

Engine Application:

• Uses engPutVariable( ) to put variables into the MATLAB Engine workspace.
• Uses engEvalString( ) to evaluate expressions in MATLAB. The function runs in the MATLAB Engine workspace.
• Uses engGetVariable( ) to get results from MATLAB Engine workspace back into Fortran.

Mex Routine:

• Uses mexCallMATLAB( ) to call a MATLAB function. No need to put variables in workspace first and get result from workspace. The function essentially runs in a temporary workspace assigned to the mex routine and result is available directly in the mex routine.

Converting variables from Fortran to MATLAB mxArray:

---------------------------------------------------------------------------------

Both Engine Application and Mex Routine use same methods:

• Create proper mxArray using one of the mxCreateEtc functions.
• Get the pointer to the MATLAB mxArray variable data via mxGetData( ).
• Copy the Fortran variable data into the mxArray data location using an mxCopyEtc function.

Converting MATLAB mxArray variables to Fortran:

------------------------------------------------------------------

Both Engine Application and Mex Routine use same methods:

• Get the pointer to the MATLAB mxArray variable data via mxGetData( ).
• Copy the MATLAB mxArray variable data into the Fortran variable using an mxCopyEtc function.

I should note at this point that it is sometimes possible to avoid the data copies noted above in mex routines, but it is a little tricky. The data copies cannot be avoided in Engine applications. This can be important if you are working with very large variables. Basically you get the address of the mxArray data memory via the mxGetData( ) function. Then you pass that address by value (e.g., using the %VAL( ) construct) into a Fortran routine with implicit interface expecting the variable to be passed by reference. The Fortran routine gets fooled into thinking it is a normal Fortran variable and will use it as such, all without making a deep copy of the data.

For your application, since you are not running your Fortran code as a MATLAB function, it probably doesn't matter much which method you use (Engine App vs Mex Routine) as long as the variables are not too large. If you already have a large Fortran project then maybe just stick with the Engine App approach for simplicity.

[Lars Schmitz]

Thank you for all those remarks! I believe they will come in really handy!