Automatic generation of .pxd from Fortran

[Ondrej Certik]

Hi,

I wanted to share the code that I wrote:

https://gist.github.com/1013569

that takes a Fortran source (that uses the iso_c_binding module), like
the following:

module c_dftatom

! C interface to dftatom

use iso_c_binding, only: c_int, c_double, c_bool
use dftatom, only: get_Vxc, dp, solve_radial_eigenproblem, &
integrate, rpoisson_inward, get_atom_orb, atom_lda, mesh_exp, mesh_hyp, &
get_Vh2, mesh_exp_deriv, spline_interpolation, atom_rlda, &
get_atom_orb_rel, solve_radial_eigenproblem_perturb, &
integrate_radial_problem_outward
implicit none

contains

subroutine c_get_vxc(n, R, rho, relat, c, V) bind(c)
integer(c_int), intent(in) :: n
real(c_double), intent(in) :: R(n), rho(n), c
logical(c_bool), intent(in) :: relat
real(c_double), intent(out) :: V(n)
! We need to retype this, so that it compiles:
logical :: relat2
real(dp) :: exc(n)
relat2 = relat
call get_Vxc(R, rho, relat2, c, exc, V)
end subroutine

subroutine c_integrate(n, x, f, s) bind(c)
integer(c_int), intent(in) :: n
real(c_double), intent(in), dimension(n) :: x
real(c_double), intent(in), dimension(n) :: f
real(c_double), intent(out) :: s
s = integrate(x, f)
end subroutine

subroutine c_integrate_radial_poisson(n, R, Rp, rho, Z, V) bind(c)
integer(c_int), intent(in) :: n
real(c_double), intent(in) :: R(n), Rp(n), rho(n)
integer(c_int), intent(in) :: Z
real(c_double), intent(out) :: V(n)
V = rpoisson_inward(R, Rp, rho, Z)
end subroutine

....

end module

That you need to write by hand. The parser then produces the following file:

# This file is automatically generated by fparser from 'src/c_dftatom.f90'.
# Do not edit by hand (rerun fparser instead).

from libcpp cimport bool

cdef extern void c_get_vxc(int *n, double *R, double *rho, bool
*relat, double *c, double *V)
cdef extern void c_integrate(int *n, double *x, double *f, double *s)
cdef extern void c_integrate_radial_poisson(int *n, double *R, double
*Rp, double *rho, int *Z, double *V)
...

And then you just write Python wrappers by hand, as usual like:

from libcpp cimport bool

from numpy cimport import_array, ndarray
from numpy import empty

cimport rdirac

class ConvergeError(Exception):
pass

def getvxc(ndarray[double, mode="c"] R not None,
ndarray[double, mode="c"] rho not None,
bool relat, double c):
cdef int n = len(R)
assert len(rho) == n
cdef ndarray[double, mode="c"] V = empty(n)
rdirac.c_get_vxc(&n, &R[0], &rho[0], &relat, &c, &V[0])
return V

def integrate_radial_poisson(ndarray[double, mode="c"] density not None,
ndarray[double, mode="c"] R not None,
ndarray[double, mode="c"] Rp not None,
int Z):
cdef int n = len(density)
assert len(R) == n
cdef ndarray[double, mode="c"] V = empty(n)
rdirac.c_integrate_radial_poisson(&n, &R[0], &Rp[0], &density[0], &Z,
&V[0])
return V

...

The parsing code is super simple (see the gist link), everything is
explicit and clear what is happening. The code removes the need of
doing manual duplication (from Fortran C bindings to C and then again
to Cython). What you need to type by hand is unique -- e.g. export
your Fortran subroutines into C by writing a simple C wrapper using
iso_c_bindings, where you decide what the best C interface is. And
then you need to write the Python wrappers by hand, again deciding
what the best Python interface is.
I like this approach the most. Maybe it will be useful for someone.

In future versions, I should add generating .h files (besides .pxd).

Ondrej