Performant Python

C/Fortran Interface

Neelofer Banglawala, Arno Proeme, Kevin Stratford, Andy Turner

Why couple python with another language?

  • Provide glue to organise complex tasks

    • Handle complex software coordination provided by Python

  • Combine performance of compiled codes with flexibility of Python

    • e.g. incorporate Python analysis and visualisation into existing codebase
    • Provide flexible way to extract results from code using Python

  • Reuse code that you already have

    • Gradually introduce new functionality using Python

Extension modules

  • Basic approach is to compile extension modules

    • Compiles native langauge source
    • Adds a wrapper which provides the interface

  • Requires

    • Appropriate compiler (e.g., gfortran, cc, ...)
    • A clear understanding of the number of types of any arguments

  • Will produce

    • Extension module (shared library .so) loaded at run time
    • import extension module as usual
  • Some care may be required with compound/opaque data types

Different approaches

  • Fortran

    • f2py is part of Numpy will handle external subroutines
    • Modern Fortran requires, e.g., f90wrap

  • C (more choice)

    • f2py can be used (a kludge)
    • C native interface
    • Cython
    • Swig
    • ctypes/CFFI

  • Various other approaches

    • Weave, Numba

Fortran and python via f2py

  • You need to provide f2py with:

    • Fortran source code
    • A signature file: a file describing the external function and its arguments
    • (f2py can help you generate a signature file)
  • Also need access to a Fortran compiler

  • f2py can:

    • Create a signature file containing argument attributes (e.g., optional) that define the Fortran interface
    • Wrap Fortran code in an extension module that can be imported from within Python

General recipe

  • Create a signature file
    • Use 
f2py <source> -m <extension_module> -h <signature>.pyf
  • Typically the signature filename stub is the same as the source filename
  • Check the signature file for correctness
    • Sequence and types of arguments to be passed from Python to Fortran function
    • Argument attributes, such as depend
  • Produce the final extension module
    • f2py -c <signature_file>.pyf <source_file>.f90
  • Import module into Python and use the external Fortran function!
import extension_module_name
    extension_module_name.function(args)
  • The source filename may not be the same as the function name

Example farray_sqrt.f90

Consider:

array_sqrt() is an external subroutine

subroutine array_sqrt(n, a_in, a_out)
  implicit none
  integer, intent(in) :: n
  real, dimension(n), intent(in)  :: a_in
  real, dimension(n), intent(out) :: a_out

  integer :: i

  do i = 1, n
    a_out(i) = sqrt(a_in(i))
  end do

  return
end subroutine array_sqrt

Create a signature file

  • f2py can try to create the signature file (farray_sqrt.pyf) automatically

    • From a terminal, issue the command:
f2py farray_sqrt.f90 -m farray -h farray_sqrt.pyf

  • The Python extension module will be called: farray

    • use the -m option

  • Signature in text file called: farray_sqrt.pyf

    • use the -h option

    • Note: will not overwrite an existing signature file: Use --overwrite-signature to overwrite.

In [ ]:
# Call from within Python to save exiting notebook...

!f2py ./code/farray_sqrt.f90 -m farray -h farray_sqrt.pyf

Check signature file

Attributes such as optional, intent and depend specify the visibility, purpose and dependencies of the arguments.

In [ ]:
!cat farray_sqrt.pyf

Note: exact result can depend on version of numpy, so it is worth checking

Compile extension module

Once you have verified that the signature file is correct

  • Use f2py to compile a module file that can be imported into Python
f2py -c farray_sqrt.pyf farray_sqrt.f90

This should produce a shared library file called: farray.so

In [ ]:
# Run f2py command from within notebook
# If you don't want to see the output, try "msg = !f2py ..."

!f2py -c farray_sqrt.pyf ./code/farray_sqrt.f90
In [ ]:
# Check we have the farray.so
!ls
In [ ]:
!cat farray_sqrt.pyf

Call external function from Python

In [ ]:
# import the extension module
import numpy as np
from farray import array_sqrt
In [ ]:
# view docsting of function (automatically produced)
array_sqrt?
In [ ]:
# let's try to use the function

ain = np.array([1.0, 4.0, 9.0, 16.0, 2.0], np.double)
print(ain)
aout = array_sqrt(ain)
print(aout)

