GSoC 2017 Report


My project was to enhance the code-generation facilities of SymPy. You can read my proposal for the motivation behind this work. The overall goals were the following:

  • Allow to render code targeting a specific precision (e.g. binary32 vs. binary64 floating point numbers). Prior to this project the printers would sometimes generate code containing a mixture of single, double and extended precision, and there were no way to change this short of subclassing the printers and overriding the methods.
  • Allow to render blocks of code and not only expressions. There was an initial effort to support this in the submodule sympy.codegen.
  • Improve the Lambdify functionality in SymEngine's python wrapper. Before this project it did not handle outputs of mixed shape, and it also had considerable overhead.

Summary of the final work product

A whole new repository with code and notebooks for code-generation was created during the first part of GSoC:

Jason Moore, Kenneth Lyons, Aaron Meurer and I (my commits) created this for the tutorial in code generation with SymPy at the SciPy 2017 conference.

The majority of the work are contained in these pull-requests:

  •, #112 (merged): Heterogeneous output in Lambdify.
  •, #171 (merged): Bug fix for heterogeneous output in Lambdify.
  • sympy, #12693 (merged): Extending the sympy.codegen.ast module with new classes (for generating ASTs).
  • sympy, #12808 & #13046 (merged): PythonCodePrinter, MpmathPrinter, SymPyPrinter NumPyPrinter, SciPyPrinter.
  • sympy, #13194 (open): Add .codegen.rewriting module.
  • sympy, #13200 (merged): Add .codegen.approximations module.
  • sympy, #13100 (open): More AST nodes. Building on #12693, this is the biggest PR. In addition to improving the AST nodes it introduces .codegen.algorithms as well as an internal testing module .utilities._compilation which allows to compile and import/run strings of C/C++/Fortran code.

In addition there were smaller pull-requests made & merged:

  • sympy_benchmarks: #37, #38, #39, #40: Benchmarks for lambidfy and common sub-expression elimination.
  • sympy: #12686 (Support for __abs__ in SymPy matrices), #12692 (subclass support for SymPy's deprecation decorator), #12762 (Fix floating point error under windows), #12805 Revert change to cse (performance regression), #12764 environment variable use, #12833 string formatting, #12944 allow relative path in autowrap, #13063 fix test timing script (and updated timings), #12833 (some of the commits) Allow custom class in autowrap & codegen, #12843 (one of the commits) allow changing compile arguments in CythonCodeWrapper.

Detailed review of the work

The first weeks of the summer was mostly spent on the code generation material presented at the SciPy conference tutorial, in parallel with that work was done to handle different choices of data types in the printers. And new AST nodes were introduced to represent type.

Code for the tutorial

During the writing of this code improvements were made to the existing code-generation facilities and SymPy (and experience with their shortcomings were gained). One of the challenges in this work was that the attendees at the conference would be using all major platforms (Linux/macOS/Windows) and different Python versions, we needed to ensure that generating code, compiling, linking and importing worked all combinations.

Lambdify in SymEngine's python wrapper

Writing the code for the tutorial provided great test cases for the code-generation capabilities of SymPy. The motivation of doing code generation is usually that of speed (but sometimes it may be motivated by wanting to work with some library written in another language). An alternative to generating high level code which then gets compiled, is to go toward assembly (or some intermediate representation). SymEnigne had support for doing this via LLVM's JIT compiler. The Python bindings however needed an overhaul (something I had included in the time-line in my proposal), and now I wanted to use Lambdify (the SymEngine version of sympy.lambdify), and together with the help of Isuru Fernando we got it to work (and benchmarks for pydy show that it is even faster than using the cython backend).

AST nodes

I had made AST nodes in my prototype for my proposal, right at the start of the project I ported those to SymPy. It took some rewriting and discussion with Aaron (both during our weekly meetings and at the conference) to get it to a point where we were confident enough to merge it into SymPy's codebase.

