This is prometeo, an experimental modeling tool for embedded high-performance computing. prometeo provides a domain specific language (DSL) based on a subset of the Python language that allows one to conveniently write scientific computing programs in a high-level language (Python itself) that can be transpiled to high-performance self-contained C code easily deployable on embedded devices.

The same program runs in the Python interpreter, where it can be debugged, checked and plotted with the whole Python ecosystem, and as transpiled C, which computes the same numbers with static memory and no Python runtime:

_{examples/python_workflow: an SQP for closed-loop MPC of a cart-pole, run in Python (the Python-only code checks the results against SciPy and IPOPT and opens an interactive figure), then transpiled to C. The animation is generated from actual runs by gifs/make_workflow_gif.py.}

• Python compatible syntax: prometeo is a DSL embedded into the Python language. prometeo programs run in the Python interpreter, where they can be debugged with pdb and combined with any Python code (NumPy, SciPy, matplotlib, CasADi, ...) that is removed when transpiling.
• NumPy-like linear algebra: matrix and vector expressions such as P = Q + A.T @ P @ A - K.T @ R @ K, K[0:n, n:] = B.T, y = A @ x + B @ u, w = u * v or c : float = x @ Q @ x, with cholesky and cho_solve, are lowered to BLASFEO calls, with temporaries allocated statically.
• checked: a static checker accepts exactly the documented subset of Python and rejects everything else with an error pointing to the source, before running or transpiling. Types and matrix sizes are checked at compile time (and at run time when they are not known).
• efficient: prometeo programs transpile to C code that calls BLASFEO, with performance comparable to hand-written code using it (see performance).
• deterministic memory usage: prometeo computes a bound on the memory used by a program, and transpiled programs allocate it once: there is no dynamic memory allocation and no garbage collection.
• self-contained and embeddable: transpiled programs are C99 code that depends only on BLASFEO and a small runtime, not on the Python runtime library.
• CasADi: models and their derivatives defined with CasADi are called from prometeo code, through C code generated by CasADi.

examples/minimal/minimal.py computes an LQR controller for a double integrator and simulates the closed loop:

from prometeo import *

# LQR for a double integrator: Riccati recursion, then a closed-loop simulation

n : dims = 2        # sizes are compile-time constants
m : dims = 1
N : dims = 100

def lqr(A: pmat[n, n], B: pmat[n, m], Q: pmat[n, n], R: pmat[m, m],
        P: pmat[n, n], K: pmat[m, n]) -> None:
    # P = Q + A'P(A + BK), K = -(R + B'PB)^-1 B'PA, iterated N times
    L : pmat = pmat(m, m)
    pmat_copy(Q, P)
    for k in range(N):
        L = cholesky(R + B.T @ P @ B)
        K = -cho_solve(L, B.T @ P @ A)
        P = Q + A.T @ P @ (A + B @ K)
        P = 0.5 * (P + P.T)
    return

def main() -> int:
    A : pmat = pmat(n, n)
    B : pmat = pmat(n, m)
    Q : pmat = pmat(n, n)
    R : pmat = pmat(m, m)
    P : pmat = pmat(n, n)
    K : pmat = pmat(m, n)
    x : pvec = pvec(n)
    h : float = 0.1
    A[0, 0] = 1.0
    A[0, 1] = h
    A[1, 1] = 1.0
    B[1, 0] = h
    Q[0, 0] = 1.0
    Q[1, 1] = 0.1
    R[0, 0] = 0.01

    lqr(A, B, Q, R, P, K)
    pmat_print(K)

    x[0] = 1.0
    print('cost-to-go: %f\n' % (0.5 * (x @ P @ x)))
    for t in range(50):
        x = (A + B @ K) @ x
    print('state after 5 s: %e\n' % x[0])
    return 0

Run it in Python, or transpile it to C, compile it and run it:

$ pmt minimal.py
 -7.60447  -4.97765

cost-to-go: 3.272843

state after 5 s: 2.385667e-06

Execution time = 0.012 sec

$ pmt minimal.py --cgen=True
 > prometeo: starting transpilation
 ...
 > prometeo: running compiled C code:

 -7.60447  -4.97765

cost-to-go: 3.272843

state after 5 s: 2.385667e-06

execution time:0.000422s