Python and C via ctypes

  • Uses only python (no additional interface file or mixed-language intermediate code)

    import ctypes

  • Must load the library (.so) file explicitly

lib = ctypes.cdll.LoadLibrary("./my_library.so")
  • Must specify the prototype for the C function
lib.my_c_function.restype = ctypes.c_int
    lib.my_c_function.argtypes = [ctypes.c_double]

Example C side

Consider the simple function:

int my_c_function(double val) {

    return (int) (val + 1.0);
}

We need to compile an extension module:

In [ ]:
!gcc -c -fPIC ./code/my_library.c
!gcc -shared -o my_library.so my_library.o
!ls

Example python side

Now:

In [ ]:
import ctypes

lib = ctypes.cdll.LoadLibrary("./my_library.so")

lib.my_c_function.restype = ctypes.c_int
lib.my_c_function.argtypes = [ctypes.c_double]

x = float(23)
result = lib.my_c_function(x)
print(result, type(result))

ctypes and numpy.ndarray

Consider again the square root example, this time in C:

#include <math.h>

void array_sqrt(int n, double * a_in, double * a_out) {

  int i;

  for (i = 0; i < n; i++) {
    a_out[i] = sqrt(a_in[i]);
  }

  return;
}
In [ ]:
!gcc -c -fPIC ./code/c_sqrt.c
!gcc -shared -o c_sqrt.so c_sqrt.o

Using numpy.ctypeslib

The corresponding ctypes code must address the two C pointers.

In [ ]:
import ctypes
import numpy
from numpy.ctypeslib import ndpointer

lib = ctypes.cdll.LoadLibrary("./c_sqrt.so")
lib.array_sqrt.restype = None
lib.array_sqrt.argtypes = [ctypes.c_int, \
            ndpointer(ctypes.c_double, flags="C_CONTIGUOUS"), \
            ndpointer(ctypes.c_double, flags="C_CONTIGUOUS")]

a_in = numpy.array([16.0, 25.0, 36.0, 49.0])
a_out = numpy.empty(4, np.double)

lib.array_sqrt(4, a_in, a_out)
print(a_out)

Alternatives

  • Native Python interface

    • Fully-flexible and portable
    • Complex and verbose
    • Option if you are interfacing a large amount of code and/or have a large software development project

  • Cython : converts Python-like code into a C library which can call other C libraries

    • Standard C-like Python (or Python-like C)

  • SWIG (or Simplified Wrapper and Interface Generator) : reads header files and generates a library Python can load

    • Very generic and feature-rich
    • Supports multiple languages other than Python (e.g. Perl, Ruby)

Alternatives contd ...

  • ctypes, cffi (C Foreign Function Interface for Python) : both provide "foreign function interfaces", or lightweight APIs, for calling C libraries from within Python

  • The goal is to provide a convenient and reliable way to call compiled C code from Python using interface declarations written in C

  • Weave : includes C/C++ code within Python code and compiles it transparently

  • Boost.python : helps write C++ libraries that Python can load and use easily

  • PyCUDA : allows you to include NVIDIA CUDA code within Python. You can also write C code by hand, that can be called by Python.

Summary

  • Calling C/Fortran allows code re-use
  • Fortran/C can give better performance than Python
  • f2py is a simple way to call Fortran code from Python
  • Modern Fortran users should consider f90wrap
  • Other options more appropriate for C

https://github.com/jameskermode/f90wrap

In [ ]:
help(ctypes)

Exercise fibonacci.f90

Use f2py to create an extension module for function fibonacci() at ./code/fibonacci.f90 and test it in Python.

fibonacci() computes the first n Fibonacci numbers: 0, 1, 1, 2, 3, 5, 8, 13,... and stores the results in the array provided.

subroutine fibonacci(n, a_out)
  implicit none
  integer, intent(in) :: n
  real*8, dimension(n) :: a_out

  integer :: i

  do i = 1, n
    if (i.eq.1) then
      a_out(i) = 0.0
    else if (i.eq.2) then
      a_out(i) = 1.0
    else
      a_out(i) = a_out(i-1) + a_out(i-2)
    end if
  end do
end subroutine fibonacci
In [ ]:
!ls

C Exercise

An equivalent C function is available to compute a Fibonacci series: ./code/fibonacci.c

Previous

Schedule

Next