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

Help
1

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
2

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'}
}
3

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)

4

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
5

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

6

@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)

7

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.)