Every matrix expression becomes a few BLASFEO calls on statically allocated temporaries; for example, the first two lines of the Riccati recursion are transpiled to:

    // L = cholesky(R + B.T @ P @ B)
    c_pmt_pmat_copy(R, ___pmt_m0);
    c_pmt_pmat_fill(___pmt_m1, 0.0);
    c_pmt_gemm_ab_tn(1.0, B, P, 1.0, ___pmt_m1, ___pmt_m1);
    c_pmt_gemm_ab_nn(1.0, ___pmt_m1, B, 1.0, ___pmt_m0, ___pmt_m0);
    c_pmt_pmat_fill(L, 0.0);
    c_pmt_potrf(___pmt_m0, L);
    // K = -cho_solve(L, B.T @ P @ A)
    c_pmt_pmat_copy(L, ___pmt_m0);
    c_pmt_pmat_fill(K, 0.0);
    c_pmt_pmat_fill(___pmt_m1, 0.0);
    c_pmt_gemm_ab_tn(1.0, B, P, 1.0, ___pmt_m1, ___pmt_m1);
    c_pmt_gemm_ab_nn(1.0, ___pmt_m1, A, 1.0, K, K);
    c_pmt_potrsm(___pmt_m0, K);
    c_pmt_pmat_scale(-1.0, K);

Sizes are declared with dims constants (and pmat[rows, cols] annotations for arguments), so mistakes are found before the program runs:

$ pmt minimal.py      # with a typo: K = -cho_solve(L, B.T @ P @ B)

 > prometeo: error at minimal.py:16:13: mismatching sizes: the result is 1 x 1, but K is 1 x 2
             K = -cho_solve(L, B.T @ P @ B)
                 ^

Code between # pure > and # pure < lines runs only in Python and is removed when transpiling: it can check intermediate results with NumPy (np.array(A) converts a matrix), compare with reference solvers, log and plot, or stop in the debugger. examples/python_workflow shows this workflow on an MPC example, including how its plots and checks found two bugs while it was written.

Since prometeo programs transpile to C code that calls the high-performance linear algebra library BLASFEO (publication: https://arxiv.org/abs/1704.02457, code: https://github.com/giaf/blasfeo), execution time can be comparable to hand-written high-performance code. The figure below shows a comparison of the CPU time necessary to carry out a Riccati factorization using highly optimized hand-written C code with calls to BLASFEO and the ones obtained with prometeo transpiled code from this example. The computation times obtained with NumPy and Julia are added too for comparison - notice however that these last two implementations of the Riccati factorization are not as easily embeddable as the C code generated by prometeo and the hand-coded C implementation. All the benchmarks have been run on a Dell XPS-9360 equipped with an i7-7560U CPU running at 2.30 GHz (to avoid frequency fluctuations due to thermal throttling).

The SQP examples compare with acados, whose solvers are hand-written C on top of BLASFEO, on the same optimal control problems (benchmarks/sqp_mpc). With the same BLASFEO build (AVX2), one real-time iteration of the prometeo SQP takes about as long as acados' (Intel i7-1370P, average over a closed loop of 30 steps):

Moreover, prometeo can largely outperform state-of-the-art Python compilers such as Nuitka. The table below shows the CPU times obtained on a Fibonacci benchmark.

prometeo needs Python 3.10 or later (it is tested with 3.10 to 3.14), a C compiler and make.

Releases are published on PyPI (pip install prometeo-dsl). The features described here may be more recent than the latest release: to use them, install from source:

git clone https://github.com/zanellia/prometeo.git
cd prometeo
git submodule update --init                     # BLASFEO
make -C prometeo/cpmt install_shared            # BLASFEO and the C runtime, into prometeo/lib
python3 -m venv .venv && source .venv/bin/activate
pip install -e .

Then run the examples in examples/ with pmt <example>.py (Python) or pmt <example>.py --cgen=True (C).

BLASFEO is built for X64_INTEL_HASWELL (x86-64 with AVX2) by default. For other CPUs, choose its target, e.g. make -C prometeo/cpmt install_shared TARGET=GENERIC (any CPU; the tests and CI use it) or TARGET=ARMV8A_ARM_CORTEX_A57; the targets are listed in external/blasfeo/Makefile.rule. On larger problems, the target matters: the AVX2 kernels are several times faster than the generic ones.

• the language: what can be written in prometeo programs, and what is checked; also on Read the Docs
• examples/python_workflow: developing in Python
• benchmarks/sqp_mpc: the comparison with acados

pip install -e '.[test]' && pytest tests runs, for every example, the transpiler against reference C code and compares the numbers printed by the Python and C versions; the conformance tests of the language (tests/conformance); random programs and random matrix expressions, run in Python and C and checked against NumPy; and the SQP examples against IPOPT.

Although translating a program written in a language into another with a comparable level of abstraction can be significantly easier than translating to one with a very different level of abstraction (especially if the target language is of much lower level), translating Python programs into C programs still involves a considerable abstraction gap it is not an easy task in general. Loosely speaking, the challenge lies in the necessity to reimplement features that are natively supported by the source language in the target language. In particular, when translating Python to C, the difficulty comes both from the different level of abstraction of the two languages and from the fact that the source and target language are of two very different types: Python is an interpreted, duck-typed and garbage-collected language and C is a compiled and statically typed language.

The task of transpiling Python to C becomes even more challenging if we add the constraint that the generated C code must be efficient (even for small to medium scale computations) and deployable on embedded hardware. In fact these two requirements directly imply that the generated code cannot make use of: i) sophisticated runtime libraries, e.g., the Python runtime library, which are generally not available on embedded hardware ii) dynamic memory allocation that would make the execution slow and unreliable (exception made for memory that is allocated in a setup phase and whose size is known a priori).

