Comparisons¶

The examples that measure themselves against other tools, on one page. Each section is copied from its example, with the table and what each port needed. The counterpart programs are on the example's page and in its compare/ folder. examples/compare.py checked that they all print the same answer before timing them.

Numerics: Compared with Numba, Codon, Mojo, C, and Rust¶

From Numerics. The same three functions, written for five other compilers, in compare/: semantics_numba.py, semantics_codon.py, semantics.mojo, semantics.c, and semantics.rs.

There is nothing to time here. The question is what each prints for may_overflow(30), floor_semantics(-7, 2), and modulo_semantics(7, -2). The answers Python gives are 265252859812191058636308480000000, -4, and -1.

	`may_overflow(30)`	`-7 // 2`	`7 % -2`	what the code says
Python, PPy on every path	265252859812191058636308480000000	-4	-1	`int` is an integer
Numba `@njit`	-8764578968847253504	-4	-1	64-bit wrap; Python's floor and sign
Mojo 1.0 `Int`	-8764578968847253504	-4	-1	64-bit wrap; Python's floor and sign
Codon	-8764578968847253504	-3	1	64-bit wrap; C's truncation and sign
C `long long` (gcc, -O0 and -O3)	-8764578968847253504	-3	1	signed overflow is undefined; truncation
Rust `i64`, release	-8764578968847253504	-3	1	wraps; truncation
Rust `i64`, debug	panics: attempt to multiply with overflow	-3	1	checked, then aborts

Numba and Mojo keep Python's division and remainder and wrap the multiply silently.
Codon, C, and Rust in release keep C's rounding and wrap.
Rust in debug is the only other one that notices the overflow, and its answer is to stop.

PPy notices and continues in Python: the same function, the same source, the number Python prints. --unsafe on ppy run or ppy build buys the wrap-semantics row, and says so.

Numba 0.67.0 on CPython 3.12.13, Codon 0.19.6, Mojo 1.0.0, gcc 13.3, rustc 1.95.0.

Value classes: Compared with Numba's `@jitclass`¶

From Value classes. A million distance2 calls from Python, and a native loop over a Ray's two fields, in compare/: vectors_bench.ppy (under ppy run and, the same file, under python) and vectors_numba.py. Milliseconds, best of five, over five processes.

PPy is a @dataclass and functions over it. Numba's nearest thing is a @jitclass with a typed field list, whose instances are Numba's own objects rather than Python ones:

@dataclass
class Ray:
    origin: float
    direction: float


@ppy.pure
@ppy.opt(3)
def travel(ray: Ray, count: int) -> float:
    position: float = ray.origin
    total: float = 0.0
    for i in range(count):
        position = position * 0.999999 + ray.direction * (i % 3)
        total += position
    return total

@jitclass([("origin", float64), ("direction", float64)])
class Ray:
    def __init__(self, origin, direction):
        self.origin = origin
        self.direction = direction


@njit
def travel(ray, count):
    position = ray.origin
    total = 0.0
    for i in range(count):
        position = position * 0.999999 + ray.direction * (i % 3)
        total += position
    return total

	PPy `ppy run`	CPython, the same file	Numba `@jitclass`
distance2, a million calls from Python	62.18 ± 0.47	63.13 ± 0.71	349.84 ± 7.47
travel, eight million steps over a Ray natively	11.16 ± 0.19	262.98 ± 2.15	10.97 ± 0.17

Inside a native loop the two are the same code: the Ray is two doubles in registers on both. The difference is the boundary.

Eight float operations are too little work to see past it. A PPy value class is still a Python dataclass: the generated boundary reads its fields and passes scalars, and a million calls cost what CPython takes to run the body itself. A @jitclass instance crosses into an @njit function through Numba's dispatcher and its own object layout, several times that per call.

The native loop is where the work is, and there the two agree. The loop here takes the class as a parameter.

Intel Core Ultra 9 386H; Numba 0.67.0 on CPython 3.12.13, PPy on CPython 3.14.5, from a checkout on a native filesystem.

Tuples: Compared with Numba¶

From Tuples. A million divmod_pair calls from Python, and a native loop that unpacks a tuple parameter and walks eight million steps, in compare/: pairs_bench.ppy (under ppy run and, the same file, under python) and pairs_numba.py. Milliseconds, best of five, over five processes.

Tuples are native in both. PPy spells the types; Numba infers them at the first call:

@ppy.pure
@ppy.opt(3)
def walk(start: tuple[float, float], count: int) -> float:
    x, y = start
    for i in range(count):
        x, y = (x + float(i)) / 2.0, (y + 1.0) / 2.0
    return x + y

@njit
def walk(start, count):
    x, y = start
    for i in range(count):
        x, y = (x + float(i)) / 2.0, (y + 1.0) / 2.0
    return x + y

	PPy `ppy run`	CPython, the same file	Numba `@njit`
divmod_pair, a million calls from Python	32.78 ± 1.40	47.86 ± 0.33	107.03 ± 5.77
walk, eight million steps from a tuple natively	11.26 ± 0.12	316.30 ± 2.14	11.25 ± 0.19

The loop is the same machine code on both: two doubles in registers, nothing allocated. The call from Python is where they differ, and the difference is the boundary. PPy's generated wrapper unpacks two ints and boxes a pair on the way out; Numba's dispatcher types the arguments on every call.

The loop here takes the tuple as a parameter and rebuilds it in place; see Limitations for why.

Intel Core Ultra 9 386H; Numba 0.67.0 on CPython 3.12.13, PPy on CPython 3.14.5, from a checkout on a native filesystem.

Buffers and JIT: Compared with Numba, Cython, NumPy, and C¶

From Buffers and JIT. The three kernels over 8,192 elements: a sum, a dot product in order and reassociated, and a modular digest. Each call is timed as the mean of 2,000 calls, best of five rounds, milliseconds per call over five processes. The programs are in compare/: kernels_bench.ppy, kernels_numba.py, kernels_cython.pyx, kernels_numpy.py, kernels.c. At this size a call is a few microseconds, so the table is as much about the call boundary as about the loop.

PPy is the loop as written, with a Buffer[float] borrowed from an array.array and @ppy.fastmath where reassociation is allowed:

@ppy.pure
@ppy.opt(3)
def dot(a: Buffer[float], b: Buffer[float]) -> float:
    result: float = 0.0
    for i in range(len(a)):
        result += a[i] * b[i]
    return result

Numba is the same loop under @njit, with no annotations, NumPy arrays, and fastmath=True for the relaxed dot. Its dispatcher types the arguments on every call:

@njit
def dot(a, b):
    result = 0.0
    for i in range(len(a)):
        result += a[i] * b[i]
    return result

Cython is typed memoryviews and cdef locals in a .pyx, built into an extension by cythonize -i. boundscheck=False is what makes it a C loop, and there is no relaxed dot to write:

