In [ ]:
from sympy import *
init_printing(use_latex='mathjax')

x = symbols('x')

Code Generation

In [ ]:
expr = R_nl(3, 1, x, 6)
expr
In [ ]:
from sympy.utilities.autowrap import ufuncify
from sympy.utilities.lambdify import lambdify

fn_numpy   = lambdify(x, expr, 'numpy')
fn_fortran = ufuncify([x], expr)

They work the same

In [ ]:
from numpy import linspace
xx = linspace(0, 1, 5)
fn_numpy(xx)
In [ ]:
fn_fortran(xx)

Plot form with matplotlib

In [ ]:
from matplotlib.pyplot import plot, show, legend
%matplotlib inline
In [ ]:
xx = linspace(0, 5, 50000)  # A bigger 
plot(xx, fn_numpy(xx))
In [ ]:
# Time fn_numpy
In [ ]:
# Time fn_fortran

Theano

In [ ]:
import theano.tensor as T
import theano
In [ ]:
tx = T.vector('x')
ty = T.vector('y')
tz = tx + 2 * ty
tz
In [ ]:
# Print z, see expression
In [ ]:
# Make function
f = theano.function([tx, ty], tz)
In [ ]:
import numpy as np
xx = np.array([1, 2, 3], dtype=np.float32)
yy = np.array([10, 20, 30], dtype=np.float32)
f(xx, yy)

Theano compiles symbolic trees down to numeric C and CUDA code

Less sophisticated mathematics but more sophisticated code generation

(Almost) all SymPy expressions can be transformed into Theano Expressions.

In [ ]:
a = Symbol('a')
from sympy.printing.theanocode import theano_code
theano_code(a)
In [ ]:
 
In [ ]:
 
In [ ]:
theano.printing.debugprint(theano_code(expr))
In [ ]:
from sympy.printing.theanocode import theano_function

fn_theano  = theano_function([x], [expr], dims={x: 1}, dtypes={x: 'float64'})
In [ ]:
theano.printing.debugprint(fn_theano)
In [ ]:
fn_theano(xx)

More complex problem

Lets compute both the expression and its derivative simultaneously.

In [ ]:
outputs = expr, simplify(expr.diff(x))
outputs
In [ ]:
fn_numpy  = lambdify([x], outputs, 'numpy')
fn_theano = theano_function([x], outputs, dims={x: 1}, dtypes={x: 'float64'})

fns_fortran = [ufuncify([x], output) for output in outputs]
fn_fortran  = lambda xx: [fn_fortran(xx) for fn_fortran in fns_fortran]
In [ ]:
xx = linspace(0, 5, 50000)

for y in fn_theano(xx):
    plot(xx, y)
legend(['$R_{31}$', "$R'_{31}$"])
In [ ]:
timeit fn_numpy(xx)
In [ ]:
timeit fn_fortran(xx)
In [ ]:
timeit fn_theano(xx)

Matrices

In [ ]:
n = Symbol('n', integer=True)
X = MatrixSymbol('X', n, n)
y = MatrixSymbol('y', n, 1)

f = theano_function([X, y], [X*y])
In [ ]:
nX = np.array([[1, 0], [0, 1]], dtype=np.float32)
ny = np.array([[10], [1]], dtype=np.float32)

f(nX, ny)
In [ ]: