Standard library · Numbers

num/grid

Imported as import "num/grid" as grid, its names are then grid.…. Every signature below is the one the checker infers.

grid — arrays of numbers with a shape, the way NumPy has them.

A grid is one flat Buf(Float) and the shape laid over it: shape says how many along each axis, strides how far to step in the buffer for one along each axis, off where element zero sits. transpose, rows, col and slice1 are views — the same buffer, a different shape — and cost nothing; copy and every operation that makes a new grid lay their result out contiguously, which is the shape every fast path here recognises.

The arithmetic (add, mul, gsqrt, sum, matmul …) is written as plain loops over the buffer rather than through a lambda per element, so that the compiler can hoist the bounds checks and vectorize; map and zip take a lambda and are for everything the named operations do not cover. A grid is storage, as a Buf is: put1/put2 write in place, and the views share it.

A comparison makes a Mask — one byte per element, as NumPy's boolean arrays are — and where, mcount, any and all read one.

What is here, by section: shapes and views; making grids; elementwise arithmetic and the rest of the maths; masks and where; reductions, including the nan-skipping ones and the order statistics; joining, sorting and searching; and solving, on an LU factorisation with partial pivoting. What is not here: dtypes other than float, more than two dimensions, broadcasting between different shapes, and the decompositions beyond LU.

Three things keep it fast. An operation writes its result into an operand's own storage when that operand is a temporary nobody else can see (unique), so a * 2.0 + 1.0 makes one grid and not two. The reductions run eight partial sums in four vector lanes, and matmul is a register-blocked kernel written in Float2. And on a --parallel build with RILL_THREADS set, anything over a quarter of a million elements is split across the workers in strands; a single-worker run never spawns.

import "num/grid" as g a = g.arange(6) |> g.reshape2(2, 3) b = g.transpose(a) println(g.matmul(a, b)) # [[5, 14], [14, 50]]

From Python it is import grid as g — bindings/python builds this library into a shared object and gives its handles NumPy's notation.

import "num/grid" as g
a = g.arange(6) |> g.reshape2(2, 3)
a                                  # => [[0, 1, 2], [3, 4, 5]]
b = g.transpose(a)
g.matmul(a, b)                     # => [[5, 14], [14, 50]]
g.sum(a)                           # => 15
g.gtk(a, 2.0)                      # => [[false, false, false], [true, true, true]]
g.where(g.gtk(a, 2.0), a, g.zeros_like(a))   # => [[0, 0, 0], [3, 4, 5]]

Types

`Grid`

Grid(shape: Buf(Int), strides: Buf(Int), off: Int, data: Buf(Float))

`Mask`

1 where a comparison held, 0 where it did not; always contiguous.

Mask(shape: Buf(Int), bits: Buf(U8))

`Acc`

Acc(c00: Float2, c01: Float2, c02: Float2, c03: Float2, c04: Float2, c10: Float2, c11: Float2, c12: Float2, c13: Float2, c14: Float2, c20: Float2, c21: Float2, c22: Float2, c23: Float2, c24: Float2, c30: Float2, c31: Float2, c32: Float2, c33: Float2, c34: Float2)

`Lu`

The LU factors of a square matrix, worked out in one pass with partial pivoting: a holds L below the diagonal (its own diagonal is all ones, and is not stored) and U on and above it, piv says which row of the original each row came from, sign is 1 or -1 by how many rows were exchanged, and ok is false when the matrix turned out to be singular. Everything below solves through this: one factorisation, then as many right-hand sides as wanted.

Lu(a: Grid, piv: Buf(Int), sign: Float, ok: Bool)

`Blk`

c -= a * b, where all three sit inside bigger buffers: coff, aoff and boff say where each starts and cld, ald and bld how wide the buffer each lives in is. c is m by nn, a is m by kk, b is kk by nn.

Four rows by eight columns at a time, the sixteen answers kept in Float2 registers while the shared dimension is walked — the same shape of kernel as matmul, and for the same reason: written as one multiply-add per sixteen bytes moved it is the memory bus that answers, not the processor. Both halves of solving go through this: the trailing update of the factorisation, and the substitutions afterwards.

