# Foreign Function Interface

``````

In [1]:

%matplotlib inline
import matplotlib.pyplot as plt
plt.style.use('ggplot')

``````
``````

In [2]:

import numpy as np

``````

## Wrapping functions written in C

### Steps

• Write the C header and implementation files
• Write the Cython `.pxd` file to declare C function signatures
• Write the Cython `.pyx` file to wrap the C functions for Python
• Write `setup.py` to automate buiding of the Python extension module
• Run `python setup.py build_ext --inplace` to build the module
• Import module in Python like any other Python module

``````

In [3]:

%%file c_math.h

#pragma once
double plus(double a, double b);
double mult(double a, double b);
double square(double a);
double acc(double *xs, int size);

``````
``````

Writing c_math.h

``````

### C implementation file

``````

In [4]:

%%file c_math.c
#include <math.h>
#include "c_math.h"

double plus(double a, double b) {
return a + b;
};

double mult(double a, double b) {
return a * b;
};

double square(double a) {
return pow(a, 2);
};

double acc(double *xs, int size) {
double s = 0;
for (int i=0; i<size; i++) {
s += xs[i];
}
return s;
};

``````
``````

Writing c_math.c

``````

The `.pxd` file is similar to a header file for Cython. In other words, we can `cimport <filename>.pxd` in the regular Cython `.pyx` files to get access to functions declared in the `.pxd` files.

``````

In [5]:

%%file cy_math.pxd

cdef extern from "c_math.h":
double plus(double a, double b)
double mult(double a, double b)
double square(double a)
double acc(double *xs, int size)

``````
``````

Writing cy_math.pxd

``````

### Cython "implementation" file

Here is whhere we actually wrap the C code for use in Python. Note especially how we handle passing in of arrays to a C function expecting a pointer to double using `typed memoryviews`.

``````

In [6]:

%%file cy_math.pyx

cimport cy_math

def py_plus(double a, double b):
return cy_math.plus(a, b)

def py_mult(double a, double b):
return cy_math.mult(a, b)

def py_square(double a):
return cy_math.square(a)

def py_sum(double[::1] xs):
cdef int size = len(xs)
return cy_math.acc(&xs[0], size)

``````
``````

Writing cy_math.pyx

``````

### Build script `setup.py`

This is a build script for Python, similar to a Makefile

``````

In [7]:

%%file setup.py

from distutils.core import setup, Extension
from Cython.Build import cythonize
import numpy as np

ext = Extension("cy_math",
sources=["cy_math.pyx", "c_math.c"],
libraries=["m"],
extra_compile_args=["-w",  "-std=c99"])

setup(name = "Math Funcs",
ext_modules = cythonize(ext))

``````
``````

Writing setup.py

``````

### Building an extension module

``````

In [8]:

! python setup.py clean
! python setup.py -q build_ext --inplace

``````
``````

Compiling cy_math.pyx because it changed.
[1/1] Cythonizing cy_math.pyx
running clean

``````
``````

In [9]:

! ls cy_math*

``````
``````

cy_math.c                     cy_math.pxd
cy_math.cpython-35m-darwin.so cy_math.pyx

``````

### Using the extension module in Python

``````

In [10]:

import cy_math
import numpy as np

print(cy_math.py_plus(3, 4))
print(cy_math.py_mult(3, 4))
print(cy_math.py_square(3))

xs = np.arange(10, dtype='float')
print(cy_math.py_sum(xs))

``````
``````

7.0
12.0
9.0
45.0

``````

### Confirm that we are getting C speedups by comparing with pure Python accumulator

``````

In [11]:

def acc(xs):
s = 0
for x in xs:
s += x
return s

``````
``````

In [12]:

import cy_math

xs = np.arange(1000000, dtype='float')
%timeit -r3 -n3 acc(xs)
%timeit -r3 -n3 cy_math.py_sum(xs)

``````
``````

3 loops, best of 3: 197 ms per loop
3 loops, best of 3: 1.95 ms per loop

``````

## C++

This is similar to C. We will use Cython to wrap a simple funciton.

``````

In [13]:

#pragma once

``````
``````

``````
``````

In [14]:

int add(int a, int b) {
return a+b;
}

``````
``````

``````
``````

In [15]:

%%file plus.pyx

def plus(a, b):

``````
``````

Writing plus.pyx

``````

#### Note that essentially the only difference from C is `language="C++"` and the flag `-std=c++11`

``````

In [16]:

%%file setup.py
from distutils.core import setup, Extension
from Cython.Build import cythonize

ext = Extension("plus",
extra_compile_args=["-w", "-std=c++11"])

setup(
ext_modules = cythonize(
ext,
language="c++",
))

``````
``````

Overwriting setup.py

``````
``````

In [17]:

%%bash
python setup.py -q build_ext --inplace

``````
``````

Compiling plus.pyx because it changed.
[1/1] Cythonizing plus.pyx

error: invalid argument '-std=c++11' not allowed with 'C/ObjC'
error: command '/usr/bin/clang' failed with exit status 1

``````
``````

In [18]:

import plus

plus.plus(3, 4)

``````
``````

---------------------------------------------------------------------------
ImportError                               Traceback (most recent call last)
<ipython-input-18-4d2911b32d13> in <module>()
----> 1 import plus
2
3 plus.plus(3, 4)

ImportError: No module named 'plus'

``````

## Wrap an R function from libRmath using `ctypes`

R comes with a standalone C library of special functions and distributions, as described in the official documentation. These functions can be wrapped for use in Python.

### Building the Rmath standalone library

