Using f2py

The purpose of the F2PY --Fortran to Python interface generator-- is to provide connection between Python and Fortran languages. Here I will show a simple example of how to use this tool and compare the its efficiency against a pure python code.

First, lets write our fortran routine:

        Subroutine test(a,b,n,c)

        Implicit none
        Integer::i
        Integer, intent(in)::n
        Real, dimension(n), intent(in)::a,b
        Real, dimension(n), intent(out)::c
!f2py   depend(n)::a,b,c       
!f2py   intent(in)::a,b
!f2py   intent(out)::c

        Do i=1,n
          c(i) = sqrt(a(i)**2.0 + b(i)**2.0)
        EndDo

        End subroutine test

This code simply calculate the distance to the origin of a bunch of points. The input is two vectors with identical dimension n. Note that you must specify what goes in and what goes out to the f2py parser, this is done by using the a comment line followed by f2py:

!f2py   depend(n)::a,b,c       
!f2py   intent(in)::a,b
!f2py   intent(out)::c

Now lets compile and link it to python using:

f2py -m func -h func.pyf func.f --overwrite-signature
f2py --fcompiler=gfortran -c func.pyf func.f

Unfortunately f2py does not easily support the Intel Fortran Compiler (ifort), although gfortran does a good job on most of the cases.

It ridiculously easy to call this routine inside python. The script bellow call the fortran routine and then compares its execution time with an identical pure python routine plus the same done using Numpy.

#!/usr/bin/python
import numpy as np
import time
import func

n = 10000
#generates a normal distribution with n points
x = np.random.randn(n)
y = np.random.randn(n)

# Using the fortran routine
t1 = time.time()
f = func.test(x,y,n)
t2 = time.time()

t = (t2-t1)*1000

print 'Fortran runtime (ms): ',t

# Using pure python
t1 = time.time()
r = np.zeros(n)
for i in range(n):
 r[i] = np.sqrt(x[i]**2 + y[i]**2)
t2 = time.time()

t = (t2-t1)*1000

print 'Python runtime (ms): ',t

# Using pure python
t1 = time.time()
r = np.sqrt(x**2 + y**2)
t2 = time.time()
t = (t2-t1)*1000

print 'Python + numpy runtime (ms): ',t

If everything goes OK the result should be something like this:

Fortran runtime (ms):  0.150918960571
Python runtime (ms):  76.4260292053
Python + numpy runtime (ms):  0.485181808472

The advantages of using f2py are obvious, the runtime of the pure python code is almost 500 times greater! When using numpy to handle the arrays we get a much better execution time, but still slower when compared to fortran.

I hope this post was useful to the people thinking of migrating to python in the near future.