F2PY examples#
Below are some examples of F2PY usage. This list is not comprehensive, but can be used as a starting point when wrapping your own code.
F2PY walkthrough: a basic extension module#
Creating source for a basic extension module#
Consider the following subroutine, contained in a file named add.f
C
SUBROUTINE ZADD(A,B,C,N)
C
DOUBLE COMPLEX A(*)
DOUBLE COMPLEX B(*)
DOUBLE COMPLEX C(*)
INTEGER N
DO 20 J = 1, N
C(J) = A(J)+B(J)
20 CONTINUE
END
This routine simply adds the elements in two contiguous arrays and places the result in a third. The memory for all three arrays must be provided by the calling routine. A very basic interface to this routine can be automatically generated by f2py:
python -m numpy.f2py -m add add.f
This command will produce an extension module named addmodule.c
in the
current directory. This extension module can now be compiled and used from
Python just like any other extension module.
Creating a compiled extension module#
Note
This usage depends heavily on numpy.distutils
, see F2PY and Build Systems
for more details.
You can also get f2py to both compile add.f
along with the produced
extension module leaving only a shared-library extension file that can
be imported from Python:
python -m numpy.f2py -c -m add add.f
This command produces a Python extension module compatible with your platform.
This module may then be imported from Python. It will contain a method for each
subroutine in add
. The docstring of each method contains information about
how the module method may be called:
>>> import add
>>> print(add.zadd.__doc__)
zadd(a,b,c,n)
Wrapper for ``zadd``.
Parameters
----------
a : input rank-1 array('D') with bounds (*)
b : input rank-1 array('D') with bounds (*)
c : input rank-1 array('D') with bounds (*)
n : input int
Improving the basic interface#
The default interface is a very literal translation of the Fortran code into
Python. The Fortran array arguments are converted to NumPy arrays and the
integer argument should be mapped to a C
integer. The interface will attempt
to convert all arguments to their required types (and shapes) and issue an error
if unsuccessful. However, because f2py
knows nothing about the semantics of
the arguments (such that C
is an output and n
should really match the
array sizes), it is possible to abuse this function in ways that can cause
Python to crash. For example:
>>> add.zadd([1, 2, 3], [1, 2], [3, 4], 1000)
will cause a program crash on most systems. Under the hood, the lists are being
converted to arrays but then the underlying add
function is told to cycle
way beyond the borders of the allocated memory.
In order to improve the interface, f2py
supports directives. This is
accomplished by constructing a signature file. It is usually best to start from
the interfaces that f2py
produces in that file, which correspond to the
default behavior. To get f2py
to generate the interface file use the -h
option:
python -m numpy.f2py -h add.pyf -m add add.f
This command creates the add.pyf
file in the current directory. The section
of this file corresponding to zadd
is:
subroutine zadd(a,b,c,n) ! in :add:add.f
double complex dimension(*) :: a
double complex dimension(*) :: b
double complex dimension(*) :: c
integer :: n
end subroutine zadd
By placing intent directives and checking code, the interface can be cleaned up quite a bit so the Python module method is both easier to use and more robust to malformed inputs.
subroutine zadd(a,b,c,n) ! in :add:add.f
double complex dimension(n) :: a
double complex dimension(n) :: b
double complex intent(out),dimension(n) :: c
integer intent(hide),depend(a) :: n=len(a)
end subroutine zadd
The intent directive, intent(out) is used to tell f2py that c
is
an output variable and should be created by the interface before being
passed to the underlying code. The intent(hide) directive tells f2py
to not allow the user to specify the variable, n
, but instead to
get it from the size of a
. The depend( a
) directive is
necessary to tell f2py that the value of n depends on the input a (so
that it won’t try to create the variable n until the variable a is
created).
After modifying add.pyf
, the new Python module file can be generated
by compiling both add.f
and add.pyf
:
python -m numpy.f2py -c add.pyf add.f
The new interface’s docstring is:
>>> import add
>>> print(add.zadd.__doc__)
c = zadd(a,b)
Wrapper for ``zadd``.
Parameters
----------
a : input rank-1 array('D') with bounds (n)
b : input rank-1 array('D') with bounds (n)
Returns
-------
c : rank-1 array('D') with bounds (n)
Now, the function can be called in a much more robust way:
>>> add.zadd([1, 2, 3], [4, 5, 6])
array([5.+0.j, 7.+0.j, 9.+0.j])
Notice the automatic conversion to the correct format that occurred.
Inserting directives in Fortran source#
The robust interface of the previous section can also be generated automatically by placing the variable directives as special comments in the original Fortran code.
Note
For projects where the Fortran code is being actively developed, this may be preferred.
Thus, if the source code is modified to contain:
C
SUBROUTINE ZADD(A,B,C,N)
C
CF2PY INTENT(OUT) :: C
CF2PY INTENT(HIDE) :: N
CF2PY DOUBLE COMPLEX :: A(N)
CF2PY DOUBLE COMPLEX :: B(N)
CF2PY DOUBLE COMPLEX :: C(N)
DOUBLE COMPLEX A(*)
DOUBLE COMPLEX B(*)
DOUBLE COMPLEX C(*)
INTEGER N
DO 20 J = 1, N
C(J) = A(J) + B(J)
20 CONTINUE
END
Then, one can compile the extension module using:
python -m numpy.f2py -c -m add add.f
The resulting signature for the function add.zadd is exactly the same
one that was created previously. If the original source code had
contained A(N)
instead of A(*)
and so forth with B
and C
,
then nearly the same interface can be obtained by placing the
INTENT(OUT) :: C
comment line in the source code. The only difference
is that N
would be an optional input that would default to the length
of A
.
A filtering example#
This example shows a function that filters a two-dimensional array of double precision floating-point numbers using a fixed averaging filter. The advantage of using Fortran to index into multi-dimensional arrays should be clear from this example.
C
SUBROUTINE DFILTER2D(A,B,M,N)
C
DOUBLE PRECISION A(M,N)
DOUBLE PRECISION B(M,N)
INTEGER N, M
CF2PY INTENT(OUT) :: B
CF2PY INTENT(HIDE) :: N
CF2PY INTENT(HIDE) :: M
DO 20 I = 2,M-1
DO 40 J = 2,N-1
B(I,J) = A(I,J) +
& (A(I-1,J)+A(I+1,J) +
& A(I,J-1)+A(I,J+1) )*0.5D0 +
& (A(I-1,J-1) + A(I-1,J+1) +
& A(I+1,J-1) + A(I+1,J+1))*0.25D0
40 CONTINUE
20 CONTINUE
END
This code can be compiled and linked into an extension module named filter using:
python -m numpy.f2py -c -m filter filter.f
This will produce an extension module in the current directory with a method
named dfilter2d
that returns a filtered version of the input.
depends
keyword example#
Consider the following code, saved in the file myroutine.f90
:
subroutine s(n, m, c, x)
implicit none
integer, intent(in) :: n, m
real(kind=8), intent(out), dimension(n,m) :: x
real(kind=8), intent(in) :: c(:)
x = 0.0d0
x(1, 1) = c(1)
end subroutine s
Wrapping this with python -m numpy.f2py -c myroutine.f90 -m myroutine
, we
can do the following in Python:
>>> import numpy as np
>>> import myroutine
>>> x = myroutine.s(2, 3, np.array([5, 6, 7]))
>>> x
array([[5., 0., 0.],
[0., 0., 0.]])
Now, instead of generating the extension module directly, we will create a signature file for this subroutine first. This is a common pattern for multi-step extension module generation. In this case, after running
python -m numpy.f2py myroutine.f90 -h myroutine.pyf
the following signature file is generated:
! -*- f90 -*-
! Note: the context of this file is case sensitive.
python module myroutine ! in
interface ! in :myroutine
subroutine s(n,m,c,x) ! in :myroutine:myroutine.f90
integer intent(in) :: n
integer intent(in) :: m
real(kind=8) dimension(:),intent(in) :: c
real(kind=8) dimension(n,m),intent(out),depend(m,n) :: x
end subroutine s
end interface
end python module myroutine
! This file was auto-generated with f2py (version:1.23.0.dev0+120.g4da01f42d).
! See:
! https://web.archive.org/web/20140822061353/http://cens.ioc.ee/projects/f2py2e
Now, if we run python -m numpy.f2py -c myroutine.pyf myroutine.f90
we see an
error; note that the signature file included a depend(m,n)
statement for
x
which is not necessary. Indeed, editing the file above to read
! -*- f90 -*-
! Note: the context of this file is case sensitive.
python module myroutine ! in
interface ! in :myroutine
subroutine s(n,m,c,x) ! in :myroutine:myroutine.f90
integer intent(in) :: n
integer intent(in) :: m
real(kind=8) dimension(:),intent(in) :: c
real(kind=8) dimension(n,m),intent(out) :: x
end subroutine s
end interface
end python module myroutine
! This file was auto-generated with f2py (version:1.23.0.dev0+120.g4da01f42d).
! See:
! https://web.archive.org/web/20140822061353/http://cens.ioc.ee/projects/f2py2e
and running f2py -c myroutine.pyf myroutine.f90
yields correct results.