Use f2py with code having real(kind=wp) declarations

Can f2py be used to call Fortran code with declarations of the form `real(kind=wp)` from Python? For example, calling the subroutines in the Fortran module

``````module foo
implicit none
private
public :: sum_absolute, twice
integer, parameter :: wp = kind(1.0d0)
contains
pure subroutine sum_absolute(n, x, sum_abs_x)
! return the sum of the absolute values of x
integer      , intent(in)  :: n
double precision, intent(in)  :: x(n)
double precision, intent(out) :: sum_abs_x
!f2py intent(in) n
!f2py intent(in) x(n)
!f2py intent(out) sum_abs_x
sum_abs_x = sum(abs(x))
end subroutine sum_absolute
!
pure subroutine twice(n, x, twice_x)
integer      , intent(in)  :: n
double precision, intent(in)  :: x(n)
double precision, intent(out) :: twice_x(n)
!f2py intent(in) n
!f2py intent(in) x(n)
!f2py intent(in,out) twice_x(n)
twice_x = 2*x
end subroutine twice
end module foo
``````

with the Python script

``````import numpy as np
from bar import foo
x = np.array([-10, 20, -30], dtype=np.float64)
twice_x = np.zeros(len(x))
foo.twice(x, twice_x)
print("x, sum_absolute(x), twice(x) =", x, foo.sum_absolute(x), twice_x)
``````

using the commands

``````f2py -c -m bar sub.f90 --compiler=mingw32
python xsub.py
``````

works on Windows, giving output

`x, sum_absolute(x), twice(x) = [-10. 20. -30.] 60.0 [-20. 40. -60.]`

However, when the Fortran module has variables declared with a parameterized kind instead of using `double precision`,

``````module foo
implicit none
private
public :: sum_absolute, twice
integer, parameter :: wp = kind(1.0d0)
contains
pure subroutine sum_absolute(n, x, sum_abs_x)
! return the sum of the absolute values of x
integer      , intent(in)  :: n
real(kind=wp), intent(in)  :: x(n)
real(kind=wp), intent(out) :: sum_abs_x
!f2py intent(in) n
!f2py intent(in) x(n)
!f2py intent(out) sum_abs_x
sum_abs_x = sum(abs(x))
end subroutine sum_absolute
!
pure subroutine twice(n, x, twice_x)
integer      , intent(in)  :: n
real(kind=wp), intent(in)  :: x(n)
real(kind=wp), intent(out) :: twice_x(n)
!f2py intent(in) n
!f2py intent(in) x(n)
!f2py intent(in,out) twice_x(n)
twice_x = 2*x
end subroutine twice
end module foo
``````

The output of
`x, sum_absolute(x), twice(x) = [-10. 20. -30.] -10.0 [0. 0. 0.]`

is wrong.

1 Like

It seems you are using the quick way rather than the smart way. Maybe you can try the latter and inspect the signature file.

Iâ€™ve only played a couple of times with f2py. I prefer f90wrap, which builds on top of f2py and has support for derived data types. f90wrap has a specific option `-k KIND_MAP` that alows the user to define a dictionary with the kind map, for example:

``````{
'real':     {'': 'float',
'8': 'double',
'dp': 'double',
'idp': 'double'},
'complex': {'': 'complex_float',
'8': 'complex_double',
'16': 'complex_long_double',
'dp': 'complex_double'},
'integer': {'': 'long_long',
'4': 'int',
'8': 'long_long',
'dp': 'long_long'},
'character': {'': 'char'}
}
``````

Is it intentional that the example of the smart way document linked above, though tagged with a copyright of 2008-2022 still uses FORTRAN? The same actually happens to all examples on this handout, which can contribute to the perception â€śFORTRAN = oldâ€ť (as in the meaning of outdated).

And what a pity - impossible for me to get either one of the three examples presented (by the numpy page) to work. A â€śthere must be a better wayâ€ť moment (I think quoting Raymond Hettinger here).

numpy_log.txt (24.4 KB)

I had the same problem. For me it works if you declare wp in each subroutine repeatedly. Not very elegant but it worksâ€¦

module foo
implicit none
private
public :: sum_absolute, twice
contains
pure subroutine sum_absolute(n, x, sum_abs_x)
! return the sum of the absolute values of x
integer, parameter :: wp = kind(1.0d0)
integer , intent(in) :: n
real(kind=wp), intent(in) :: x(n)
real(kind=wp), intent(out) :: sum_abs_x
!f2py intent(in) n
!f2py intent(in) x(n)
!f2py intent(out) sum_abs_x
sum_abs_x = sum(abs(x))
end subroutine sum_absolute
!
pure subroutine twice(n, x, twice_x)
integer, parameter :: wp = kind(1.0d0)
integer , intent(in) :: n
real(kind=wp), intent(in) :: x(n)
real(kind=wp), intent(out) :: twice_x(n)
!f2py intent(in) n
!f2py intent(in) x(n)
!f2py intent(in,out) twice_x(n)
twice_x = 2*x
end subroutine twice
end module foo

2 Likes

The obstacle was reported as an issue relevant to numpyâ€™s documentation (ticket on GitHub).

@Beliavsky I spent some time today reading more about f2py. It turns out that (just like f90wrap) it also relies on a user-defined dictionary to deal with kind specifiers. By default, the dictionary file name is `.f2py_f2cmap`.

You can find a working example here: HugoMVale/fortran-in-python: Examples of how to build/call fortran modules in python (github.com)

Update: My ticket on GitHub was resolved successfully, the examples of all three sections on numpyâ€™s page now can be replicated 1:1.

There are Python header files provided by an optional additional package `python3-dev` (name in Linux Debian 12/bookworm, currently version 3.11.2-1+1b). An installation may require update and installation of additional dependencies synaptic (or an other package manager of your preference) resolves on the fly. Yet then, the obstacle previously reported is absent. (For now, numpyâ€™s documentation doesnâ€™t explicitly state this dependency.)