Since source-to-source code transformation, or transpilation, and in particular transpilation of Python code into C code is not an unexplored realm, in the following, we mention a few existing projects that address it. In doing so, we highlight where and how they do not satisfy one of the two requirements outlined above, namely (small scale) efficiency and embeddability.

Several software packages exist that address Python-to-C translation in various forms.

In the context of high-performance computing, Numba is a just-in-time compiler for numerical functions written in Python. As such, its aim is to convert properly annotated Python functions, not entire programs, into high-performance LLVM code such that their execution can be sped up. Numba uses an internal representation of the code to be translated and performs a (potentially partial) type inference on the variables involved in order to generate LLVM code that can be called either from Python or from C/C++. In some cases, namely the ones where a complete type inference can be carried out successfully, code that does not rely on the C API can be generated (using the nopython flag). However, the emitted LLVM code would still rely on Numpy for BLAS and LAPACK operations.

Nuitka is a source-to-source compiler that can translate every Python construct into C code that links against thelibpython library and it is therefore able to transpile a large class of Python programs. In order to do so, it relies on the fact that one of the most used implementations of the Python language, namely CPython, is written in C. In fact, Nuitka generates C code that contains calls to CPython that would normally be carried out by the Python parser. Despite its attractive and general transpilation approach, it cannot be easily deployed on embedded hardware due to its intrinsic dependency on libpython. At the same time, since it maps rather closely Python constructs to their CPython implementation, a number of performance issues can be expected when it comes to small to medium scale high-performance computing. This is particularly due to the fact that operations associated with, for example, type checking, memory allocation and garbage collection that can slow down the execution are carried out by the transpiled program too.

Cython is a programming language whose goal is to facilitate writing C extensions for the Python language. In particular, it can translate (optionally) statically typed Python-like code into C code that relies on CPython. Similarly to the considerations made for Nuitka, this makes it a powerful tool whenever it is possible to rely on libpython (and when its overhead is negligible, i.e., when dealing with sufficiently large scale computations), but not in the context of interest here.

Finally, although it does not use Python as source language, we should mention that Julia too is just-in-time (and partially ahead-of-time) compiled into LLVM code. The emitted LLVM code relies however on the Julia runtime library such that considerations similar to the one made for Cython and Nuitka apply.

Transpilation of programs written using a restricted subset of the Python language into C programs is carried out using prometeo's transpiler. This source-to-source transformation tool analyzes abstract syntax trees (AST) associated with the source files to be transpiled in order to emit high-performance and embeddable C code. In order to do so, special rules need to be imposed on the Python code. This makes the otherwise extremely challenging task of transpiling an interpreted high-level duck-typed language into a compiled low-level statically typed one possible. In doing so, we define what is sometimes referred to as an embedded DSL in the sense the resulting language uses the syntax of a host language (Python itself) and, in prometeo's case, it can also be executed by the standard Python interpreter.