Search code examples
matlabfortranmexmatlab-deployment

Portable declaration of REAL variables in mex gateway for Fortran


I am writing a mex gateway for a piece of Fortran code.

In the Fortran code, for portability, the floating-point variables are declared as

REAL(kind(0.0D0)) :: x, y, etc

(BTW, I am aware that there are better ways to do it, as discussed at Fortran: integer*4 vs integer(4) vs integer(kind=4), What does "real*8" mean?, and https://software.intel.com/en-us/blogs/2017/03/27/doctor-fortran-in-it-takes-all-kinds )

However, it seems to me that mex supports only REAL*8 and REAL*4, the former being Double, the latter being Single. I got this impression from the following functions/subroutines:

mxIsDouble, mxIsSingle, mxCopyPtrToReal8, mxCopyReal8ToPtr, mxCopyPtrToReal4, mxCopyReal4ToPtr

My questions are as follows.

  1. Is it true that mex supports only REAL*8 and REAL*4?

  2. Does it improve the portability of the mex gateway if I declare double-precision floating-point variables as

    REAL(kind(0.0D0)) :: x, y, etc

or even

integer, parameter :: dp = selected_real_kind(15, 307)
real(kind=dp) :: x, y, etc

Or should I simply declare

REAL*8 :: x, y, etc 
  1. Are REAL*8 and/or REAL*4 supported on all platforms? If no, does this mean that MATLAB mex is intrinsically unportable?

  2. What is the best way to specify the kind of floating-point variables in mex gateways for Fortran code?

The following code is an example. See the declaration of x, y, and xs.

#include "fintrf.h"

      subroutine mexFunction(nlhs, plhs, nrhs, prhs)
C     y = square (x) 
C     x: a floating point scalar
C     y: x^2

      implicit none

C     mexFunction arguments 
      integer, intent(in) :: nlhs, nrhs 
      mwPointer, intent(in) :: prhs(nrhs)
      mwPointer, intent(inout) :: plhs(nlhs)

C     function declarations:
      mwPointer, external :: mxCreateDoubleScalar, mxGetPr
      mwSize, external :: mxGetM, mxGetN
      integer*4, external :: mxIsDouble, mxIsSingle 

C     variables
      mwSize, parameter :: mwOne = 1
      integer, parameter :: dKind = kind(0.0D0) 
      integer, parameter :: sKind = kind(0.0) 
      real(kind=dKind) :: x, y ! Does this improve the portablity?
      real(kind=sKind) :: xs ! Does this improve the portablity?

C     validate number of arguments
      if (nrhs .ne. 1) then
         call mexErrMsgIdAndTxt ('mex:nInput', '1 input required.')
      endif
      if (nlhs .gt. 1) then
         call mexErrMsgIdAndTxt ('mex:nOutput', 'At most 1 output.')
      endif

C     validate input
      if (mxIsDouble(prhs(1)) .ne. 1 .and. mxIsSingle(prhs(1)) .ne. 1) 
      ! What if the input is a floating point number but neither Double nor Single? 
     + then 
          call mexErrMsgIdAndTxt ('mex:Input', 'Input a real number.')
      endif
      if (mxGetM(prhs(1)) .ne. 1 .or. mxGetN(prhs(1)) .ne. 1) then
          call mexErrMsgIdAndTxt ('mex:Input', 'Input a scalar.')
      endif

C     read input
      if (mxIsDouble(prhs(1)) .eq. 1) then
         call mxCopyPtrToReal8(mxGetPr(prhs(1)), x, mwOne) 
      else
         call mxCopyPtrToReal4(mxGetPr(prhs(1)), xs, mwOne)
         x = real(xs, dKind) 
         ! What if the input is a floating point number but neither REAL*8 nor REAL*4
      endif

C     do the calculation
      y = x**2

C     write output
      plhs(1) = mxCreateDoubleScalar(y)

      return
      end subroutine mexFunction

The code runs correctly. Yet I am not sure whether it is portable.


Solution

  • REAL*4 and REAL*8 are non-standard and non-portable. REAL(KIND(0.0D0) gets you DOUBLE PRECISION on every platform, as this is required by the Fortran standard.

    I can't speak to MEX gateways, but you should avoid obvious non-standard features.

    A popular choice is to define a module that declares named (PARAMETER) constants for the kinds in use. For example:

    module kinds integer, parameter :: SP = KIND(0.0) integer, parameter :: DP = KIND(0.0D0) end module kinds

    Then you can use SP and DP as kind values. If you ever need to change these, just edit the module.