cpdef double dot(double[::1] a, double[::1] b):
    cdef double result = 0.0
    cdef Py_ssize_t i
    for i in range(a.shape[0]):
        result += a[i] * b[i]
    return result

NumPy has no loop: np.dot is BLAS, np.sum is pairwise, and (counts % m).sum() builds a temporary.

C is the loop with no Python around it, timed inside the process: the floor.

	PPy `ppy run`	Numba `@njit`	Cython	NumPy	C (the loop alone)
total	0.0035 ± 0.0000	0.0036 ± 0.0000	0.0036 ± 0.0001	0.0021 ± 0.0001	0.0034 ± 0.0001
dot	0.0035 ± 0.0000	0.0036 ± 0.0000	0.0036 ± 0.0000	0.0012 ± 0.0001	0.0035 ± 0.0001
dot_relaxed	0.0011 ± 0.0001	0.0011 ± 0.0000	0.0036 ± 0.0001	0.0012 ± 0.0001	0.0034 ± 0.0000
digest	0.0072 ± 0.0000	0.0112 ± 0.0001	0.0078 ± 0.0001	0.0199 ± 0.0002	0.0045 ± 0.0001

The ordered sums sit on the C loop in PPy, Numba, and Cython alike, a microsecond of call above it. The reassociated dot is where @ppy.fastmath and Numba's fastmath=True let the vectorizer in, and Cython, with no such switch, stays scalar.

NumPy's dot is BLAS and the fastest row. Its digest pays for an 8,192-element temporary and is the slowest. The C digest is the loop's floor; what stands between it and the others is the price of a call through the interpreter.

Intel Core Ultra 9 386H; Numba 0.67.0, Cython 3.3.0, NumPy 2.5.3 on CPython 3.12.13, gcc 13.3, PPy on CPython 3.13.13, from a checkout on a native filesystem.

Algorithms: The eight kernels against C¶

From Algorithms. Same machine, same session, kernel wall time: mean ± standard deviation over 7 runs, each a fresh process. Every row prints the same answer in every column.

kernel	plain	`ppy run`	`ppy build --unsafe`	C (`gcc -O3`)	C (`clang -O3`)
sieve 2e6	173.1 ± 2.9	9.2 ± 0.5	8.8 ± 0.4	13.9 ± 0.6	13.8 ± 0.3
collatz 3e5	1111.0 ± 13.1	39.9 ± 0.8	30.3 ± 0.4	41.4 ± 0.5	32.2 ± 0.3
knapsack 400×2e4	423.1 ± 2.7	5.0 ± 0.1	3.6 ± 0.2	2.4 ± 0.1	2.1 ± 0.3
edit 2000×2000	471.4 ± 6.5	2.9 ± 0.1	2.9 ± 0.3	3.6 ± 0.1	3.5 ± 0.2
floyd 220	436.3 ± 5.3	3.4 ± 0.5	3.1 ± 0.2	5.0 ± 0.1	2.8 ± 0.0
matmul 220	473.9 ± 5.7	3.6 ± 0.4	3.5 ± 0.1	2.1 ± 0.1	3.4 ± 0.1
union-find 5e5	149.6 ± 3.8	3.2 ± 0.2	3.2 ± 0.4	3.6 ± 0.2	3.8 ± 0.3
fermat 6e4	23.9 ± 0.3	2.3 ± 0.1	2.6 ± 0.2	1.8 ± 0.1	1.5 ± 0.1

Bold is the fastest cell in the row. The native path is 9× to 160× faster than plain CPython on these eight, and all three paths print identical answers.

The two PPy columns differ in integer semantics:

ppy run keeps Python-integer semantics. Overflow is guarded and falls back to arbitrary precision.
ppy build --unsafe is the wrap-semantics artifact. That is where collatz picks up its remaining 10 ms and knapsack its 1.4 ms. The other kernels are already guard-free in the loop and do not move.

--host-cpu is within the noise on all eight.

Guard hoisting ([tool.ppy.llvm] safeguards) closed most of the old gaps. A multiplied index like i * n + k proves its extreme cases once in a guard block ahead of the loop. The body then runs plain mul nsw with no side exits, which lets LLVM strength-reduce and vectorize.

Against both compilers:

PPy wins sieve, edit distance, and union-find, and collatz once the artifact wraps.
floyd beats gcc and loses to clang by a third of a millisecond.
gcc and clang keep knapsack, matmul, and fermat, whose remaining guards live on data values no range can prove.

gcc 13.3 and clang 22.1, both -O3, no -march. Writes go through borrowed buffers, so the caller sees them.

The eight kernels against Numba, Mojo, and Codon¶

The same eight kernels are ported to three tools that compile Python-shaped code, in compare/: algorithms_numba.py, algorithms.mojo, algorithms_codon.py. The table shows kernel wall time inside the program, in milliseconds, mean and spread over seven fresh processes. Numba is timed after one warm call, so its compile is not in the number.

Here is Collatz as each of them spells it.

PPy is annotated Python, and ppy run keeps Python's integers: the multiply is overflow-checked and falls back to arbitrary precision.

def collatz_longest(limit: int) -> int:
    best: int = 0
    for start in range(1, limit):
        n: int = start
        steps: int = 0
        while n != 1:
            n = n // 2 if n % 2 == 0 else 3 * n + 1
            steps += 1
        best = max(best, steps)
    return best

Numba is the same function with @njit in place of the annotations and a NumPy array in place of array.array. Codon is the same source with List[int] annotations, run with codon run -release. Both wrap at 64 bits without a word.

@njit
def collatz_longest(limit):
    best = 0
    for start in range(1, limit):
        n = start
        ...

Mojo is typed by hand: Int for every index, Int64 for every element, var on every local. This is the plain List version without UnsafePointer:

def collatz_longest(limit: Int) -> Int64:
    var best: Int64 = 0
    for start in range(1, limit):
        var n = Int64(start)
        var steps: Int64 = 0
        while n != 1:
            n = n // 2 if n % 2 == 0 else 3 * n + 1
            steps += 1
        best = max(best, steps)
    return best

	`ppy run`	`ppy build --unsafe`	Numba `@njit`	Mojo	Codon
sieve 2e6	9.14 ± 0.67	9.04 ± 0.67	8.04 ± 0.25	11.42 ± 0.41	8.28 ± 0.38
collatz 3e5	41.56 ± 1.13	31.40 ± 0.55	33.14 ± 0.73	70.72 ± 0.97	33.38 ± 0.56
knapsack 400x2e4	5.08 ± 0.04	3.66 ± 0.09	4.48 ± 0.08	2.26 ± 0.11	5.90 ± 0.10
edit 2000x2000	3.02 ± 0.16	2.86 ± 0.05	1.76 ± 0.05	8.32 ± 0.20	2.46 ± 0.05
floyd 220	3.22 ± 0.04	3.10 ± 0.00	2.48 ± 0.08	3.15 ± 0.18	1.30 ± 0.34
matmul 220	3.64 ± 0.05	3.72 ± 0.08	3.48 ± 0.11	8.25 ± 0.16	3.90 ± 0.10
union-find 5e5	3.10 ± 0.07	3.18 ± 0.41	2.46 ± 0.05	6.68 ± 0.51	3.98 ± 0.32
fermat 6e4	2.36 ± 0.05	2.62 ± 0.04	1.74 ± 0.05	2.20 ± 0.03	1.52 ± 0.04

