Parallelism

Deterministic reductions

A par reduction combines its elements in one fixed order that depends only on how many elements there are, never on the number of threads or the device. So a float par sum prints the same bits on 1 thread, 32 threads or a GPU. The order is a tree, not left to right, so it can differ from a sequential sum in the last bits, usually by being closer to the exact result.

The order

For par sum, par product, par min, par max, par any and par all:

Leaves. The elements are split into leaves of 16 consecutive elements: elements 0-15 are leaf 0, 16-31 are leaf 1, and so on. The last leaf may be shorter.
Inside a leaf, left to right. Each leaf is folded from its first element to its last, starting from the reducer’s start value (0 for sum, 1 for product): ((((0 + x0) + x1) + x2) + ... + x15).
Between leaves, a balanced tree. The leaf results are combined level by level: leaves 0 and 1, 2 and 3, 4 and 5, …; then those results in pairs; and so on until one value is left. A node without a partner at the end of a level moves up a level unchanged.
Filters. With if, elements that don’t pass are skipped inside their leaf. A leaf whose elements were all filtered out keeps its place in the tree but has no value, so its partner moves up unchanged.

Every step depends only on positions, so the result depends only on the elements and their count.

A picture

par sum(x for x in xs) with 80 elements makes 5 leaves:

elements   0..15     16..31    32..47    48..63    64..79
             |         |         |         |         |
leaves      L0        L1        L2        L3        L4       each folded left to right
             \        /           \       /          |
level 1     (L0 + L1)            (L2 + L3)          L4       L4 has no partner: moves up
                  \                 /                |
level 2       ((L0 + L1) + (L2 + L3))               L4
                          \                         /
level 3       (((L0 + L1) + (L2 + L3)) + L4)                 the result

A sequential sum(x for x in xs) is one long left fold instead: (((0 + x0) + x1) + ...) + x79.

Why a fixed order

Floating-point addition isn’t associative: (a + b) + c and a + (b + c) can round differently. A reduction that let each thread add up its own share and then combined the shares as they finished would print different bits depending on the thread count, the machine, and timing. Tests that compare output would fail for no reason, and a bug would be hard to tell from noise.

Brine fixes the tree instead. Threads (and GPU blocks) only decide who computes each subtree, never which elements go in it, so:

the result is the same on 1 thread, 32 threads, a 128-core server or a GPU;
the result is the same from run to run;
dev, --release and native and C backends agree.

The tree also loses less precision. In a left-to-right sum of n terms the rounding errors can grow with n; in a balanced tree they grow with the height, about log2(n). Adding small terms to a large running total, the usual source of lost bits, mostly turns into adding numbers of similar size.

Same bits everywhere: a demonstration

This program adds up the first ten million terms of the harmonic series (1 + 1/2 + 1/3 + …), once with par sum and once with a sequential sum:

fn main():
    let n = 10_000_000
    let tree = par sum(1.0 / float(i) for i in 1..=n)
    let left = sum(1.0 / float(i) for i in 1..=n)
    print(f"par sum: {tree}")
    print(f"sum:     {left}")

Built with --release and run on a Ryzen 9 9950X3D and an RTX 4090 with different thread counts, with the GPU off and on:

brine build harmonic.br --release -o harmonic
BRINE_THREADS=1 BRINE_GPU=off ./harmonic
BRINE_THREADS=2 BRINE_GPU=off ./harmonic
BRINE_THREADS=8 BRINE_GPU=off ./harmonic
BRINE_THREADS=32 BRINE_GPU=off ./harmonic
BRINE_GPU=on BRINE_GPU_LOG=1 ./harmonic

Every run prints the same two lines:

par sum: 16.69531136585985
sum:     16.695311365857272

The GPU run also logs where the reduction ran:

brine gpu: NVIDIA GeForce RTX 4090 ready in 88 ms
brine gpu: harmonic.br:3:16 par reduction n=10000000 -> gpu 1.252 ms (kernels 0.321 ms, copies 0.930 ms, 0.0 MB; cpu estimate 0.000 ms)

The exact sum is 16.6953113658598518153991… The par sum result is the closest float to it (it is off by 7.2 × 10⁻¹⁶, a fifth of a unit in the last place). The sequential sum is off by 2.6 × 10⁻¹², about 730 units in the last place, because by the end it adds terms near 10⁻⁷ to a total near 16.

The two lines differ, and they always will: they’re different orders. What never changes is each line.

Integers and overflow

Integer par sum and par product check every addition and multiplication, like all integer arithmetic in Brine, and the checks happen in the tree order. Whether a reduction overflows is decided by the partial results of that order, so it is the same on every machine, but it can differ from a sequential sum.

The tree can overflow where the left fold doesn’t. Here leaf 1 holds two copies of the largest int, and folding that leaf overflows, while the left fold’s running total never leaves the range:

fn main():
    let big = 9_223_372_036_854_775_807
    var xs = [0] * 32
    xs[0] = -big
    xs[16] = big
    xs[17] = big
    print(sum(x for x in xs))
    print(par sum(x for x in xs))

9223372036854775807

panic: integer overflow in addition
  --> ov.br:8:11

And the other way round: here the left fold passes through big + 5, while the tree combines big - 5 with 10 - 10 and never leaves the range:

fn main():
    let big = 9_223_372_036_854_775_807
    var xs = [0] * 64
    xs[0] = big
    xs[16] = -5
    xs[32] = 10
    xs[48] = -10
    print(par sum(x for x in xs))
    print(sum(x for x in xs))

9223372036854775802

panic: integer overflow in addition
  --> ov2.br:9:11

When the exact total fits, both orders give the same integer, since integer addition is associative as long as nothing overflows. If partial results may come near the limits, use a wider type (i128) or float.

min, max, any and all never overflow. min and max compare the way the sequential versions do: a later value replaces the current one only if it is strictly smaller (or larger).

par sum and sum

	`sum(...)`	`par sum(...)`
order	left to right, from `0`	leaves of 16, then a balanced tree
threads	one	the pool, and the GPU when it pays
float result	the left fold’s bits	the tree’s bits, the same on any machine
integer overflow	of the running total	of any partial result in the tree
accuracy for long float sums	error can grow with `n`	error grows with `log2(n)`

Both are deterministic. They’re just different orders, so pick one and keep it. If you need a result that matches another language’s plain loop bit for bit, use the sequential sum (or a for loop).

The 0.4.0 change

Up to Brine 0.3.0, a par reduction split its elements into at most 256 chunks, summed each chunk left to right and added the chunk results in order. That was also independent of the thread count, but a GPU needs thousands of independent partial results, so 0.4.0 changed to the leaves and tree on this page. It is a breaking change (see the changelog):

par sum and par product of floats can differ from 0.3.0 in the last bits, usually closer to the exact result;
whether an integer par sum or par product overflows is decided by the new order.

Programs that need the old bits should sum sequentially.

Parallelism overview: the reducers and filters
GPU offload: how the GPU runs the same tree
Multiple processes: world: world.sum combines ranks in rank order