One of the major challanges when designing the new classes for sympy.codegen.ast was dealing with optional arguments in our subclasses of symyp.core.basic.Basic. The solutions which worked best was to have a subclass sympy.codegen.Node which stored such optinoal information as instances in a SymPy Tuple as its last argument (accessible as .attrs`). This allowed the code printers for Python, C and Fortran to support the same ``Variable class for instance, where the C printer would also look for attributes "value_const", "volatile" etc. and the Fortran printer would look for e.g. "intent".

Language specific nodes have been added under their own submodules in sympy.codegen (e.g. sympy.codegen.fnodes for Fortran and sympy.codegen.cnodes for C). The most common statements are now implmeneted, but the nodes are by far not exhaustive. There are now also helper functions for generating e.g. modules in sympy.codegen.pyutils & sympy.codegen.futils (for Python and Fortran respectively).

Code printers

Dealing with floating point types is tricky since one want to be pragmatic in order for the types to be helpful (IEEE 754 conformance is assumed), but general enough that people targeting hardware with non-standard conformance can still generate useful code using SymPy. For example, one can now choose the targeted precision:

>>> from sympy import ccode, symbols, Rational
>>> x, tau = symbols("x, tau")
>>> expr = (2*tau)**Rational(7, 2)
>>> from sympy.codegen.ast import real, float80
>>> ccode(expr, type_aliases={real: float80})
'8*M_SQRT2l*powl(tau, 7.0L/2.0L)'

Here we have assumed that the targeted architechture has x87 FPU (long double is a 10 byte extended precision floating point data type). But it is fully possible to generate code for some other targeted precision, e.g. GCC's software implemented float128:

>>> from sympy.printing.ccode import C99CodePrinter
>>> from sympy.codegen.ast import FloatType
>>> f128 = FloatType('_Float128', 128, nmant=112, nexp=15)
>>> p128 = C99CodePrinter(dict(
...     type_aliases={real: f128},
...     type_literal_suffixes={f128: 'Q'},
...     type_func_suffixes={f128: 'f128'},
...     type_math_macro_suffixes={
...         real: 'f128',
...         f128: 'f128'
...     },
...     type_macros={
...         f128: ('__STDC_WANT_IEC_60559_TYPES_EXT__',)
...     },
...     math_macros={}
... ))
>>> p128.doprint(tau**Rational(7, 2))
'powf128(tau, 7.0Q/2.0Q)'

For generating Python code there was previosuly one function (sympy.printing.python) which generated code dependent on SymPy. During the project a proper code printer for Python was introduced (an example of its output is shown later). The much used function lambdify was also changed to use this new printer. Introducing such a big change without breaking backward compatibility was certainly a challenge, but the benefit is that the user may now subclass the printers to override their default behaviour and use their custom printer in lambdify.


One usual challenge when working with symbolic expressions is that there are many ways to write the same expresisons. For code-generation purposes we want to write it in a manner which maximizes performance and minimizes significance loss (or let the user make that choice when the two are at odds). Since SymPy already has a great tools for traversing the expression tree and applying quite advanced pattern matching based replacements using Wild it was reasonably straightforward to implement rewriting rules for transforming e.g. 2**x to exp2(x) etc. Using the same structure, rules for rewriting expressions to drop small elements in sums (based on a user-predefined bounds).


One of the great benefitst from being able to represent abstract syntax trees as (largetly) language agnostic SymPy obejcts is that we can create functions for building these trees. Simpler numerical algorithms (which are ubiquitous in scientific codes) can be collected under sympy.codegen.algorithms. As a first case Newton's algortihm was implemented:

>>> from sympy import cos
>>> from sympy.codegen.algorithms import newtons_method_function
>>> ast = newtons_method_function(cos(x) - x**3, x)
>>> print(ccode(ast))
double newton(double x){
   double d_x = INFINITY;
   while (fabs(d_x) > 9.9999999999999998e-13) {
      d_x = (pow(x, 3) - cos(x))/(-3*pow(x, 2) - sin(x));
      x += d_x;
   return x;

once we have the AST we can print it using the python code printer as well:

>>> from sympy.printing import pycode
>>> print(pycode(ast))
def newton(x):
    d_x = float('inf')
    while abs(d_x) > 1.0e-12:
        d_x = (x**3 - math.cos(x))/(-3*x**2 - math.sin(x))
        x += d_x
    return x

or the Fortran code printer:

>>> from sympy.printing import fcode
>>> print(fcode(ast, source_format='free', standard=2003))
real*8 function newton(x)
real*8 :: x
real*8 :: d_x = (huge(0d0) + 1)
do while (abs(d_x) > 1.0d-12)
   d_x = (x**3 - cos(x))/(-3*x**2 - sin(x))
   x = x + d_x
end do
newton = x
end function

Newton's method is quite simple, but what makes SymPy suitable for this is that it needs the ratio between the function and its derivative.


I think that I managed to address all parts of my proposal. That being said, there is still a lot of potential to expand the sympy.codegen module. But now there are purposefully made base classes for creating AST node classes (sympy.codegen.ast.Token & sympy.codegen.ast.Node), the language agnostic ones are general enough that an algorithm represented as a single AST can be printed as Python/C/Fortran. At some level code will still be needed to be written manually (presumably as templates), but the amount of template rendering logic can be significantly reduced. Having algorithm AST factories such as the one for Newton's method in sympy.codegen.ast.algorithms is also exciting since those algorithms can be unit-tested as part of SymPy. Ideas for furthor work on code-generation with SymPy have been added to the list of potential ideas for next years GSoC.


I plan to continue to contribute to the SymPy project, and start using the new resources in my own research. Working with the new classes should also allow us to refine them if needed (preferably before the next release is tagged in order to avoid having to introduce deprecation cycles). SymPy is an amazing project with a great community. I'm really grateful to Google for funding me (and others) to do a full summers work on this project.

Status update week 13 GSoC

This was the last week of work on GSoC. I have been hard at work improving documentation and examples for the code.

Adding much needed documentaiton

I've spent the weekend adding examples and writing up documentation for my big PR #13100 which is not yet merged. I am quite excited how this PR turned out and I am happy with the design of the underlying AST nodes.


A new submodule .codegen.rewriting was added (in #13194), this allows a user to rewrite expressions using special math functions. The provided rules are those to rewrite to C99's special math functions (expm1, log1p etc.). I think it will be a useful addition (I have myself had the need for exactly this in my own research). The design is quite simple thanks to the excellt replace function in SymPy. There are still some corner cases (I have an "xfailed" test checked in for example).

Status update week 12 GSoC

Working on new AST nodes

#13100 is shaping up to be the largest PR of my GSoC project. The design of the new AST nodes especially (Token) is really helpful. But there is still a design issue: some nodes would naturally take different arguments depending on what language is being targeted. So I came to the conclusion that I needed some way of representing attributes. The solution I came up with would be to have a slightly more capable Node class (subclassing Token) which would in turn be subclassed from for nodes that need attributes.

I also enhanced the printing of both of these classes and introduced a String class, which in contrast to Symbol does not accept assumptions in its constructor, and does not have implied printing rules of sub- & superscript etc.

Adding .codegen.algorithms

A new submodule .codegen.algorithms was added, containing a AST generating function for Newton's method. This makes a nice design target for both the printers and AST nodes: being able to express the same AST in differnt languages is definitely an indication that we have a versatile printing system.

Compiling code on the fly

Introducing the .codegen.algorithms module also made the need to test generated code during CI runs clear. Jason Moore has previously mentioned that he thinks one of my python packages (pycompilation) would fit nicely into SymPy. I've been a bit relucatant to port it over since I have felt that it has not seen enough testing (and only under Linux). But now there was a need and we could start by making it an internal package only used by our own tests. That way it will get to mature without having to worry about deprecation cycles. And once more platforms are added to SymPy's CI configuration it would also see testing on other platforms (using AppVeyor for SymPy has been discussed for a long while now).

Status update week 11 GSoC

Checking in the new AST node types

#12693 got merged 🎉. It took a few rewrites essentially but I fell that the design of the new nodes will allow us to scale with reasonable maintance cost when adding new language specific nodes. The base class for new AST nodes (Token) to sublcass from allows one to implement nodes in an expressive manner by setting __slots__. The constructor of Token then sets the .args of Basic based on __slots__ this has the benefit that you need not write setters and getters using @property decorators (which quickly becomes tiresome when you have many classes).

Checking in the refactored PythonPrinter and changes to lambdify

Finally the challenging work of refactoring lambdify got merged into SymPy's master branch: #13046. We eventually decided to drop the contents of the old translations dictionaries but leave them be (in an empty state) in case users were modifying those in their code. Hopefully this approach doesn't break any code out there. Given how popular lambdify is among SymPy's users, it is a bit worrying that the test suite is not that extensive. I do remember a google engineer mentioning that the follow the "Beyonce principle": "I you liked it you should have put a test on it". Funny at is may be I hope I don't need to defend these changes with that arguement.

Status update week 10 GSoC

Refactoring lambdify

In my work to refactor lambdify I had come up with a solution where I would dynamically subclass the CodePrinters in lambdify to add translations from the old translation dictionaries. I was not happy with the solution and I don't think Aaron was either, we decided to keep the old import mechanism of lambdify which populated the namespace (instead of trying to generate code for the used imports which I had been trying).

So the work on refactoring lambdify has continued in #13046. And with this <> commit I hope we are close to getting the new version of lambdify out the door.

New FloatType representation

With some underlying assumptions about floating point representation (two's complement etc.) I have now a new representation of FloatType. I'm much happier with this representation and I think with it #12693 is much closer to getting merged.

Status update week 8 GSoC

I finished up the new PythonCodePrinter in #12808. And I have started work locally on having the NumPy- and Mpmath-printers subclass this new CodePrinter. Both of these printers were previously subclassing LambdaPrinter, so the change is a bit intrusive with lots of failing tests in the SymPy test suite. But this is a prerequisite for changing sympy.utilities.lambdify to accept custom CodePrinter subclasses as user input.

Status update week 9 GSoC

This week has been mostly spent on two refactoring lambdify and enhancing the representation of types in our AST nodes under .codegen.ast

Refactoring lambdify

The tricky business of chaning lambdify to use our new CodePrinter (or more importantly: letting the user pass their own subclasses continues).

There are currently dictionaries in sympy.utilities.lambdify which contains translations of SymPy entities to functions in corresponding targeted python module (e.g. NumPy functions if the numpy module is targeted). This week I've worked toward removing dependence on these, and instead rely on the user providing print methods instead.

Work on types in .codegen.ast

In my WIP PR for AST nodes I have made more type related changes. I've also changed the Fortran printer to use the new type information. I am not 100% satisfied with the way the floating point types are being represented at the moment (I'm using floating point numbers to store the max, min and epsilon values for example). I will probably rewrite those classes to use the underlying number of bits for mantissa and exponent as the canonical data for the class instances.

I finished up the new PythonCodePrinter in #12808. And I have started work locally on having the NumPy- and Mpmath-printers subclass this new CodePrinter. Both of these printers were previously subclassing LambdaPrinter, so the change is a bit intrusive with lots of failing tests in the SymPy test suite. But this is a prerequisite for changing sympy.utilities.lambdify to accept custom CodePrinter subclasses as user input.

Status update week 7 GSoC

This was an exciting week. I flew to Autstin, TX to attend the SciPy 2017 conference. On Monday morning Aaron Meurer, Jason Moore, Kenneth Lyons and I gave our new tutorial on code-generation. You can view the tutorial here:

You can even follow along the tutorial video in a jupyter notebook in your browser (thanks to the generous hosting by the folks at

It was nice to see that quite a few people were using SymPy for code-generation in the community. And coincidentally a proper Python code printer was requested by some users (this is on my todo-list for this summers project).

At the end of the conference there were also sprints, and SymPy did sprint as well. For my project, the greatest benefit from the sprints was that I got to discuss design choices wiht Aaron. We now have a plan for how to deal with types in the code-printers, so now it is "only" a matter of writing up the code to get it to work.

Status update week 6 GSoC

Finishing up tutorial material

We have been busy adding the final touches to the upcoming tutorial. (you can read about the tutorial in ealier blog posts).

I thought I would share a technique for doing exercises in Jupyter notebooks.

The %exercise magic

There are quite a few ways to provide exercises in jupyter notebooks (eventually I think there will be a de facto standard, either as a separate package or included with the notebook itself). But I essentially wanted a special version of %load which would load a file into a cell but replace some parts of the expressions with gaps for the tutorial attendee to fill in.

Thanks to IPython being a popular project there was already a solution posted on StackOverflow which did nearly exactly this.

After some tweaking I came up with this

import IPython.core.magic as ipym

class ExerciseMagic(ipym.Magics):

    def exercise(self, line, cell=None):
        token = '  # EXERCISE: '
        out_lst = []
        for ln in open(line).readlines():
            if token in ln:
                pre, post = ln.split(token)
                out_lst.append(pre.replace(post.rstrip('\n'), '???') + '\n')
        out_str = '# %exercise {0}\n{1}'.format(line, ''.join(out_lst)), replace=True)

def load_ipython_extension(ipython):

def unload_ipython_extension(ipython):

The downside of using this magic is that reloading the cell requires the user to uncomment the # %exercise ... line, delete the rest of the cell and rerun it. This will only feel natural if you've already worked with the %load magic before (which in my humble opinion is suboptimal since you often want bi-directional editing whenever you pull out that magic).

Status update week 5 GSoC

Continued work on tutorial material

Cross-platform work is always a bit tricky. Python does in general a great job of allowing programs to access the filesystem, network etc in an OS agnostic manner, but sometimes there is no way around having to write special code for each platform. In Python that often means that we will do runtime inspection:

>>> sys.platform
'linux'  # or e.g. darwin
'posix'  # or e.g. nt

To add another dimension of complexity we sometimes deal with 32-bit and sometimes 64-bit systems (most headaches are caused by pointer sizes, size of int and long), in python we can get the size of int of the platform quite easily:

>>> sys.maxsize - 2**63 + 1

when you, in addition to this, are doing code-generation, you are most likely relying on a compiler. That is our third dimension of complexity, on linux that (usually) translates to gcc & llvm/clang (but Intel's icc/ifort & portland group compilers and not uncommon). This third "dimension" has been the number one struggle the past week (and in particular the combination mingw/windows). At this point I feel that it is not worthwhile to continue pursuing mingw support on Windows.

It is quite likely that a few of our attendees will have some technical difficulties with system installed compilers. Therefore, I'm particularly happy that we now have got our binder environment to work flawlessly. It took some work to setup a Dockerfile, but now it feels reassuring to have this as a fall-back.

Work on SymPy for version 1.1

The work on the PythonCodePrinter has progressed since last week and should hopefully soon be ready for merging. The purpose of this class is to be stepping stone moving away from having e.g. lambdify relying on the StrPrinter, and instead use the CodePrinter class in .printing.codeprinter. Changing the super class showed to be trickier than first anticipiated (I thought it would be an afternoons work at most). I originally planned to rename the original PythonPrinter to something more appropriate (SymPyCodePrinter) together with its convenience functions. However, since this would incur a deprecation for only name change, it was decided that we will have to live with the old naming (to avoid raising tons of deprecation warnings in downstream projects).

We had identified two extensions we wanted to introduce to the codegeneration/autowrap factilities before the v1.1 release:

  • Kenneth Lyons (@ixjlyons) spear-headed the effort to allow custom CodeGen subclasses to be passed to codegen
  • ...and Jason the changes to the CythonWrapper to allow compiler arguments to be passed to autowrap.

Plans for the upcoming week

With a release candidate of SymPy v1.1 out the door (an impressive effort led by @asmeurer) and most cross-platform issues out of the way, it is about time for the final iterations of the tutorial material (and hopefully iron out any related issues in the rc).

At the end of the week it's time for me to fly to the US to attend SciPy 2017. I have high expectations and believe that it will be a great catalyst for the rest of the summers work on code-generation.