Read ppy build --unsafe against the other three: it is the wrap-semantics column. ppy run carries the overflow guards, which is where collatz pays its 10 ms.

Numba leads on six of the eight by ten to thirty percent. It uses the same LLVM, a for over a NumPy array, and no guard of any kind.
Mojo's plain List loop is the slowest on five kernels and the fastest on knapsack, where its bounds-free table indexing shows.

Numba 0.67.0 on CPython 3.12.13, Mojo 1.0.0 (-O3), Codon 0.19.6, PPy on CPython 3.13.13; Intel Core Ultra 9 386H.

The six problems¶

Each subfolder is one competitive-programming problem, the shapes a judge sets. It reads from standard input and answers on stdout, with its own input and output format, so nothing depends on an outside site being up. Nothing inside the programs is instrumented. The times below are wall time of the whole process, measured from outside, input and interpreter startup included. The C reference reads the same input with scanf.

	problem	plain	`ppy build --unsafe`	`--standalone --unsafe`	C (`gcc`)	C (`clang`)
15a	N-Queens	142.0 ms	48.7 ms	5.7 ms	4.8 ms	5.4 ms
15b	shortest path	4858.5 ms	4426.9 ms	111.7 ms	142.2 ms	134.3 ms
15c	substring search	292.4 ms	53.3 ms	—	9.3 ms	9.0 ms
15d	range sums	1571.7 ms	1524.4 ms	30.1 ms	50.8 ms	50.9 ms
15e	longest increasing subsequence	541.6 ms	111.5 ms	41.5 ms	57.7 ms	53.9 ms
15f	counting inversions	603.7 ms	98.4 ms	38.6 ms	43.5 ms	44.2 ms

Every cell is the mean of five runs, recorded in measurements.json with the machine it was measured on. bench.py reproduces it, and scripts/refresh.py says when a number here has drifted. Both fail if the paths stop agreeing on the answer. Bold is faster than both C references.

ppy run is left out because it compiles before it runs. That adds a flat two seconds or so to every row, and it is the development path rather than the one to submit.

`ppy build --unsafe` and `--standalone --unsafe`¶

ppy build --unsafe is the hybrid. The kernels are native, but the glue around them (main, the buffers, print of a Python int) is the optimized Python the build wrote. So the binary embeds an interpreter and imports the runtime before the program begins. That is ~35 ms, and on the smaller problems it is most of what separates the column from C. Where main reads its input line by line (1.2 million edge lines in 15b, 400 thousand command lines in 15d), each ppy.input[...]() is a call into the runtime from the interpreter, about a microsecond apiece. Those reads are most of the plain and ppy build --unsafe columns there.
--standalone is a binary with no CPython in it at all. ldd shows libc and nothing else, and ppy.scan[int]() lowers to the same buffered scan of standard input that scanf does. The whole N-Queens executable is 17.0 KB against the C one's 16.1 KB.

Four of the five standalone rows beat both C references for one reason: ppy.scan reads into memory faster than scanf parses. N-Queens reads a single integer and then computes, so there is nothing for a faster reader to win back, and it stays behind.

Reading input¶

n = ppy.input[int]()                 # one line, as int(input()) reads it
a, b = ppy.input[tuple[int, int]]()  # one line, exactly two fields
row = ppy.input[Buffer[int]]()       # one line of integers, into a buffer
values = ppy.scan[Buffer[int]](n)    # n integer tokens, whatever lines they are on

ppy.input[T]() reads one line the way the builtin input() does and types the result from T. ppy.scan[T]() reads tokens across lines. Both buffer reads go straight into memory rather than building a Python object per field: 12.6 ms for 500k integers against 49.8 ms for sys.stdin.read().split() and 20.6 ms for C's scanf.

The conversion writes the line reads for you, and each reads the line the original read:

int(input()) becomes ppy.input[int]().
a, b = map(int, input().split()) becomes the tuple read.
array.array("q", map(int, input().split())) becomes the buffer line read.

A module that also touches sys.stdin keeps the input it has, since the typed reader owns the file descriptor.

Parallel ranges: Compared with Numba, Taichi, Mojo, and NumPy¶

