when using f2py, function scope within fortran module different than when compiled for fortran program?

Question

My problem is that when compiling with f2py, some module variables are not recognized by functions defined within the module. The errors are raised where variable types of arguments passed to the function are declared (such as a variable describing the type of real, or dimension elements). I do not get this error when compiling with gfortran. What is the difference and how can I remedy these errors when compiling with f2py?

My example file, moddata.f90, contains the following code:

module mod

  implicit none

  integer, parameter :: dp = selected_real_kind(15)
  integer, parameter :: nelem = 3
  real(kind=dp), dimension(nelem) :: b

  parameter (b=(/3,1,2/))

contains
  function foo(x,y) result(z)
!  dp, nelem are defined in module above
    real(kind=dp), intent(in) :: x !scalar
    integer, dimension(nelem), intent(in) :: y
!  local variable
    real(kind=dp) :: z

    z = sum(b*y*x)

  end function foo
end module mod

and I compile with

f2py -c -m moddata moddata.f90

and I get these errors:

  y_Dims[0]=(nelem);
             ^
1 warning and 1 error generated.reduce to a constant expression

if I redefine integer, parameter :: nelem=3 before integer, dimension(nelem), intent(in) :: y and recompile, I get

      real(kind=dp) foof2pywrap
                1
Error: Parameter 'dp' at (1) has not been declared or is a variable, which does not reduce to a constant expression

with same errors for each of the real(kind=dp) declarations, and

      foof2pywrap = foo(x, y)
                    1
Warning: Possible change of value in conversion from REAL(8) to REAL(4) at (1)

So I have to redefine dp by integer, parameter :: dp = selected_real_kind(15) in the function. And then it works.

I don't get these errors when I compile this module with a fortran wrapper. I wonder why nelem and dp in the function is not properly scoped with f2py?

Answer 1

Perhaps I am wrong, but I don't think f2py can handle the module+contains feature of Fortran 90. If you turn your code into

function foo(x,y) result(z)
   integer, parameter :: dp = selected_real_kind(15)
   real(kind=dp), intent(in) :: x
   integer, parameter :: nelem = 3
   integer, dimension(3), parameter :  = (/3, 1, 2/)
   integer, dimension(nelem), intent(in) :: y
   real(kind=dp) :: z

   z = sum(b*y*x)
end function

and compile it as before, it works:

 >>> x = 1.0000000000
 >>> y = [2, 3, 4]
 >>> moddata.foo(x,y)
 17.0

---

EDIT

The answer to this question on SO says that f2py doesn't understand how to turn a Fortran function into a python function. So I changed the function foo to subroutine foo2 and then compiled as f2py moddata.f90 -m moddata and got as output

Reading fortran codes...
        Reading file 'moddata.f90' (format:free)
Post-processing...
        Block: moddata
                        Block: moddata
In: :moddata:moddata.f90:moddata
get_parameters: got "invalid syntax (<string>, line 1)" on '(/3, 1, 2/)'
                                Block: foo2
Post-processing (stage 2)...
        Block: moddata
                Block: unknown_interface
                    Block: moddata
                            Block: foo2
Building modules...
        Building module "moddata"...
                Constructing F90 module support for "moddata"...
                    Variables: nelem b dp
                    Constructing wrapper function "moddata.foo2"...
getctype: "real(kind=dp)" is mapped to C "float" (to override define dict(real = dict(dp="<C typespec>")) in /home/jdwood/Documents/Physics/Fortran/tests/f2py/.f2py_f2cmap file).
getctype: "real(kind=dp)" is mapped to C "float" (to override define dict(real = dict(dp="<C typespec>")) in /home/jdwood/Documents/Physics/Fortran/tests/f2py/.f2py_f2cmap file).
getctype: "real(kind=dp)" is mapped to C "float" (to override define dict(real = dict(dp="<C typespec>")) in /home/jdwood/Documents/Physics/Fortran/tests/f2py/.f2py_f2cmap file).
getctype: "real(kind=dp)" is mapped to C "float" (to override define dict(real = dict(dp="<C typespec>")) in /home/jdwood/Documents/Physics/Fortran/tests/f2py/.f2py_f2cmap file).
getctype: "real(kind=dp)" is mapped to C "float" (to override define dict(real = dict(dp="<C typespec>")) in /home/jdwood/Documents/Physics/Fortran/tests/f2py/.f2py_f2cmap file).
getctype: "real(kind=dp)" is mapped to C "float" (to override define dict(real = dict(dp="<C typespec>")) in /home/jdwood/Documents/Physics/Fortran/tests/f2py/.f2py_f2cmap file).
                          z = foo2(x,y)
Wrote C/API module "moddata" to file "./moddatamodule.c"
Fortran 90 wrappers are saved to "./moddata-f2pywrappers2.f90"

So it does look like the double-precision was being lost, so following the suggestion to edit a file called .f2py_f2cmap, I did so and got no errors when it came to dp. However, it still gives an error for nelem, so I can think of two solutions to this:

1. Stick with using 3 in place of nelem
2. Pass nelem and b as variables to the subroutine

I also found, when using the parameter(b = (/3.d0, 1.d0, 2.d0/) ) line, that I was getting a warning that read

analyzeline: Failed to evaluate '/3.e0+1j*( 1.e0+1j*( 2.e0/)'. Ignoring: invalid syntax (<string>, line 1)

I am not sure what to make of this. However, I do get 18 as an answer when I use x=1.0 and y=(/ 3, 6, 2/) (both using python and using the module in a Fortran program).

In short, completely avoid functions when you are using f2py.