Blk(c00: Float2, c01: Float2, c02: Float2, c03: Float2, c10: Float2, c11: Float2, c12: Float2, c13: Float2, c20: Float2, c21: Float2, c22: Float2, c23: Float2, c30: Float2, c31: Float2, c32: Float2, c33: Float2)

Implementations

impl Show(Grid)
impl Show(Mask)

Functions

`fn ndim(g: Grid) -> Int`

a = g.arange(6) |> g.reshape2(2, 3)
g.ndim(a)              # => 2
g.ndim(g.arange(3))    # => 1

`fn size(g: Grid) -> Int`

a = g.arange(6) |> g.reshape2(2, 3)
g.size(a)   # => 6

`fn count(shape: Buf(Int)) -> Int`

g.count(g.shape2(2, 3))   # => 6

`fn rows_of(g: Grid) -> Int`

a = g.arange(6) |> g.reshape2(2, 3)
g.rows_of(a)   # => 2

`fn cols_of(g: Grid) -> Int`

a = g.arange(6) |> g.reshape2(2, 3)
g.cols_of(a)   # => 3

`fn shape1(n: Int) -> Buf(Int)`

#g.shape1(4)      # => 1
g.shape1(4)[0]    # => 4

`fn shape2(r: Int, c: Int) -> Buf(Int)`

s = g.shape2(2, 3)
s[0]   # => 2
s[1]   # => 3

`fn same_shape(a: Grid, b: Grid) -> Bool`

g.same_shape(g.zeros2(2, 3), g.ones2(2, 3))   # => true
g.same_shape(g.zeros2(2, 3), g.ones2(3, 2))   # => false

`fn contiguous(g: Grid) -> Bool`

Whether the elements sit in the buffer in order, starting at its front.

a = g.arange(6) |> g.reshape2(2, 3)
g.contiguous(a)                # => true
g.contiguous(g.transpose(a))   # => false

`fn zeros(shape: Buf(Int)) -> Grid`

g.zeros(g.shape2(2, 2))   # => [[0, 0], [0, 0]]

`fn zeros1(n: Int) -> Grid`

g.zeros1(3)   # => [0, 0, 0]

`fn zeros2(r: Int, c: Int) -> Grid`

g.zeros2(1, 2)   # => [[0, 0]]

`fn full(shape: Buf(Int), v: Float) -> Grid`

g.full(g.shape1(3), 7.5)   # => [7.5, 7.5, 7.5]

`fn ones(shape: Buf(Int)) -> Grid`

g.ones(g.shape2(2, 2))   # => [[1, 1], [1, 1]]

`fn ones1(n: Int) -> Grid`

g.ones1(2)   # => [1, 1]

`fn ones2(r: Int, c: Int) -> Grid`

g.ones2(2, 1)   # => [[1], [1]]

`fn full2(r: Int, c: Int, v: Float) -> Grid`

g.full2(2, 2, 3.0)   # => [[3, 3], [3, 3]]

`fn arange(n: Int) -> Grid`

0, 1, …, n - 1

g.arange(5)   # => [0, 1, 2, 3, 4]

`fn linspace(a: Float, b: Float, n: Int) -> Grid`

g.linspace(0.0, 1.0, 5)   # => [0, 0.25, 0.5, 0.75, 1]

`fn of(xs: List(Float)) -> Grid`

A 1-D grid from a list of floats.

g.of(Cons(1.5, Cons(2.5, Nil)))   # => [1.5, 2.5]

`fn from_fn2(r: Int, c: Int, f: (Int, Int) -> Float) -> Grid`

A grid whose element at (i, j) is f(i, j).

g.from_fn2(2, 3, \i, j -> to_float(i * 10 + j))   # => [[0, 1, 2], [10, 11, 12]]

`fn from_fn1(n: Int, f: (Int) -> Float) -> Grid`

g.from_fn1(4, \i -> to_float(i * i))   # => [0, 1, 4, 9]

`fn at1(g: Grid, i: Int) -> Float`

g.at1(g.arange(5), 3)   # => 3

`fn at2(g: Grid, i: Int, j: Int) -> Float`

a = g.arange(6) |> g.reshape2(2, 3)
g.at2(a, 1, 2)   # => 5

