The Scaly Programming Language

Chapter 6

Lanes

Four Numbers in One

A core can do more than one thing per instruction. Its vector registers hold several numbers side by side — two, four or more doubles, sixteen or thirty-two bytes —, and one instruction adds or multiplies all of them at once, each in its own lane. The cores of the previous chapters multiply the work by the number of cores; lanes multiply it once more, on each of them.

A compiler can sometimes find lanes in a loop by itself, but not in a loop that stops early or that takes a different path for each number. So Scaly lets a program say it: a vector type holds a fixed number of values as one value, and its operators work lane by lane.

let a f64x4(1.0, 2.0, 3.0, 4.0)
let b f64x4.splat(0.5)
let c f64x4.iota(10.0, 0.25)
let d a * b + c
print "`d[0]` `d[1]` `d[2]` `d[3]`"
let e a.with(2, 9.0) * 2.0
print "`e[0]` `e[1]` `e[2]` `e[3]`"
10.5 11.25 12 12.75
2 4 18 8

f64x4(…) gives each lane its value, splat the same value to every lane, and iota(base, step) counts up from base. a * b + c is four multiplications and four additions, done as one of each. A scalar on one side of an operator stands for itself in every lane, as the 2.0 does. d[i] reads a lane and stops the program with exit code 15 when there is no lane i; a.with(i, x) is a copy of a with one lane replaced. A vector is a value like a double: it is copied, passed and returned in registers, and it never lives on a page.

The types are named by their element and their number of lanes: f64x2, f64x4, f64x8; f32x4, f32x8, f32x16; i64x2, i64x4, i32x4, i32x8; u8x16, u8x32, u32x4, u32x8. A vector wider than the hardware is fine: on a core whose registers hold two doubles, an f64x4 is two registers, and its operations are two instructions each.

Masks

A comparison of two vectors compares lane by lane, and its answer is a mask — mask2 to mask32, one flag per lane. A mask answers questions about its lanes, combines with and, or and not lane by lane, and chooses between two vectors with select:

let a f64x4(1.0, 2.0, 3.0, 4.0)
let big a > 2.5
print "`big.count()` lanes above 2.5, any: `big.any()`, all: `big.all()`"
let clipped select(big, f64x4.splat(2.5), a)
print "`clipped[0]` `clipped[1]` `clipped[2]` `clipped[3]`"
print "sum `a.sum()`, largest `a.max_lane()`"
2 lanes above 2.5, any: true, all: false
1 2 2.5 2.5
sum 10, largest 4

Masks are how a vector program takes a different path per lane: it computes both paths for all lanes and keeps, with select, what each lane should have. And any and all are how it leaves a loop: when no lane has anything left to do. sum, product, min_lane and max_lane fold the lanes into one number; a sum of floating-point lanes adds them in order, lane 0 first, so it gives the digits of a plain loop over them.

Lanes from Memory

f64x4.load(xs, i) reads the elements xs[i] to xs[i + 3] of a Slice, Vector or Array as one vector, and v.store(ys, i) writes them. Each checks once that all four elements exist — a scalar loop checks every one. A dot product four elements at a time, with a fused multiply-add per step and the rest done one by one:

function dot(a: Vector[double], b: Vector[double]) returns double
{
    var acc f64x4.splat(0.0)
    var i: size_t 0
    while i + 4 <= a.length
    {
        acc := f64x4.load(a, i).fma(f64x4.load(b, i), acc)
        i := i + 4
    }
    var sum acc.sum()
    while i < a.length
    {
        sum := sum + a[i] * b[i]
        i := i + 1
    }
    sum
}

let n 1001
var a Vector[double](n as size_t)
var b Vector[double](n as size_t)
for i in n
{
    a.put(i as size_t, (i as double) * 0.5)
    b.put(i as size_t, 2.0)
}
print "a · b = `dot(a, b)`"
a · b = 500500

v.fma(w, acc) is v * w + acc with one rounding instead of two; sqrt, abs, floor, min and max work lane by lane as well. The four lanes of acc are four running sums, added up at the end. That is a different order of additions than a loop over one element at a time, and for numbers that are not whole it can change the last digits; the order is written in the program, so the result is the same on every run and every machine.

Mandelbrot, Eight Points at a Time

The program of A Plain for counts one point after the other, and each point leaves its loop when it escapes. With lanes, a row is counted eight points at a time, as two f64x4 — two, because while one multiplication is still in flight the core can already start the other. The eight points go through their iterations together; every eight iterations the program asks whether any lane is still inside, and stops when none is. A point that escaped keeps iterating with the others; after fifty iterations its value is so large that it has become infinite or not a number, and the comparison with 4.0 sorts it out at the end:

function count8(x: int, y: int, n: int) returns int
{
    let step 2.0 / (n as double)
    let ci f64x4.splat(step * (y as double) - 1.0)
    let cr0 f64x4.iota(step * (x as double) - 1.5, step)
    let cr1 cr0 + 4.0 * step
    var zr0 cr0
    var zi0 ci
    var zr1 cr1
    var zi1 ci
    var s0 f64x4.splat(0.0)
    var s1 s0
    var block 0
    while block < 6
    {
        var k 0
        while k < 8
        {
            let r0 zr0 * zr0
            let i0 zi0 * zi0
            let ri0 zr0 * zi0
            s0 := r0 + i0
            zr0 := r0 - i0 + cr0
            zi0 := ri0 + ri0 + ci
            let r1 zr1 * zr1
            let i1 zi1 * zi1
            let ri1 zr1 * zi1
            s1 := r1 + i1
            zr1 := r1 - i1 + cr1
            zi1 := ri1 + ri1 + ci
            k := k + 1
        }
        if ((s0 <= 4.0) or (s1 <= 4.0)).any() = false
            return 0
        block := block + 1
    }
    var k 0
    while k < 2
    {
        let r0 zr0 * zr0
        let i0 zi0 * zi0
        let ri0 zr0 * zi0
        s0 := r0 + i0
        zr0 := r0 - i0 + cr0
        zi0 := ri0 + ri0 + ci
        let r1 zr1 * zr1
        let i1 zi1 * zi1
        let ri1 zr1 * zi1
        s1 := r1 + i1
        zr1 := r1 - i1 + cr1
        zi1 := ri1 + ri1 + ci
        k := k + 1
    }
    (s0 <= 4.0).count() + (s1 <= 4.0).count()
}

function row_count(y: int, n: int) returns int
{
    var count 0
    var x 0
    while x < n
    {
        count := count + count8(x, y, n)
        x := x + 8
    }
    count
}

let n 2000
var counts Vector[int](n as size_t)
for y in n
    counts.put(y as size_t, row_count(y, n))

var total 0
for c in counts
    total := total + c
print "mandelbrot `n`: `total` points inside"
$ scaly run mandelbrot.scaly
mandelbrot 2000: 1588053 points inside
$ scalyc --task-plan mandelbrot.scaly
task-plan: 2 for-loops, 2 parallelizable
task-plan: mandelbrot.scaly:71:1: parallel
task-plan: mandelbrot.scaly:75:1: parallel
$ SCALY_WORKERS=1 scaly run mandelbrot.scaly
mandelbrot 2000: 1588053 points inside

The count is the one of the scalar program. The row loop is the same plain for as before, and it still runs on every core: a vector is a value like a number, and count8 is a function like row_count was. Built with scaly build --release and n at 8000, the scalar program took 4.1 seconds on one core and 0.54 on all ten; this one takes 0.76 on one core and 0.11 on all ten — the lanes bring five times, the cores another seven, together thirty-seven.

Lanes and Cores

Lanes and cores combine because they divide the work at different places: the lanes inside one round, the cores across the rounds. What the compiler checks for a for does not change. One thing a round may not do with lanes is store a vector into memory outside the round: the compiler proves that two rounds write different elements only for put at the loop's index, not for the four elements a store writes.

let n 1000
var xs Vector[double](n as size_t)
for i in n
    xs.put(i as size_t, i as double)

var ys Vector[double](n as size_t)
for k in n / 4
{
    let v f64x4.load(xs, (4 * k) as size_t)
    (v * 2.0).store(ys, (4 * k) as size_t)
}
print "ys[999] = `ys[999]`"
$ scaly run double.scaly
ys[999] = 1998
$ scalyc --task-plan double.scaly
task-plan: 2 for-loops, 1 parallelizable
task-plan: double.scaly:3:1: parallel
task-plan: double.scaly:7:1: blocked: SIMD store

Such a loop runs in order, on one core, with its lanes. A loop whose rounds load and compute and hand back numbers — a sum, a count, a slot of their own — runs on all cores with its lanes, as the rows of the Mandelbrot set do.

Which instructions the lanes become depends on the processor the program is compiled for, and without being told, the compiler compiles for a generic one of the target's family: on x86-64 that means no AVX, even where the machine has it. scaly build --release -mcpu=native compiles for the machine at hand (or -mcpu= with a processor's name, for another one). On the laptop of these measurements it changed nothing: its vector registers are as wide as the generic ones.