Passing variables to Fortran from python causes wrong values

Question

Based on this minimal example, I want to manipulate variables of a Fortran module with Python and Fortran subroutines. Please take a look at the following example:

vars.f90

module vars
implicit none
real(kind=selected_real_kind(p=15))      :: fk(10)
end module vars

sub.f90

subroutine sub
    use vars
    print *, "sub: fk(1) = ", fk(1)
    print *, "adding 1 to fk(1)"
    fk(1) = fk(1) + 1
    print *, "fk(1) = ", fk(1)
end

mytest.f90

include "vars.f90"
include "sub.f90"

The compilation is done with the following command and gfortran:

f2py -c -m mytest mytest.f90

Finally, here is the testcase to reproduce the problem in a Python 3.6.5 console:

>>> import mytest
>>> mytest.vars.fk[1]
0.0
>>> mytest.vars.fk[1]=1.5000
>>> mytest.vars.fk[1]
1.5
>>> mytest.sub()
 sub: fk(1) =   0.12500000000000000     
 adding 1 to fk(1)
 fk(1) =    1.1250000000000000     
>>> mytest.vars.fk[1]
1.890625
>>> 

From my point of view, fk(1) should be 2.5 in the end. But unfortunately Fortran reads the variable in the subroutine wrong although the direct access to the variable via python console displays the correct value. After modifying the variable the python console displays a wrong number too.

Any advice or proposal to resolve/reason this behavior? I appreciate any help! Thank you!

P.S.: First I thought it depends on the way Python and Fortran handle the arrays differently(Python starts at index 0 and Fortran at index 1), but this is not the mistake. mytest.vars.fk[0] is not modified after calling the subroutine (initial value 0.0).

Answer 1

There are two issues with the code:

1. The indexing isssue, that you have noted. You still need to access fk[0] in Python to obtain the desired behaviour.

2. The use of selected_real_kind(p=15) for the kind specification. See the FAQ of F2PY and this related SO answer.

Currently, f2py sees your varibale as a default real variable (here, 32 bit) instead of a double precision. You can fix it by creating a file .f2py_f2cmap (in the current working directory) with the content:

{'real':{'selected_real_kind(p=15)': 'double'}}

that will select "double" as the corresponding C type for your parameter. F2PY has no way of detecting modern Fortran types, unfortunately.

I flagged this as duplicate of f2py: Specifying real precision in fortran when interfacing with python? but there is also the indexing issue here. Let's see what others think.