`fn put1(g: Grid, i: Int, v: Float) -> Unit`

v = g.zeros1(3)
g.put1(v, 1, 9.0)
v   # => [0, 9, 0]

`fn put2(g: Grid, i: Int, j: Int, v: Float) -> Unit`

m = g.zeros2(2, 2)
g.put2(m, 1, 0, 9.0)
m   # => [[0, 0], [9, 0]]

`fn transpose(g: Grid) -> Grid`

Two-dimensional transpose: the same buffer read the other way round.

a = g.arange(6) |> g.reshape2(2, 3)
g.transpose(a)   # => [[0, 3], [1, 4], [2, 5]]

`fn rows(g: Grid, lo: Int, hi: Int) -> Grid`

Rows lo up to hi of a 2-D grid, or elements lo up to hi of a 1-D one.

m = g.arange(6) |> g.reshape2(3, 2)
g.rows(m, 1, 3)   # => [[2, 3], [4, 5]]
g.rows(g.arange(5), 1, 3)   # => [1, 2]

`fn slice1(g: Grid, lo: Int, hi: Int) -> Grid`

g.slice1(g.arange(5), 2, 4)   # => [2, 3]

`fn row(g: Grid, i: Int) -> Grid`

a = g.arange(6) |> g.reshape2(2, 3)
g.row(a, 1)   # => [3, 4, 5]

`fn col(g: Grid, j: Int) -> Grid`

a = g.arange(6) |> g.reshape2(2, 3)
g.col(a, 2)   # => [2, 5]

`fn copy(g: Grid) -> Grid`

A contiguous grid with the same elements.

a = g.arange(6) |> g.reshape2(2, 3)
t = g.transpose(a)
g.contiguous(g.copy(t))   # => true
g.copy(t)                 # => [[0, 3], [1, 4], [2, 5]]

`fn reshape(g: Grid, shape: Buf(Int)) -> Grid`

g.reshape(g.arange(4), g.shape2(2, 2))   # => [[0, 1], [2, 3]]

`fn reshape2(g: Grid, r: Int, c: Int) -> Grid`

g.reshape2(g.arange(6), 3, 2)   # => [[0, 1], [2, 3], [4, 5]]

`fn flatten(g: Grid) -> Grid`

a = g.arange(6) |> g.reshape2(2, 3)
g.flatten(a)   # => [0, 1, 2, 3, 4, 5]

`fn empty_like(g: Grid) -> Grid`

A grid of the same shape, filled with nothing in particular, with zeros, with ones, or with one number.

a = g.arange(6) |> g.reshape2(2, 3)
g.same_shape(g.empty_like(a), a)   # => true

`fn zeros_like(g: Grid) -> Grid`

a = g.arange(6) |> g.reshape2(2, 3)
g.zeros_like(a)   # => [[0, 0, 0], [0, 0, 0]]

`fn ones_like(g: Grid) -> Grid`

g.ones_like(g.arange(3))   # => [1, 1, 1]

`fn full_like(g: Grid, v: Float) -> Grid`

g.full_like(g.arange(2), 4.0)   # => [4, 4]

`fn squeeze(g: Grid) -> Grid`

A 2-D grid with one row or one column becomes the 1-D grid of its elements; anything else is unchanged.

g.squeeze(g.ones2(1, 3))   # => [1, 1, 1]
g.squeeze(g.ones2(2, 2))   # => [[1, 1], [1, 1]]

`fn ravel(g: Grid) -> Grid`

a = g.arange(6) |> g.reshape2(2, 3)
g.ravel(g.transpose(a))   # => [0, 3, 1, 4, 2, 5]

`fn logspace(a: Float, b: Float, n: Int, base: Float) -> Grid`

n numbers from a to b, evenly spaced in the exponent: logspace takes the ends as powers of the base, geomspace as the values.

g.logspace(0.0, 2.0, 3, 10.0)   # => [1, 10, 100]

`fn geomspace(a: Float, b: Float, n: Int) -> Grid`

g.geomspace(1.0, 8.0, 4)   # => [1, 2, 4, 8]

`fn concat(a: Grid, b: Grid, axis: Int) -> Grid`