```git clone https://github.com/JuliaLang/Rmath-julia.git
cd Rmath-julia/src
make
cd ../..
```

#### Functions to wrap

``````

In [ ]:

! grep "\s.norm(" Rmath-julia/include/Rmath.h

``````
``````

In [ ]:

from ctypes import CDLL, c_int, c_double

``````
``````

In [ ]:

%%bash
ls Rmath-julia/src/*so

``````
``````

In [ ]:

lib = CDLL('Rmath-julia/src/libRmath-julia.so')

def rnorm(mu=0, sigma=1):
lib.rnorm.argtypes = [c_double, c_double]
lib.rnorm.restype  = c_double
return lib.rnorm(mu, sigma)

def dnorm(x, mean=0, sd=1, log=0):
lib.dnorm4.argtypes = [c_double, c_double, c_double, c_int]
lib.dnorm4.restype  = c_double
return lib.dnorm4(x, mean, sd, log)

def pnorm(q, mu=0, sd=1, lower_tail=1, log_p=0):
lib.pnorm5.argtypes = [c_double, c_double, c_double, c_int, c_int]
lib.pnorm5.restype  = c_double
return lib.pnorm5(q, mu, sd, lower_tail, log_p)

def qnorm(p, mu=0, sd=1, lower_tail=1, log_p=0):
lib.qnorm5.argtypes = [c_double, c_double, c_double, c_int, c_int]
lib.qnorm5.restype  = c_double
return lib.qnorm5(p, mu, sd, lower_tail, log_p)

``````
``````

In [ ]:

pnorm(0, mu=2)

``````
``````

In [ ]:

qnorm(0.022750131948179212, mu=2)

``````
``````

In [ ]:

plt.hist([rnorm() for i in range(100)])
pass

``````
``````

In [ ]:

xs = np.linspace(-3,3,100)
plt.plot(xs, list(map(dnorm, xs)))
pass

``````

### Using Cython to wrap standalone library

``````

In [ ]:

%%file rmath.pxd

cdef extern from "Rmath-julia/include/Rmath.h":
double dnorm(double, double, double, int)
double pnorm(double, double, double, int, int)
double qnorm(double, double, double, int, int)
double rnorm(double, double)

``````
``````

In [ ]:

%%file rmath.pyx

cimport rmath

def rnorm_(mu=0, sigma=1):
return rmath.rnorm(mu, sigma)

def dnorm_(x, mean=0, sd=1, log=0):
return rmath.dnorm(x, mean, sd, log)

def pnorm_(q, mu=0, sd=1, lower_tail=1, log_p=0):
return rmath.pnorm(q, mu, sd, lower_tail, log_p)

def qnorm_(p, mu=0, sd=1, lower_tail=1, log_p=0):
return rmath.qnorm(p, mu, sd, lower_tail, log_p)

``````
``````

In [ ]:

%%file setup.py
from distutils.core import setup, Extension
from Cython.Build import cythonize

ext = Extension("rmath",
sources=["rmath.pyx"],
include_dirs=["Rmath-julia/include"],
library_dirs=["Rmath-julia/src"],
libraries=["Rmath-julia"],
runtime_library_dirs=["Rmath-julia/src"],
extra_compile_args=["-w",  "-std=c99", "-DMATHLIB_STANDALONE"],
)

setup(
ext_modules = cythonize(
ext
))

``````
``````

In [ ]:

! python setup.py build_ext --inplace

``````
``````

In [ ]:

import rmath

plt.hist([rmath.rnorm_() for i in range(100)])
pass

``````
``````

In [ ]:

xs = np.linspace(-3,3,100)
plt.plot(xs, list(map(rmath.dnorm_, xs)))
pass

``````

### `Cython` wrappers are faster than `ctypes`

``````

In [ ]:

%timeit pnorm(0, mu=2)
%timeit rmath.pnorm_(0, mu=2)

``````

## Fortran

``````

In [ ]:

! pip install fortran-magic

``````
``````

In [ ]:

``````
``````

In [ ]:

%%fortran

subroutine fort_sum(N, s)
integer*8, intent(in) :: N
integer*8, intent(out) :: s
integer*8 i
s = 0
do i = 1, N
s = s + i*i
end do
end

``````
``````

In [ ]:

fort_sum(10)

``````

#### Another example from the documentation

``````

In [ ]:

subroutine solve(A, b, x, n)
! solve the matrix equation A*x=b using LAPACK
implicit none

real*8, dimension(n,n), intent(in) :: A
real*8, dimension(n), intent(in) :: b
real*8, dimension(n), intent(out) :: x

integer :: pivot(n), ok

integer, intent(in) :: n
x = b

! find the solution using the LAPACK routine SGESV
call DGESV(n, 1, A, n, pivot, x, n, ok)

end subroutine

``````
``````

In [ ]:

A = np.array([[1, 2.5], [-3, 4]])
b = np.array([1, 2.5])

solve(A, b)

``````

## Interfacing with R

``````

In [ ]:

``````
``````

In [ ]:

%%R
library(ggplot2)
suppressPackageStartupMessages(
ggplot(mtcars, aes(x=wt, y=mpg)) + geom_point() + geom_smooth(method=loess)
)

``````

#### Converting between Python and R

``````

In [ ]:

%R -o mtcars

``````

#### `mtcars` is now a Python dataframe

``````

In [ ]:

``````

#### We can also pass data from Python to R

``````

In [ ]:

x = np.linspace(0, 2*np.pi, 100)
y = np.sin(x)

``````
``````

In [ ]:

%%R -i x,y
plot(x, y, main="Sine curve in R base graphics")

``````
``````

In [ ]:

``````