From Parallel ranges. The same four kernels over eight million elements, each written the way its tool wants it, in compare/: ranges_bench.ppy (this example's kernels, with a relaxed dot beside the ordered one), ranges_numba.py, ranges_taichi.py, ranges.mojo, ranges_numpy.py. Milliseconds, best of five warm calls, over five processes.

PPy: parallel.range, typed parameters, the same file on CPython. dot keeps a serial sum's order by design, and @ppy.fastmath is the relaxed one:

def dot(a: Buffer[float], b: Buffer[float]) -> float:
    total = 0.0
    for i in parallel.range(len(a)):
        total += a[i] * b[i]
    return total

Numba: @njit(parallel=True) and prange. A prange sum reassociates on its own, so there is no ordered dot to write, and integers wrap silently:

@njit(parallel=True)
def dot(a, b):
    total = 0.0
    for i in prange(len(a)):
        total += a[i] * b[i]
    return total

Taichi: arrays are ti.fields declared at module level, and the outermost loop of a kernel is parallel by itself. The defaults overflow in 32 bits and accumulate in f32 silently; ti.init(default_ip=ti.i64, default_fp=ti.f64) is what makes the answers match:

@ti.kernel
def dot() -> ti.f64:
    total = 0.0
    for i in range(N):
        total += a_f[i] * b_f[i]
    return total

Mojo: a different language. It uses List[Float64], mut parameters, Int64(i) where an index meets an element, and parallelize from max.algorithm taking a function of one index, so a reduction is partial sums per worker written by hand.

NumPy has no loop to write and no way to keep one: dot is BLAS, and count_odd builds two temporaries to count.

	PPy `ppy run`	Numba `prange`	Taichi	Mojo	NumPy
squares	1.34 ± 0.03	0.93 ± 0.06	1.16 ± 0.04	1.05 ± 0.05	13.82 ± 0.52
dot, in order	6.39 ± 0.40	—	—	14.93 ± 1.30	—
dot, reassociated	1.80 ± 0.13	1.34 ± 0.04	1.50 ± 0.04	1.52 ± 0.10	1.33 ± 0.05
count_odd	0.99 ± 0.04	0.70 ± 0.04	1.23 ± 0.06	1.37 ± 0.19	27.64 ± 1.35
fill, serial	6.17 ± 0.11	7.51 ± 0.04	4.02 ± 0.38	19.57 ± 1.09	15.43 ± 0.28

On the elementwise kernels the four compiled tools are within a few tenths of a millisecond, which is the memory bandwidth of the machine. NumPy's temporaries are ten times behind.

The reductions are where the models differ. PPy's ordered dot is the only one that promises a serial sum's answer, and it is a serial sum. Its relaxed dot, Numba's and Taichi's are the same vectorized tree.

The serial fill is the row that measures a plain loop with no threads, and Mojo's plain List loop is the slow one there.

Intel Core Ultra 9 386H (16 threads), threads backend; Numba 0.67.0, Taichi 1.7.4, NumPy 2.5.3 on CPython 3.12.13, Mojo 1.0.0, PPy on CPython 3.13.13, from a checkout on a native filesystem.

Derivatives: Compared with JAX and PyTorch¶

From Derivatives. Newton's method on g from 100,000 starting points near 0.4, six steps each, with the derivative from each framework's own autodiff. The programs are in compare/: gradients_bench.ppy, gradients_jax.py, gradients_torch.py. Milliseconds for the whole batch, best of five, over five processes.

PPy: ppy.grad(g) is a function. newton calls it in a loop, and a loop over the starting points calls newton. Everything is scalar and native; nothing is batched:

dg = ppy.grad(g)


def newton(x: float) -> float:
    for _ in range(6):
        x = x - g(x) / dg(x)
    return x


def newton_all(count: int) -> float:
    total = 0.0
    for i in range(count):
        total += newton(0.4 + i * 1e-6)
    return total / count

JAX: jax.grad over jnp functions. To run 100,000 solves it wants them batched: vmap under jit, the loop as lax.fori_loop so it traces, and jax_enable_x64 so the digits match:

def newton(x):
    def step(_, x):
        return x - g(x) / dg(x)

    return jax.lax.fori_loop(0, 6, step, x)


newton_all = jax.jit(jax.vmap(newton))

PyTorch: torch.func.grad and vmap. The six steps stay a Python loop over a 100,000-element tensor in float64. The thread count is pinned to the performance cores', since its default of one per logical core spins on this hybrid CPU and the batch takes hundreds of milliseconds some runs:

dg = grad(g)


def newton(x):
    for _ in range(6):
        x = x - g(x) / dg(x)
    return x


newton_all = vmap(newton)

	PPy	JAX `vmap` + `jit`	PyTorch `torch.func`
newton, 100k starts	11.89 ± 0.07	1.23 ± 0.10	6.52 ± 1.04

The derivative is the same nine digits in all three: reverse mode over the same rule table. What differs is the shape of the program.

JAX and PyTorch are fast here because the problem batches, and XLA and ATen vectorize across the batch. A scalar newton in a Python loop would cost them tens of microseconds per call. PPy's scalar loop pays nothing per call and is not vectorized across the batch.

Which one is faster depends on whether your problem comes as 100,000 independent solves or as one.

Intel Core Ultra 9 386H; JAX 0.11.1 and PyTorch 2.14.0 (CPU) on CPython 3.13.13, PPy on CPython 3.13.13, from a checkout on a native filesystem.

Coroutines: Compared with asyncio and uvloop¶

From Coroutines. Twenty thousand eight-byte round trips over one loopback connection (a client and an echo server in one process), in compare/: echo_bench.ppy, echo_asyncio.py, echo_uvloop.py. Milliseconds for the whole run, over five processes.

PPy is the coroutines as the example writes them, lowered to state machines over the native loop.
asyncio is the standard library's streams, as an echo server is usually written.
uvloop is the same program on the libuv-backed loop.

async def serve(listening: int, rounds: int) -> int:
    client = await aio.accept(listening)
    room = native.stack_alloc[ppy.u8](8)
    total = 0
    for _ in range(rounds):
        got = await aio.read(client, room, 8)
        total += await aio.write(client, room, got)
    aio.close(client)
    return total

async def serve(reader, writer):
    total = 0
    for _ in range(ROUNDS):
        data = await reader.readexactly(8)
        writer.write(data)
        await writer.drain()
        total += len(data)
    writer.close()
    await writer.wait_closed()
    return total

	PPy `ppy.aio`	asyncio streams	uvloop
20000 round trips	102.52 ± 2.81	682.49 ± 26.46	153.19 ± 3.17

A round trip is four socket operations and four resumptions of a coroutine.

Under asyncio each of those is a Python frame, a Future, and a trip through the event loop's callback queue.
uvloop moves the loop and the transports into C and leaves the coroutines in Python.
PPy's state machines resume in native code and call read and write directly, so what remains per round trip is the syscalls.

The program pays for that in vocabulary. Sockets are integers, bytes go through pointers, and an error is a negative errno, where asyncio's streams give you bytes and exceptions.

Intel Core Ultra 9 386H; uvloop 0.22.1 on CPython 3.12.13, asyncio and PPy on CPython 3.14.5, from a checkout on a native filesystem.

Generics: Compared with Numba¶

From Generics. sweep over eight million values, in compare/:

generic_bench.ppy, under ppy run and, the same file, under python
generic_numba.py

Times are milliseconds, best of five, over five processes.

PPy is the generics as Python 3.12 spells them, with bounds, called from a plain function. Numba is the same three functions under @njit, generic by dispatch: each is compiled once per tuple of argument types it meets. That is the same instance-per-type-tuple rule without the declaration.

def largest[T: int | float](a: T, b: T) -> T:
    return a if a > b else b


def clamp[T: int | float](x: T, lo: T, hi: T) -> T:
    return largest(lo, x) if x < hi else hi

@njit
def largest(a, b):
    return a if a > b else b


@njit
def clamp(x, lo, hi):
    return largest(lo, x) if x < hi else hi

	PPy `ppy run`	CPython, the same file	Numba `@njit`
sweep, eight million clamps and comparisons	5.87 ± 0.03	428.92 ± 9.67	5.09 ± 0.06

Both compile largest twice (once for ints, once for floats) and call the instances directly from the loop, and the loop is the same code either way.

PPy adds a check at the source. The bound is written, so largest("a", 1) is refused before anything runs, and the file is still the file python runs.
Numba asks you to write nothing. The types are whatever arrives first.

Intel Core Ultra 9 386H; Numba 0.67.0 on CPython 3.12.13, PPy on CPython 3.14.5, from a checkout on a native filesystem.

Regular expressions: Compared with CPython's `re` and Rust's `regex`¶

From Regular expressions. The benchmark makes four passes over 400,000 generated lines (7.9 MB): count the words, find the longest, sum the name = value pairs, count the hex numbers. The programs are in compare/: patterns_bench.ppy, patterns_re.py, and patterns.rs with its Cargo.toml. The site shows them whole. Times are in milliseconds, best of five passes, mean and spread over five processes.

PPy: a pattern compiled from a bytes literal, search from a position, m.end(). The same file runs on CPython.

def sum_values(text: Buffer[ppy.u8]) -> int:
    total = 0
    pos = 0
    while True:
        m = PAIR.search(text, pos)
        if m is None:
            return total
        value = 0
        for i in range(m.start(2), m.end(2)):
            value = value * 10 + (text[i] - 48)
        total += value
        pos = m.end()

CPython re: the idiomatic spelling, finditer and int(m.group(2)).

def sum_values(text: bytes) -> int:
    return sum(int(m.group(2)) for m in PAIR.finditer(text))

Rust regex: captures_iter over &[u8], the number parsed from the captured bytes, built with cargo build --release.

fn sum_values(pair: &Regex, text: &[u8]) -> u64 {
    pair.captures_iter(text)
        .map(|c| std::str::from_utf8(&c[2]).unwrap().parse::<u64>().unwrap())
        .sum()
}

	PPy `ppy run`	CPython `re`	Rust `regex`
count_words	6.57 ± 0.05	107.07 ± 1.27	18.23 ± 0.13
longest_word	6.53 ± 0.05	133.05 ± 8.69	36.85 ± 2.99
sum_values	56.70 ± 0.33	198.20 ± 1.86	65.62 ± 0.79
count_hex	6.59 ± 0.03	10.22 ± 0.10	2.33 ± 0.10

PPy compiles each pattern into a function of core operations: a backtracker in which a byte-class repeat scans its run instead of stepping byte by byte. A match is a few locals rather than an object, so a search loop is a native loop. That is the whole distance to re, whose matcher is a bytecode interpreter allocating a match object per hit.

Against Rust, the pattern decides:

[A-Za-z]+ is a scan in both engines.
\w+\s*=\s*(\d+) makes regex fall back to its slower capturing engine.
0x[0-9a-f]+ under (?i) is where its lazy DFA wins outright. A literal-prefix prefilter is the obvious next step for the PPy matcher there.

Intel Core Ultra 9 386H; rustc 1.95.0 with regex 1.x, CPython 3.12.13 for re, PPy on CPython 3.13.13, from a checkout on a native filesystem.

GPU kernels: Compared with CuPy, Numba, Mojo, and CUDA C¶

From GPU kernels. The same two kernels over sixteen million doubles, written as thread-level kernels the way each tool spells them, in compare/: saxpy_bench.ppy, saxpy_cupy.py, saxpy_numba.py, saxpy.mojo, saxpy.cu. Times are milliseconds, best of warm launches with the device synchronized, over five processes.

Tile-level tools such as Triton and Taichi program a block as one vector and never write the shared-memory exchange, which makes it a different kernel. They are compared with PPy's tile kernels instead of these.

How each tool writes the block max¶

PPy uses a Python function with cuda.global_id(), shared memory, a barrier, and a shuffle. The same file runs on CPython through the reference launch.

@cuda.kernel
def block_max(x: native.const_ptr[float], out: native.ptr[float]) -> None:
    parked = cuda.shared[float, 64]()
    tid = cuda.thread_id()
    native.store(native.offset(parked, tid), native.load(native.offset(x, cuda.global_id())))
    cuda.syncthreads()
    mine = native.load(native.offset(parked, tid))
    other = cuda.shfl_xor(mine, 1)
    ...

Numba CUDA reads the same way, with @cuda.jit, cuda.grid(1), cuda.shared.array, cuda.syncthreads(), and cuda.shfl_xor_sync:

@cuda.jit
def block_max(x, out):
    parked = cuda.shared.array(64, dtype=np.float64)
    tid = cuda.threadIdx.x
    parked[tid] = x[cuda.grid(1)]
    cuda.syncthreads()
    mine = parked[tid]
    other = cuda.shfl_xor_sync(0xFFFFFFFF, mine, 1)
    ...

CuPy writes saxpy in one line, an ElementwiseKernel, and the block max as CUDA C in a string handed to RawKernel.

CUDA C is the kernel the others approximate, timed with events around the launch alone.

Mojo writes the kernel as a def with global_idx, thread_idx, a stack_allocation in AddressSpace.SHARED, barrier() from MAX's max.gpu.sync, and shuffle_xor. Some details differ from the others:

shuffle_xor has no Float64 form, so the double crosses as its bits.
Arguments must be fixed-width (Int64, not Int).
out is a reserved parameter name.
The launch is DeviceContext.enqueue_function with grid_dim and block_dim.

def block_max(x: UnsafePointer[Float64, MutAnyOrigin], result: UnsafePointer[Float64, MutAnyOrigin]):
    var parked = stack_allocation[64, Float64, address_space = AddressSpace.SHARED]()
    var tid = Int(thread_idx.x)
    parked[tid] = x[Int(global_idx.x)]
    barrier()
    var mine = parked[tid]
    var other = bitcast[DType.float64, 1](shuffle_xor(mine.to_bits[DType.uint64](), 1))
    ...

Results¶

	PPy `cuda.launch`	CuPy	Numba CUDA	Mojo	CUDA C
saxpy, arrays on the device	0.77 ± 0.08	0.78 ± 0.08	0.72 ± 0.06	0.55 ± 0.05	0.50 ± 0.00
block max, arrays on the device	2.03 ± 0.05	1.98 ± 0.02	2.02 ± 0.03	1.98 ± 0.02	1.88 ± 0.00
saxpy, arrays copied in and out per launch	39.19 ± 1.60	48.92 ± 2.05	36.79 ± 1.50	28.43 ± 1.45	30.46 ± 0.74
block max, array copied in per launch	15.66 ± 0.43	12.33 ± 0.34	14.11 ± 0.42	11.21 ± 0.33	11.26 ± 0.20

With the arrays on the device, a launch is the kernel. The five tools run the same block max in the same time. On saxpy the Python-hosted ones sit a tenth of a millisecond above Mojo and CUDA C, the cost of a launch through the interpreter.

The copying rows use the other memory model: a native.stack_alloc array sent in and brought back on every launch. There the driver reads the host array in place, and the traffic costs what it costs every tool. A program that launches more than once keeps its data on the device by allocating it there, with cuda.device_alloc.

NVIDIA GeForce RTX 5080 Laptop GPU, driver 610.71, CUDA 13.3; CuPy 14.2.0, Numba 0.67.0 on CPython 3.12.13; Mojo 1.0.0 with MAX 26.5; nvcc 13.3; PPy on CPython 3.13.13.

Tile kernels: Compared with Triton and Taichi¶

From Tile kernels. The same two kernels over sixteen million doubles, written as programs over tiles in the three tools that have them. They are in compare/: tiles_bench.ppy, tiles_triton.py, tiles_taichi.py. Times are in milliseconds, best of warm launches with the device synchronized, over five processes. The thread-level ports of the same kernels (CuPy, Numba, Mojo, CUDA C) are in the CUDA example's comparison.

Here is the block max as each tool spells it.

PPy: a program loads its tile and reduces it. The same file runs on CPython.

@tile.kernel
def block_max(x: native.const_ptr[float], out: native.ptr[float]) -> None:
    pid = tile.program_id()
    values = tile.load(x, pid * 64 + tile.arange(64))
    tile.store(out, pid, tile.max(values))

Triton: the same shape, with tl.program_id, tl.arange, tl.load, tl.max. It is launched over CuPy memory through a six-line data_ptr() wrapper, and its compiler tiles the work across a warp group.

@triton.jit
def block_max(x_ptr, out_ptr, BLOCK: tl.constexpr):
    offsets = tl.program_id(0) * BLOCK + tl.arange(0, BLOCK)
    tl.store(out_ptr + tl.program_id(0), tl.max(tl.load(x_ptr + offsets), axis=0))

Taichi: no program either. Arrays are ti.fields, the outer loop of a kernel is parallel over the blocks, and the max of each block is an inner serial loop. ti.init(arch=ti.cuda, default_fp=ti.f64) makes the answers match, and on WSL the driver goes on the library path by hand.

@ti.kernel
def block_max():
    for b in out:
        best = values[b * 64]
        for k in range(1, 64):
            candidate = values[b * 64 + k]
            best = candidate if candidate > best else best
        out[b] = best

	PPy `tile.launch`	Triton	Taichi
saxpy, arrays on the device	0.72 ± 0.09	0.77 ± 0.10	0.89 ± 0.11
block max, arrays on the device	0.41 ± 0.02	0.81 ± 0.02	0.31 ± 0.03
saxpy, arrays copied in and out per launch	41.11 ± 0.67	51.73 ± 6.92	163.66 ± 7.49
block max, array copied in per launch	14.87 ± 0.40	11.38 ± 0.36	29.90 ± 2.21

On saxpy, the three run at the memory bandwidth of the device.

On the block max they are three lowerings of one idea:

PPy gives the 64-lane tile a block of 64 threads and reduces it with a shuffle tree.
Triton gives it a warp group and its own reduction.
Taichi gives each block one thread that loops over the 64, with no exchange at all. That is why it is the fastest here, and why it would not be on a wider tile.

The copying rows are the other memory model: an array sent in and brought back on every launch. PPy's driver reads the host array in place, Triton goes through CuPy's asarray, and Taichi's from_numpy copies twice.

NVIDIA GeForce RTX 5080 Laptop GPU, driver 610.71, CUDA 13.3; Triton 3.8.0, Taichi 1.7.4 on CPython 3.12.13; PPy on CPython 3.13.13.

NumPy fusion: Compared with NumPy, numexpr, Numba, and JAX¶

From NumPy fusion. The two expressions over eight million doubles, in compare/: fusion_bench.ppy, fusion_numpy.py, fusion_numexpr.py, fusion_numba.py, fusion_jax.py. Milliseconds, best of five calls, over five processes.

PPy is the NumPy expression as written, in a function:

@ppy.pure
def blend(a: np.ndarray, b: np.ndarray) -> np.ndarray:
    return np.sin(a) * 2.0 + np.cos(b)

NumPy is the same line with no function around it: three calls, two 64 MB temporaries.

numexpr takes the expression as a string and evaluates it in chunks across threads:

ne.evaluate("sin(a) * 2.0 + cos(b)", local_dict={"a": a, "b": b})

Numba is an explicit loop under @njit writing into np.empty_like:

@njit
def blend(a, b):
    out = np.empty_like(a)
    for i in range(a.shape[0]):
        out[i] = math.sin(a[i]) * 2.0 + math.cos(b[i])
    return out

JAX is the expression under jax.jit on the CPU, with jax_enable_x64:

@jax.jit
def blend(a, b):
    return jnp.sin(a) * 2.0 + jnp.cos(b)

	PPy	NumPy	numexpr	Numba `@njit`	JAX `jit`
normalize	22.25 ± 0.49	21.68 ± 0.54	13.03 ± 0.28	25.63 ± 0.79	8.03 ± 0.14
blend	73.52 ± 0.72	88.02 ± 0.73	10.07 ± 0.20	84.53 ± 0.84	22.84 ± 0.22

blend is two transcendentals per element, and on one thread that is what the time is. PPy, NumPy, and Numba each call sin and cos once per element, and fusing away NumPy's two temporaries takes off the fifth of the time that was memory traffic. numexpr and JAX evaluate sin and cos across vector lanes and across threads, which is where their rows come from.

normalize is a reduction in NumPy's own order followed by a division pass, on every tool that keeps the order. PPy's sum is NumPy's sum, so the row is NumPy's time. @ppy.fastmath is the permission to reassociate it (Parallel shows what that buys).

Intel Core Ultra 9 386H (16 threads); Numba 0.67.0, numexpr 2.14.2, NumPy 2.5.3 on CPython 3.12.13, JAX 0.11.1 on CPython 3.13.13, PPy on CPython 3.14.5, from a checkout on a native filesystem.

Parallel fused kernels: Compared with NumPy, numexpr, Numba, and JAX¶

From Parallel fused kernels. The fused expression and the sum of squares over eight million doubles, in compare/: fused_bench.ppy, fused_numpy.py, fused_numexpr.py, fused_numba.py, fused_jax.py. Milliseconds, best of five calls, over five processes.

PPy is the NumPy expression in a function, with @ppy.parallel to split it. The serial row is the same expression fused into one loop with no decorator, and strict_total is NumPy's sum in NumPy's order:

@ppy.pure
@ppy.parallel
@ppy.opt(3)
def parallel(a: np.ndarray, b: np.ndarray) -> np.ndarray:
    return (a * b + a) * (b - a) + a * 0.5 - b * 0.25

NumPy is the expression as written: five operations, four 64 MB temporaries, one thread.

numexpr is the expression as a string, compiled to its own virtual machine and evaluated in chunks across threads:

ne.evaluate("(x * y + x) * (y - x) + x * 0.5 - y * 0.25", local_dict={"x": x, "y": y})

Numba is an explicit loop under @njit(parallel=True) writing into np.empty_like:

@njit(parallel=True)
def fused(a, b):
    out = np.empty_like(a)
    for i in prange(a.shape[0]):
        out[i] = (a[i] * b[i] + a[i]) * (b[i] - a[i]) + a[i] * 0.5 - b[i] * 0.25
    return out

JAX is the expression under jax.jit on the CPU, fused by XLA and run on its thread pool, with jax_enable_x64:

@jax.jit
def fused(a, b):
    return (a * b + a) * (b - a) + a * 0.5 - b * 0.25

	PPy	NumPy	numexpr	Numba `prange`	JAX `jit`
fused, serial	12.09 ± 0.41	—	—	—	—
fused	4.10 ± 0.18	56.26 ± 0.39	6.80 ± 0.31	3.37 ± 0.41	5.91 ± 0.46
sum of squares	13.35 ± 0.50	12.35 ± 0.27	6.90 ± 0.13	3.99 ± 0.08	0.65 ± 0.01
sum of squares, relaxed	4.47 ± 0.21	—	—	—	—

Fusing the expression is what removes NumPy's four temporaries. The serial fused loop reads the two inputs once and writes the output once. Splitting that loop across the cores is a memory-bandwidth problem that PPy, Numba, JAX, and numexpr solve the same way, within a couple of milliseconds of each other.

The fused kernel also checks its own result in the same loop: one add-reduction of non-finite elements the vectorizer keeps in a register, and one guard after. That is how it keeps NumPy's floating-point reporting without a second pass over 64 MB.

The ordered sum is NumPy's order and NumPy's time. The relaxed one vectorizes on one thread. JAX's reduction, reassociated across its pool, is the row that shows what the same permission buys with threads.

Intel Core Ultra 9 386H (16 threads), threads backend; NumPy 2.5.3, numexpr 2.14.2, Numba 0.67.0 on CPython 3.12.13, JAX 0.11.1 on CPython 3.13.13, PPy on CPython 3.13.13, from a checkout on a native filesystem.

PyTorch ATen regions: Compared with PyTorch eager and `torch.compile`¶

From PyTorch ATen regions. The same layer on an 8×32 input, twenty thousand calls, in compare/: layer_bench.ppy, layer_eager.py, layer_compile.py. Milliseconds per call, to four places, best of five rounds, over five processes; one PyTorch thread.

PPy is the function as written, and ppy run compiles it into one ATen region.
PyTorch eager is the same function with no decorator: three dispatches from Python.
torch.compile is the same function under the decorator, traced by Dynamo and written by Inductor.

@ppy.opt(3)
def layer(x: torch.Tensor, weight: torch.Tensor, bias: torch.Tensor) -> torch.Tensor:
    return torch.relu(torch.add(torch.matmul(x, weight), bias))

@torch.compile
def layer(x, weight, bias):
    return torch.relu(torch.add(torch.matmul(x, weight), bias))

	PPy ATen region	PyTorch eager	`torch.compile`
layer, per call	0.0019 ± 0.0000	0.0020 ± 0.0000	0.0087 ± 0.0002

Three operators on a tensor this small cost about two microseconds either way. PyTorch's eager dispatch is cheap enough that the Python round trips the region removes are within the noise of the ATen calls themselves. The region's value is not this number. It is that the function keeps its source, its autograd, and its dispatcher, and that a built artifact carries it with no compiler in the process. torch.compile pays for its guards on every call, which on an 8×32 input is more than the work.

One layer of one shape is not the question anyone asks about torch.compile, and this table cannot answer it: the work is two microseconds, so every column is measuring its own overhead. 46_gpt2 asks it at model scale instead: GPT-2 XL, forty-eight blocks, one region each, on a datacenter GPU.

Intel Core Ultra 9 386H; PyTorch 2.14.0 (CPU) on CPython 3.13.13, PPy against the same PyTorch on CPython 3.14.5, from a checkout on a native filesystem.

GPT-2 XL: one region per block: Compared with PyTorch eager and `torch.compile` at model scale¶

From GPT-2 XL: one region per block.

	PPy ATen regions	PyTorch eager	`torch.compile`	`torch.compile` reduce-overhead
forward b1 s32	8.78 ± 0.12	10.77 ± 0.08	7.46 ± 0.28	13.36 ± 0.02
forward b1 s128	9.24 ± 0.18	11.06 ± 0.08	7.88 ± 0.28	14.11 ± 0.01
forward b1 s512	14.87 ± 0.01	14.90 ± 0.01	14.74 ± 0.07	22.86 ± 0.02
forward b8 s512	95.33 ± 0.30	95.50 ± 0.18	95.55 ± 0.26	102.96 ± 0.66
decode 128 tokens	1204.24 ± 28.34	1484.26 ± 17.63	1057.01 ± 41.66	1658.16 ± 1.84
training step b4 s512	136.62 ± 0.23	136.73 ± 0.21	135.58 ± 0.25	143.53 ± 0.67

Milliseconds per pass of the whole model, bfloat16, best of ten in each of five processes. The decode row is a 128-token prompt and then 128 tokens one at a time out of the cache, timed without the prefill.

Every column is the same forty-eight blocks on the same weights, and all four agree on the answer to three decimals. That is what the agreement is checked on, because the greedy token is not stable under any of this: on weights drawn from a generator the top two logits of a position are usually a tie, and fusing a reduction differently flips it. For the same reason the decode loop feeds the next token of a fixed sequence rather than the model's own argmax. Sampling would send the four columns down four different sequences.

Where the time goes¶

The rows vary one thing: how much of the wall clock is the Python between operators.

At b1 s32 the GEMMs are small and reading the weights is about three milliseconds, so the rest is dispatch. Removing it is worth 18% to the region and 31% to torch.compile.
At b8 s512 the same forty-eight blocks are 95 ms of cuBLAS and the Python is a rounding error. The first three columns are the same number, and they should be.
The training step is the same story with a backward pass on it.

Decode is the far end of that axis. One token through forty-eight blocks is 6144 region calls for 128 tokens, with almost no arithmetic per call. The region takes 19% off eager, torch.compile 29%, and the two are within 14% of each other. It is the row where removing the interpreter is worth the most, and it is still the row where a fusing compiler wins.

What a region leaves unfused¶

A region removes the interpreter between the operators and leaves the kernels where they were. Where Inductor's own win is also interpreter overhead, the two land close. Inductor could fuse the two layer_norms, the gelu, and the residual adds, but there is not enough of that in a model that is 98% GEMM for the fusion to show at this size.

The columns do not show what each cost to get. The region is one C++ translation unit, 19 s to build once and cached on disk for every process after. torch.compile spends about 18 s on the first shape's forward graph in every new process, and several minutes to reach all five shapes and the backward with a cold Inductor cache.

The reduce-overhead column¶

reduce-overhead (the CUDA-graphs mode) is slower here in every row than plain torch.compile, and slower than eager too. It reproduces: the same numbers in its own process with a single shape rather than five, and no warning from Inductor that it declined to capture.

It also does not run a KV cache as written. A cache is a tensor the caller keeps across invocations, which is the memory CUDA graphs reuse. The decode loop fails outright (accessing tensor output of CUDAGraphs that has been overwritten by a subsequent run) until the counterpart announces each invocation with cudagraph_mark_step_begin() and sets cudagraph_trees_generation_cloning = "user_visible", the two things its own error message asks for. The cloning that second switch turns on is part of what the last column costs.

The training step needed the same care for the same reason: gradients are set to None before each backward rather than accumulated, which is what an optimizer does anyway.

Setup¶

The PPy column is ppy <file> (the staged Python that loads the compiled regions), not ppy run. Compiling the driver natively has nothing to do here: outside the regions the program is forty-eight calls and a loop.

NVIDIA GeForce RTX 4090 24 GB, driver 580.126.20, CUDA 13.0; PyTorch 2.14.0+cu130 on CPython 3.13.15, on a rented Pod. This is not the machine the other tables come from: scripts/compare_docs.py skips this comparison unless PPY_TORCH_CUDA_PYTHON names a PyTorch with a device, so the bench runner leaves it alone, and scripts/cloud/runpod_bench.py is what re-measures it.

Training a torch MLP: Compared with PyTorch as it is usually written¶

From Training a torch MLP. The benchmark covers the standardization of 20,000×16 rows and one forward pass of the step's operators. The programs are in compare/:

train_bench.ppy, under ppy run and, the same file, under python
train_eager.py
train_compile.py

Milliseconds, best of five (the forward pass best of two hundred), over five processes; eight PyTorch threads.

PPy standardizes in the loop the script was written with, over Buffer[float]. PyTorch is how the same preprocessing is written for PyTorch, vectorized over the batch, with roll for the interaction term:

def standardize(raw: Buffer[float], out: Buffer[float], rows: int, cols: int) -> float:
    total: float = 0.0
    for row in range(rows):
        base: int = row * cols
        target: int = row * cols * 2
        sum_: float = 0.0
        for i in range(cols):
            sum_ += raw[base + i]
        mean: float = sum_ / cols
        ...

def standardize(raw: torch.Tensor) -> tuple[torch.Tensor, float]:
    mean = raw.mean(dim=1, keepdim=True)
    deviation = torch.sqrt(((raw - mean) ** 2).mean(dim=1, keepdim=True)) + 1e-8
    z = (raw - mean) / deviation
    interaction = z * z.roll(-1, dims=1)
    out = torch.cat([z, interaction], dim=1)
    return out, float(out.sum())

forward_loss is the same six operators on every side: one ATen region under PPy, eager under PyTorch, and under torch.compile in the third column.

	PPy `ppy run`	CPython, the same file	PyTorch, vectorized	`torch.compile`
standardize, 20000 rows	0.83 ± 0.02	62.74 ± 0.20	0.72 ± 0.04	0.72 ± 0.08
forward pass, per call	0.11 ± 0.01	0.11 ± 0.00	0.11 ± 0.00	0.11 ± 0.00

Standardize. The loop as the script wrote it, compiled, is level with the vectorized rewrite. Both are one pass over the rows, and the rewrite costs seven tensor temporaries the loop never makes. Under python the same loop is the cost the conversion removed.
Forward pass. This is PyTorch's on every side: the matmul over 20,000 rows is the time, and the Python round trips between six operators are not. The region neither gains nor loses there, and torch.compile's guards are a few microseconds on top.

Intel Core Ultra 9 386H; PyTorch 2.14.0 (CPU) on CPython 3.13.13, PPy against the same PyTorch on CPython 3.14.5, from a checkout on a native filesystem.

Columnar expressions: Compared with pandas and polars¶

From Columnar expressions. The two expressions over eight million rows, a NaN in every seventh row of s and every eleventh of t, in compare/: frames_bench.ppy, frames_pandas.py, frames_polars.py. Times are milliseconds, best of five calls, over five processes.

PPy is the pandas expression in a function.
pandas is the same expression with no function around it, one pass and one temporary per operator (numexpr under it, where installed).
polars is the expression rewritten in its own vocabulary, with the NaNs as nulls, planned and run across threads.

@ppy.pure
def blend(s: pd.Series, t: pd.Series) -> pd.Series:
    return s * t + s.fillna(0.0)

def blend(frame):
    mixed = pl.col("s") * pl.col("t") + pl.col("s").fill_null(0.0)
    return frame.select(mixed.alias("out"))["out"]

	PPy	pandas	polars
blend	22.47 ± 1.63	31.09 ± 1.11	13.24 ± 0.42
above	4.65 ± 0.37	6.78 ± 0.19	4.91 ± 0.19

The fused loop reads s and t once and writes the answer once, on one thread, into the array the result Series wraps. pandas makes a temporary per operator, three passes over 64 MB each, and numexpr under it splits the arithmetic across threads to make that up. polars plans the expression and runs it across every core, which is what its row is.

The point of the PPy column is the source: it is the pandas line, unchanged, with pandas' own NaN convention. The polars column is a different program.

Intel Core Ultra 9 386H (16 threads); pandas 3.0.5 with numexpr 2.14.2, polars 1.44.2 on CPython 3.12.13, PPy with pandas 3.0.5 on CPython 3.14.5, from a checkout on a native filesystem.

Comparisons¶

Numerics: Compared with Numba, Codon, Mojo, C, and Rust¶

Value classes: Compared with Numba's @jitclass¶

Tuples: Compared with Numba¶

Buffers and JIT: Compared with Numba, Cython, NumPy, and C¶

Algorithms: The eight kernels against C¶

The eight kernels against Numba, Mojo, and Codon¶

The six problems¶

ppy build --unsafe and --standalone --unsafe¶

Reading input¶

Parallel ranges: Compared with Numba, Taichi, Mojo, and NumPy¶

Derivatives: Compared with JAX and PyTorch¶

Coroutines: Compared with asyncio and uvloop¶

Generics: Compared with Numba¶

Regular expressions: Compared with CPython's re and Rust's regex¶

GPU kernels: Compared with CuPy, Numba, Mojo, and CUDA C¶

How each tool writes the block max¶

Results¶

Tile kernels: Compared with Triton and Taichi¶

NumPy fusion: Compared with NumPy, numexpr, Numba, and JAX¶

Parallel fused kernels: Compared with NumPy, numexpr, Numba, and JAX¶

PyTorch ATen regions: Compared with PyTorch eager and torch.compile¶

GPT-2 XL: one region per block: Compared with PyTorch eager and torch.compile at model scale¶

Where the time goes¶

What a region leaves unfused¶

The reduce-overhead column¶

Setup¶

Training a torch MLP: Compared with PyTorch as it is usually written¶

Columnar expressions: Compared with pandas and polars¶

Value classes: Compared with Numba's `@jitclass`¶

`ppy build --unsafe` and `--standalone --unsafe`¶

Regular expressions: Compared with CPython's `re` and Rust's `regex`¶

PyTorch ATen regions: Compared with PyTorch eager and `torch.compile`¶

GPT-2 XL: one region per block: Compared with PyTorch eager and `torch.compile` at model scale¶