Joining: axis 0 adds rows, 1 adds columns. Two 1-D grids join end to end.

g.concat(g.ones2(1, 2), g.zeros2(1, 2), 0)   # => [[1, 1], [0, 0]]
g.concat(g.ones2(1, 2), g.zeros2(1, 2), 1)   # => [[1, 1, 0, 0]]
g.concat(g.arange(2), g.arange(2), 0)        # => [0, 1, 0, 1]

`fn as_rows(g: Grid) -> Grid`

g.as_rows(g.arange(3))   # => [[0, 1, 2]]

`fn vstack(a: Grid, b: Grid) -> Grid`

g.vstack(g.arange(2), g.ones1(2))   # => [[0, 1], [1, 1]]

`fn hstack(a: Grid, b: Grid) -> Grid`

g.hstack(g.ones2(2, 1), g.zeros2(2, 1))   # => [[1, 0], [1, 0]]

`fn stack(a: Grid, b: Grid) -> Grid`

Two 1-D grids as the two rows of a 2-D one.

g.stack(g.arange(2), g.ones1(2))   # => [[0, 1], [1, 1]]

`fn tile(g: Grid, k: Int) -> Grid`

The whole grid laid down k times, and each element k times over.

g.tile(g.arange(2), 2)   # => [0, 1, 0, 1]

`fn repeat(g: Grid, k: Int) -> Grid`

g.repeat(g.arange(2), 2)   # => [0, 0, 1, 1]

`fn diag(g: Grid) -> Grid`

From a 1-D grid, the square matrix with it down the diagonal; from a 2-D one, the diagonal itself.

g.diag(g.arange(2))                          # => [[0, 0], [0, 1]]
g.diag(g.arange(4) |> g.reshape2(2, 2))      # => [0, 3]

`fn diagonal(g: Grid) -> Grid`

g.diagonal(g.arange(4) |> g.reshape2(2, 2))   # => [0, 3]

`fn tril(g: Grid, k: Int) -> Grid`

The lower and upper triangles, keeping everything on or beyond the k-th diagonal and zeroing the rest.

g.tril(g.ones2(2, 2), 0)   # => [[1, 0], [1, 1]]

`fn triu(g: Grid, k: Int) -> Grid`

g.triu(g.ones2(2, 2), 0)   # => [[1, 1], [0, 1]]

`fn meshgrid(x: Grid, y: Grid) -> (Grid, Grid)`

Every pair from x and y, as NumPy lays them out: rows of x repeated down, columns of y repeated across.

(xs, ys) = g.meshgrid(g.arange(2), g.arange(3))
xs   # => [[0, 1], [0, 1], [0, 1]]
ys   # => [[0, 0], [1, 1], [2, 2]]

`fn fliplr(g: Grid) -> Grid`

Reversed: along the last axis, or both axes of a 2-D grid.

g.fliplr(g.arange(3))                          # => [2, 1, 0]
g.fliplr(g.arange(4) |> g.reshape2(2, 2))      # => [[1, 0], [3, 2]]

`fn flipud(g: Grid) -> Grid`

g.flipud(g.arange(4) |> g.reshape2(2, 2))   # => [[2, 3], [0, 1]]

`fn flip(g: Grid) -> Grid`

g.flip(g.arange(4) |> g.reshape2(2, 2))   # => [[3, 2], [1, 0]]

`fn roll(g: Grid, k: Int) -> Grid`

Everything moved along by k, what falls off one end coming back at the other.

g.roll(g.arange(4), 1)    # => [3, 0, 1, 2]
g.roll(g.arange(4), -1)   # => [1, 2, 3, 0]

`fn take(g: Grid, idx: Grid) -> Grid`

The elements at the given positions, which are a grid of numbers.

g.take(g.of(Cons(10.0, Cons(20.0, Cons(30.0, Nil)))), g.of(Cons(2.0, Cons(0.0, Nil))))   # => [30, 10]

`fn dense(g: Grid) -> Grid`

The contiguous copy every operation works from when a view is handed in.

a = g.arange(6) |> g.reshape2(2, 3)
g.dense(g.transpose(a))   # => [[0, 3], [1, 4], [2, 5]]