# Viv: the full reference for language models

> The language, every builtin and which engine runs it, every documented example with its output, the CLI, and the post format, in one file. The compact version is https://viv-lang.org/viv-prompt-pack.md.

Viv is a small language for numerical programs. A number is a 64-bit
float (`Num`) unless it is declared otherwise: `n: Int` is a 64-bit integer,
`Float32` an IEEE binary32 number, and `2.0 + 3.0i` a complex one. Arrays are
rank 1 to 4, and programs end in `test` blocks that check their claims. One
program runs in three engines: `viv run f.viv` interprets it (the site runs the
same interpreter compiled to wasm), `viv f.viv -o f.c` compiles it to C99, and
`viv f.viv -o f.wasm` compiles it to wasm. All three are meant to print the
same bytes. Figures are records on stdout, so a plot is part of the output and
needs no library.

How to use this pack: write Viv, not Python. The grammar is short and differs
from Python in the ways the "Python habits" table lists; the examples show
every feature with the exact output it prints.

## Reference

### Lexical structure

- **Numbers**: A number is a 64-bit float. `10` lexes as `10.0`, scientific notation works (`1e-10`), and a 64-bit integer is an `Int`, declared with an annotation (`n: Int`). `%` and `//` follow Python's rules for floats (see operators (#operators)); `mod(a, b)` is C's fmod, whose result takes the sign of `a`.

- **No implicit multiplication**: A name cannot follow a number directly: `2x` is an error that says to write `2 * x`. The `i` of an imaginary literal (`2.5i`) and the `e` of an exponent (`1e-10`) are part of the number.

- **Names**: A name starts with a letter or `_` and continues with letters, digits and `_`. Letters are Unicode letters, so `Δx` and `σ` are names. Names that start with `viv__` are used by the parser for the hidden variables of `a[i] += x` and the longer `for` forms.

- **The dot rule**: A `.` is part of a number only with digits on both sides. So `1.5` is a float, `x.y` is field access, and `1.` is neither.

- **Statements end at the newline**: No semicolons. Inside parentheses or brackets the lexer suppresses newlines, so long expressions wrap freely.

- **Comments**: `#` to end of line, dropped by the lexer. The annotation markers — `#|` what, `#||` how, `#|||` why, `#~` tags, `#:` variable glossary, `#[`/`#]` regions — are comments the compiler ignores and `viv doc` extracts. Above a `def` they form its structured doc.

- **Strings**: `"text"`, with the escapes `\n` `\t` `\r` `\0` `\\` `\"` `\'` and `\u{1F600}` (any other backslash is an error), and f-strings, `f"x = {x:.3f}"`. A string is a value: see Strings (#strings).

- **Docstrings**: A bare string as the first statement of a `def` documents it and compiles to nothing.

### Statements

| Form | Notes |
|---|---|
| `x = expr` | Assignment. At top level it creates a global; inside a `def`, a local hoisted to function scope. |
| `x += e` `-=` `*=` `/=` | Compound assignment; desugars to the plain form, `x = x + e`. |
| `a[i, j] += e` `-=` `*=` `/=` | Compound assignment on an array element: `a[i, j] = a[i, j] + e`, with each index evaluated once. An index that is not a number or a variable is computed into a hidden variable first, so `c[floor(10 * rand())] += 1` draws one random number, not two. |
| `a[i, j] = expr` | Indexed store, up to rank 4. The base may be any array-valued expression, including a record field. |
| `a[lo:hi] = v` `A[:, j] = v` | A write into a slice (#slices): `v` is a number, which every selected element takes, or an array of exactly the slice's shape. The value is copied before the write, so `s[1:5] = s[0:4]` shifts. `a[lo:hi] += v` and the other compound forms read the slice, combine and write it back, with each bound evaluated once. |
| `return expr` | A `def` that ends without a `return` gives nothing: it is called as a statement, and the compiler refuses a use of its value (`show f(x)`, `f(x) + 1.0`). |
| `return a, b` | Several results. Every `return` of the def gives the same number of values, each with one type (a number, a boolean, a complex number, an array, a record); `return g(x)` passes on the values of a def `g` that returns as many. The caller unpacks them, `a, b = f(x)`, on a line of its own, or calls `f(x)` alone and drops them. `f(x) + 1`, `show f(x)` and `f` as a function value are errors that name `f`. The def returns a record the compiler declares for it (Records (#types)), so every engine runs it. |
| `a, b = e` | Unpacking. `e` is a call of a def that returns as many values, a list of values (`a, b = b, a` swaps: every value on the right is computed first, then the targets are assigned left to right), or an array with exactly one element per name, such as `B, C, H, W = shape(X)`; an array of another length stops the run, located at the line. A target is a name, an element `v[i]` or a slice `v[lo:hi]`, and its indices are computed after the right side, as in Python: `v[i], v[j] = v[j], v[i]` swaps two elements. `loss, g = value_and_grad(...)` is this form too. |
| `if c { } else { }` | With `else if` chains. |
| `while c { }` | With `break` and `continue`. |
| `for i to n { }` | `i` runs 0 … n−1. |
| `for i from a to b { }` | `i` runs a, a+1, … while it is below b. As in `for i to n`, `to` names the end and the end is left out. |
| `for i from a to b by s { }` | `i` runs a, a+s, a+2s, … and stops before reaching b: below b for a positive s, above b for a negative one (`for i from 10 to 0 by -3` is 10, 7, 4, 1). The loop makes ceil((b − a) / s) passes, none when that is 0 or less, and pass k has `i = a + k·s`. a, b and s are evaluated once, in that order. A step of 0 is an error: at compile time when it is a number, at run time otherwise. |
| `for i to n reversed { }` | The same values in the opposite order: `for t to T reversed` is T−1, …, 0 and `for i from 1 to n reversed` is n−1, …, 1. `reversed` goes last, after any `from` and `by`. These longer forms are rewritten by the parser into `for k to count` over a hidden counter, with `i` assigned at the top of each pass, so every engine runs them. |
| `parallel for i to n { }` | Sequential in the interpreter; the C backend emits an OpenMP pragma. Only this form: no `from`, `by` or `reversed`. Each pass may write only the elements of a shared array indexed by `i` (`a[i]`, `m[i, j]`, `m[:, i] = col`, `setslab(W, i, mine)`, or a def given `i` as the index it writes at), and may read a shared array it writes only at those elements; a write two passes could share, a write of a global through a def, a `break` or a `return` is refused before the run, in every engine. Arrays a pass makes are its own, and are freed when the pass ends. |
| `show expr` | Prints the value. A variable, or `fmt` of one, is labeled with its name (`x: 2`); any other expression prints bare, as `print` would. |
| `show a, b, c` | Several labeled values on one line, two spaces apart: `a: 1.5 b: 0.3 v: [1, 2]`. A variable, or `fmt` of one, is labeled with its name, so `show fmt(x, 17)` prints `x: 0.30000000000000004`. Arrays print on one line here. It is `print("a: ", a, " b: ", b, ..., sep="")`; see Output (#output). |
| `print(...)` | Text and values on one line, a space between them; see Output (#output). |
| `for x in e { }` | Each element of `e` in order: a string's characters, a list's strings, a vector's numbers. `for i, x in enumerate(e)` gives the index too, and `enumerate(e, 1)` (or `start=1`) counts from 1. `for a, b in zip(u, v)` walks several sequences together and stops at the end of the shortest. `for row in rows(M)` gives the rows of a matrix, `M[k, :]` (an array's rank is known only when the program runs, so `for x in M` reads numbers and stops on a matrix). `e` is evaluated once; the parser writes the loop as `for k to len(e)` with `x = e[k]` at the top of each pass. |
| `stop("text")` | Ends the run with `error: text` and the line it is on, exit status 1, in every engine: the error a failed check should give. The argument is one string, a literal or a value (`stop(f"no atom {name}")`). Output printed before it stays printed (the wasm engine reports the error alone, as for every run-time error). |
| `narrate "text"` | A progress line. A non-string argument is commentary and compiles to nothing. |
| `test "name" { }` | Runs the body and reports under the name: one `PASS name` or `FAIL name` line, however many asserts are inside. A failure makes the process exit nonzero. |
| `def f(a, b, c = 1.0) { }` | A parameter with a default comes after those without and is passed by name only: `f(1, 2, c=0.5)`; a default may read the parameters before it (`scale = 1.0 / sqrt(rows)`) and is computed at each call that leaves it out. The front end fills every call out, so it costs nothing at run time; a call through a function value passes its arguments in order. Optional annotations: `def f(a: Num, v: Arr) -> Num`. Recursion works to a depth of 2500 nested calls, the same in every engine; a call past that stops the run with `recursion deeper than 2500 nested calls`, located at the def. A def and a top-level variable cannot share a name, even when the def comes from a `#! use` module: the compiled program has one namespace for both, so every engine refuses it and names the two lines. |
| `step name "text" { }` | A named phase of a def (or of the top level), with optional text. A step is a block, not a scope: its variables are the def's. The engines record entering and leaving it, and the explorer, the story and the intent and debug lenses address it as `def/name`, so its name is unique in the def. Not inside a `parallel for`. |
| `step name: statement` | A step of one statement, on that statement's line: `step attn_residual: X1 = X + Y` is `step attn_residual { X1 = X + Y }` to every engine and tool, hashes included. The statement starts after the `:` on the same line and ends where a statement ends, so `step bump: for i to n {` runs to the loop's closing brace, and an `else` on the next line continues an `if`. `;` is a line break, so `step a: x = 1; y = 2` holds only `x = 1`. With text: `step name "text": statement`. |
| `### … ###` | A fenced block of text: everything between two lines that are exactly `###` is skipped, so prose needs no `#` in front of it. A post (https://viv-lang.org/posts) keeps its title and its notes in one. |
| `rung "title" { }` | One step of a ladder (https://viv-lang.org/ladders): defs, tests and statements that build on the rungs before it. A def in a later rung replaces an earlier one of the same name. `viv run` runs the top rung; `viv ladder` runs every rung and checks the ladder's rules. |
| `record P { x: Num, v: Arr }` | Record type (a user-defined struct), in every engine. Build one with every field named, `p = P { x: 1, v: zeros(3) }`; read a field with `p.x`; assign one with `p.x = 2`; write into an array field with `p.v[i] = e`. A record is a value: `q = p` gives `q` its own, and `q.x = 3` leaves `p` as it was (an array inside is shared until it is written into, as arrays are). A write is into a variable of the def it is in: a def does not change a record it was given, or a global. |
| `record M { blocks: [B] }` | An array of records: `[b1, b2]`, `append(bs, b)`, `bs[i]`, `len(bs)`, `bs[i] = b`, and a place through it, `m.blocks[l].mlp.up = W` or `m.blocks[l].mlp.up[:] = W`. A value, as each record in it is. |
| `a + b`, `zeros_like(r)`, `for_leaves(f, r1, r2)` | Over a record of arrays (every field an array, a record of arrays, or an array of them): `+ - * /` with another of its type or a number, `-r`, `sqrt(r)` and the other elementwise functions, act on each array; `zeros_like(r)` is one of the same shape; `sum(r)` adds its arrays' sums; `for_leaves(f, r1, ..., x)` calls `f` on each array of the records at once, with `x` after them, and `f` may write them in place. `loss, g = value_and_grad(f, model, ...)` differentiates a record of arrays and gives `g` of its type. A number field is refused, with its name: keep sizes in a record of their own. |
| `fn(x) -> expr` | A lambda, and a value: assign it, pass it, call it through the parameter. |

### Output with `print` and `fmt`

`print` writes its arguments on one line with a space between them, as Python's does, then a newline. `sep=", "` puts other text between them and `sep=""` none; `sep` is a string literal or a variable. There is no `end=`: every `print` ends its line. An argument is a string (written as is), `fmt(x, ...)` (a value in a format), or any other value (as `show` writes it: six significant digits, the way C's `%g` does). Every value is computed before anything is written. The interpreter, C and wasm print the same bytes, checked over 10^6 doubles (`tools/print-fmt-sweep/`).

| Form | Writes | Python |
|---|---|---|
| `print("loss", loss)` | `loss 0.693147` | `print(f"loss {loss:g}")` |
| `print("loss", fmt(loss))` | `loss 0.6931471805599453`: the shortest digits that read back as the same double | `print("loss", loss)` (`repr`) |
| `print(a, b, c, sep=", ")` | `1, 2, 3` | the same |
| `print("[", x, "]", sep="")` | `[0.5]` | the same |
| `fmt(x, 17)` | 17 significant digits, C's `%.17g`; any count 1 to 17 | `{x:.17g}` |
| `fmt(x, ".6g")`, `".3e"`, `".2f"` | C's `%.6g`, `%.3e`, `%.2f` | the same spec |
| `fmt(x, "10.4f")` | at least 10 characters, padded on the left | the same spec |
| `fmt(n, "d")`, `"5d"` | a whole number, no decimal point (`%.0f`) | `{n}`, `{n:5d}` for an int |
| `print(v)`, `fmt(v, ".2f")` | an array on one line, `[1, 2.5]`, `[[1, 2], [3, 4]]`, each element in the format | — |
| `print(flag, z)` | `true 1+2i`: booleans and complex numbers, bare only | — |
| `print(x)` | one number: the one-argument `print`, unchanged | — |
| `print()` | an empty line | `print()` |

A format is `[width][.precision]` then `e`, `f`, `g` or `d`, or empty for `repr`; precision and width are at most 40 and 100. The rest of Python's spec (fill, alignment, `+`, `,`, `%`, a precision with no type) is refused with its location and what to write instead. The digit count and the format are written in the source, not computed. `nan` prints as `nan` whatever its sign, and `fmt` keeps the sign of a negative zero (`-0.0`, `-0`) where a bare value writes `0`. A program's own `def fmt` or `def print` replaces the builtin, and no separator is added to its calls.

Until 2026-09-30 `print` put nothing between its arguments, and programs wrote `print(a, " ", b)`. `tools/migrate_print_sep.py` rewrote every program in the repository to the new form, checked call by call to reach the engines unchanged, so each prints the bytes it printed before.

### Operator precedence

| Level | Operators | Notes |
|---|---|---|
| 1 | `x ^ y` | power, the same as `pow(x, y)` (so an array base works elementwise). Right-associative, and tighter than a minus on its left: `-x ^ 2` is `-(x ^ 2)`, `2 ^ 3 ^ 2` is 512, `2 ^ -1` is 0.5 |
| 2 | `-x` `not x` | unary |
| 3 | `*` `/` `@` `%` `//` | `*` on two arrays is elementwise; `@` is always the matrix product (NumPy's `matmul`: a stack of matrices multiplies matrix by matrix, under Axes (#arrays)). `a % b` is the floored remainder and `a // b` the floor of the quotient, Python's rules for floats: the remainder takes the sign of `b`, so `-7 % 3` is 2 and `(0 - 1) // 2` is -1, and `a == (a // b) * b + a % b`. They take numbers; a zero divisor stops the program with `division by zero in a % b` |
| 4 | `+` `-` | elementwise on arrays; all four arithmetic ops broadcast between an array and a scalar, either side, and between two (not complex) arrays by NumPy's rule (below (#arrays)) |
| 5 | `<` `<=` `>` `>=` `==` `!=` `in` `not in` | `x in a`: a number in a real array of any rank (some element `==` x, so a NaN is in nothing, as in NumPy), a key in a map, a substring in a string, a string in a list of strings. On numbers, the comparisons give a Bool. With a real array on either side they compare elementwise and give an array of 1.0 where the comparison holds and 0.0 where it does not, so `dY * (Y > 0.0)` is ReLU's adjoint on any rank. `all`, `any` and `where` read such a mask; `if a == b` on arrays is refused, with a note that says to write `all(a == b)` |
| 6 | `and` `or` | short-circuit in every engine: the right side runs only when the left does not decide, so `j > 0 and a[j - 1] > v` never reads `a[-1]` |
| 7 | `a if c else b` | a when c is true, else b; only the one chosen is computed. Loosest of all, as in Python: `a if c else b if d else e` is `a if c else (b if d else e)` |

Measured: `2.0 + 3.0 * 4.0` is 14, `-7 % 3` is 2 and `7 // -2` is -4. Parentheses override everything.

### Comprehensions, `if` expressions and keys

A comprehension builds an array, or a list of strings, from loops, as Python's does. Its `for` takes any head a `for` statement takes, and its `if` keeps the elements it holds for:

```viv
xs = [3, 1, 4, 1, 5, 9]
evens = [x for x in xs if x % 2 == 0]           # [4]
squares = [i * i for i to 5]                    # [0, 1, 4, 9, 16]
table = [[i * j for j from 1 to 4] for i from 1 to 4]
words = [upper(w) for w in split("a bb ccc")]   # ['A', 'BB', 'CCC']
print(" ".join(f"{x:.2f}" for x in xs[:3]))     # 3.00 1.00 4.00
print(all(x > 0 for x in xs), max(words, key=len))
sign = "neg" if xs[-1] < 0 else "pos"
print(sorted(words, key=len, reverse=true))
```

| Form | Gives |
|---|---|
| `[e for x in xs]`, `[e for i to n if c]` | a vector of numbers, a matrix when every `e` is a vector of one length (a row each), or a list of strings. More `for`s nest, the first outermost; `xs` may be an array, a string, a list of strings, `enumerate(...)`, `zip(...)` or `rows(M)`. An empty result is `[]` |
| `{k: v for x in xs}` | a map (#maps), Python's dict comprehension: `stoi = {c: i for i, c in enumerate(vocab)}`. A key seen again keeps its first place and takes the later value. The same heads and filters as a list's. A brace in an expression is only this: a map literal is `map()` and stores, and a set (`{x, y}`) is refused with a pointer to `unique` |
| `f(e for x in xs)` | a generator as the only argument of a call: `f([e for x in xs])`, so `sum`, `max` and `" ".join` take one |
| `all(c for x in xs)`, `any(...)` | a Bool, and the loop stops at the first element that decides it, as Python's does |
| `a if c else b` | a when c holds, else b; only the one chosen is computed |
| `max(xs, key=f)`, `min(xs, key=f)` | the first element whose key is largest (smallest), as Python's; a NaN key is never chosen after the first element; an empty `xs` stops with `max() arg is an empty sequence` |
| `sorted(xs, key=f, reverse=r)`, `argsort(xs, key=f, reverse=r)` | the elements (or their positions) in the order of their keys, numbers or strings: stable, so equal keys keep their order, with `reverse=true` too, as Python's; a NaN key goes after every number |

A key is a def's name, a builtin's (`len`, `abs`), a variable holding a function, or a lambda `fn(w) -> ...`, which may read the program's variables. A comprehension's variables are its own: `x` outside it is neither read nor changed, as in Python 3. The parser writes each of these forms as the loop it stands for, computed just before the statement that uses it, so every engine runs code it already runs; what a statement computes before it (a call of a def, a `print`) still runs first, and in a `while` condition it is computed again each time. For that reason one cannot be the body of a `fn(...) -> ...` (except a key's): write a def.

### Axes, broadcasting, and agreement with NumPy

The array builtins take NumPy's names and give NumPy's results: reshapes and transposes, `eye`, `diag`, `arange`, reductions along an axis, gathers and masks, `cumsum`, `argsort`, `einsum`, and solvers that return what they compute. Their tables are under Builtins (#builtins); the rules they share are here.

```viv
X = reshape(arange(12), 3, 4)
Xc = X - mean(X, 0)                 # the column means, broadcast down the rows
C = transpose(Xc) @ Xc / 3.0
vals, vecs = eig_sym(C)
w = lin_solve(C + eye(4), ones(4))
top = argsort(vals)
picked = take(X, [2.0, 0.0], 0)     # X[[2, 0], :]
big = extract(X > 5.0, X)           # X[X > 5]
```

- **Broadcasting**: `+ - * /`, the comparisons and two-array `max` and `min` between two arrays of different shapes (not complex ones) line up their last axes; each pair of lengths is equal, or one of them is 1 and repeats. So a row vector against a matrix repeats down its rows, `reshape(v, n, 1)` is a column that repeats across, and shapes that do not fit stop with the words they always did: `array shape mismatch: [2, 3] vs [3, 2]`. Complex arrays do not broadcast. `grad` sums a broadcast operand's gradient back to its shape.

- **Matrix products**: `A @ B` is NumPy's `matmul`. Two matrices give their product; a vector on the left is a row and on the right a column, and the axis it adds is dropped. An array of rank 3 or 4 is a stack of matrices in its last two axes, and the axes before them broadcast, so `X @ W` with `X` of shape (4, 8, 16) and `W` (16, 8) multiplies each of the four matrices by `W`. Two vectors are refused (`dot(u, v)` is their product). Each entry adds its terms in order of the inner index from 0.0, skipping a term whose left factor is zero, with a separate multiply and add (no fused multiply-add); the loop is blocked for the cache without changing that order, so every engine gives the same bits (`lang/viv/src/interp/matmul.rs`). NumPy calls BLAS, which sums in its own order and agrees to rounding.

- **Axes**: An axis is 0 to rank − 1, or −rank to −1 counted from the end, as in NumPy and as an index is. A reduction along an axis (`axis=`) drops it; with `keepdims=true` it keeps it with length 1, as NumPy's does, which is what lets `A - mean(A, axis=1, keepdims=true)` subtract each row's mean. A builtin's options are always given by name, with NumPy's names (`axis=`, `keepdims=`, `k=` for `triu`/`tril`/`diag`, `offset=` for `trace`, `endpoint=`, `axes=`); written by position they are refused. A reduction along the only axis of a vector would give a number, and is refused with the whole-array form to write instead (`mean(v)`).

- **Orders of summation**: `sum`, `mean`, `var` and `std`, along an axis or of a whole array, and `trace` add in NumPy's order: from 0.0, a contiguous run in blocks of 8192 entries, each block by NumPy's pairwise sum; along an axis that is not the last, the slabs in order. They agree with NumPy to the bit (`tools/py2viv/probe_numpy_order.py` states the order), so `mean(a)` is `sum(a) / numel(a)` to the bit. The order binds the exact engines; fast mode may add in another. `prod`, `cumsum` and `cumprod` go left to right, as NumPy does. `dot` and `norm` go left to right; NumPy calls BLAS, which agrees to rounding. `einsum` loops over the output entries in row-major order and over the summed indices in their order of first appearance, adding each product to 0.0; NumPy's einsum has its own loops, so the two agree to rounding.

- **NaN and signed zeros**: `amax`, `amin`, `argmax` and `argmin` are NumPy's: a NaN wins, and `amax` puts 0.0 above −0.0. (`max(a)` of one array is older and skips NaN, as C's `fmax` does.) `argsort` puts NaN last and keeps −0.0 and 0.0 in their order. `clip` keeps a NaN entry and treats a NaN bound as no bound.

- **Masks and indices**: A mask is an array of 1.0 and 0.0, as a comparison gives, so an array cannot say whether it holds truth values or indices. NumPy's `a[mask]` is `extract(mask, a)`, `a[mask] = v` is `putmask(a, mask, v)`, `a[idx]` along an axis is `take(a, idx, axis=k)`, and `np.flatnonzero` is `flatnonzero`. The indices `take` and `put` read are not negative: `a[-1]` counts from the end, `take(a, [-1.0], 0)` stops.

- **Solvers that return**: `x = lin_solve(A, b)` (a matrix `B` solves for every column), `inv(A)`, `L = cholesky(A)`, `vals, vecs = eig_sym(A)`, `w, V = eig_herm(A)` and `U, s, Vt = svd(A)` run the same code as the forms that fill arrays you pass, which still work, and give the same bits. The compiler writes an unpacking of `eig_sym`, `eig_herm` or `svd` as that call with arrays it allocates. `svd` is one-sided Jacobi (`lib/svd.viv`'s algorithm); its vectors can have other signs than LAPACK's.

`tools/arrays/numpy_cases.py` writes 4,812 cases whose results NumPy computes in an order Viv follows, and 89 that agree with NumPy within a tolerance (`einsum`, the solvers, `dot`); `cargo test --test arrays` runs them in the interpreter, C and wasm, which print NumPy's results byte for byte.

### Types

- **Num**: f64, and the default: a literal is a Num. The other scalars are declared: `Int` (a 64-bit integer, `n: Int`), `Float32` (IEEE binary32) and `Cplx` (below). A double holds every integer up to 2^53 exactly, so counters, indices and `n == 1.0` are exact.

- **Bool**: `true` or `false`: what a comparison yields and what `if` and `while` take. It can be stored, passed and returned (`def f(x: Num) -> Bool`), combined with `and`, `or`, `not`, and compared with `==` and `!=`. It is **not a number**: `true + 1.0` and `if count {` are refused, identically, by the interpreter and the compiler. `show` prints it as `true` or `false`. Arrays hold numbers, so a mask is still an array of 1.0 and 0.0: a comparison with an array on either side makes one, and `all(m)` and `any(m)` turn one back into a Bool. See booleans (https://viv-lang.org/docs/booleans).

- **Str**: An immutable UTF-8 string, indexed by code point (Strings (#strings)). A list of them is `[Str]`.

- **Arr**: Rank 1 to 4, row-major, f64 elements. Arrays pass into functions and return from them, and scalars broadcast against them from either side, and two arrays (not complex ones) broadcast as NumPy's do (Axes, broadcasting, NumPy (#arrays)). **Assignment and parameter passing alias:** `b = a` shares storage, and a function that writes into an array parameter writes into the caller's array. A real copy is `copy(a)`, and a slice (#slices) is a copy too — see scope & mutation (https://viv-lang.org/docs/scope). Indices count from 0, and an index outside its dimension stops the run in the interpreter and in C with the same message, `index 2 out of bounds for dimension 0 (size 2)`, and the line it happened on. A negative index counts from the end, as in Python: `a[-1]` is the last element, and `a[-3]` of two elements stops the run. A fractional index is cut toward zero, so `a[1.7]` is `a[1]`.

- **Map**: Numbers under number keys, in the order the keys were first stored, made by `map()`. A handle: assignment, parameters and record fields share one map. A parameter may be annotated `m: Map`. See Maps (#maps).

- **Cplx**: A first-class complex scalar: the `i` suffix on a numeric literal (`2.0 + 3.0i`), full arithmetic with real promotion, and complex overloads of `abs/exp/log/sqrt/sin/cos/tan`. Compiled to C99 `double complex`. Every engine computes these the same way, to the bit, with Viv's own formulas rather than C's operators or `<complex.h>`: A real operand is lifted to `(x, 0.0)` first, so `2.0 * z` is the complex product with its zero terms. `+` and `-` work part by part. `a * b` is `(ar br − ai bi) + i (ar bi + ai br)`, each product rounded on its own (no fused multiply-add). `a / b` is Smith's algorithm (1962): divide through by the larger part of `b`. If `|br| ≥ |bi|`, then with `r = bi/br` and `d = br + bi r` the quotient is `((ar + ai r)/d, (ai − ar r)/d)`; otherwise `r = br/bi`, `d = bi + br r`, and the quotient is `((ar r + ai)/d, (ai r − ar)/d)`. It never squares `b`, so it overflows only where the quotient does, not once `|b|` passes 1e154 as the textbook `a conj(b) / |b|²` does. A real divisor divides each part once, so `(1 + 2i) / 3` is exactly `1/3 + (2/3)i`. Dividing by zero gives `nan + nan i`. `abs(z)` is `hypot(re, im)`. `exp` is `e^re (cos im + i sin im)`. `log` is `ln|z| + i atan2(im, re)`. `sqrt` is `√|z|` at half the argument. `sin` and `cos` use the `cosh`/`sinh` formulas, and `tan` is `sin z / cos z` by the division above. All of them use Viv's libm. A complex array divided by a scalar divides each element the same way.

- **CArr**: An array of complex numbers, rank 1 to 4: `czeros(n)`, or a literal with a complex element in it (`[1.0 + 2.0i, 3.0]` — the reals promote). Indexing gives a `Cplx`; storing a real promotes it. It is **packed**: interleaved (re, im) doubles, which is exactly how C99 lays out `double complex`, so both engines walk the same memory and `fft` agrees between them to the byte. A complex and a real array are different types and do not mix silently — `sin(a)` on a complex array is an error. Annotate a parameter as `CArr` when you want to say so; otherwise it is inferred from the call site.

- **Records**: User-defined structs with typed fields, copied by value: the C backend compiles them to C structs, the wasm backend to a block of memory, and the interpreter holds the fields. Array fields refer to the original array, so bundling a mesh into one value copies no array data. A def that takes or returns one annotates it (`def f(p: P) -> P`). Records have no operators and no `show`: compare and print their fields. The several results of `return a, b` are a record the compiler declares for the def, with one field per value and the field types inferred from the returns; the program never names it, because the caller unpacks it.

- **Function values**: A `def` name or a lambda is a value. The compiled backends implement higher-order calls by monomorphization — `run_it(sq, x)` specializes to `run_it__sq` — with `viv_fn` function pointers for the cases where a function is held in a variable. Arrays of functions work in the interpreter and are a compile error in C.

### What a test prints

A `test` block is the unit that passes or fails, so it prints exactly one line under the name you gave it — not one line per assertion:

```viv
test "simpson integrates exp correctly" {
    assert_near(simpson(fn(x) -> exp(x), 0.0, 1.0, 60), exp(1.0) - 1.0, 1.0e-9)
}

$ viv run integrate.viv
PASS simpson integrates exp correctly
```

Several asserts in one block still give one line. A failure names the test and carries the first failing assertion's detail, and a block with no assertions says so instead of printing PASS:

```
PASS two asserts, one verdict
FAIL a failing one -- assert_near: |1 - 2| = 1 > 1e-09
PASS nothing asserted (no assertions)
assert_near PASS
1 assertion FAILED
```

That last `assert_near PASS` is an assertion at top level, outside any `test` — there is no name to report under, so it prints its own line as before. The trailing count is the summary, and it is what makes the process exit nonzero, so a build script can fail on the exit code. It goes to stderr, after everything the program printed, in the interpreter and in a compiled program alike, so stdout holds only the program's own output.

### A note on contracts

`integrate` and `find_root` are not builtins. They are **contracts**: functions named for the mathematical problem rather than the method, backed by several methods from different families of algorithms. By default they run every applicable method and check that the answers agree, and a disagreement stops the program with an error. `viv run --release` runs only the one method the contract selects (Gauss quadrature for `integrate`, a safeguarded secant for `find_root`).

The methods are ordinary Viv `def`s in `lang/viv/std/contracts.viv`, brought in with `viv --use lang/viv/std/contracts.viv program.viv`. The methods a contract chooses among are Viv code, not code in the compiler.

**Interpreter only, today.** A contract picks its method out of an array of functions, and the C backend implements higher-order calls by monomorphization — it specializes the callee per call site — which cannot express a callee chosen at runtime. So the contract library runs in the interpreter, and the C backend stops with a compile error naming the feature instead of producing wrong code.

### The backend matrix

The interpreter and the two compilers share one parser. Which features each one runs:

| Feature | Interpreter | C backend | Direct wasm |
|---|---|---|---|
| Scalars, control flow, `def`, recursion | ✓ | ✓ | ✓ |
| Arrays | ✓ | ✓ | ✓ |
| Lambdas, `simpson`, `deriv` | ✓ | ✓ | ✓ |
| `grad`, `value_and_grad` (rewritten before any engine) | ✓ | ✓ | ✓ |
| Records | ✓ | ✓ | ✓ |
| `return a, b`, `a, b = e` (rewritten before any engine) | ✓ | ✓ | ✓ |
| Complex scalars | ✓ | C99 `double complex` | two f64s, the interpreter's formulas |
| Arrays of functions (for example, a list of methods to try) | ✓ | ✓ | a table, `call_indirect` |
| `parallel for` | sequential | OpenMP pragma | compiles, sequential |
| `narrate` | `narrate:` line | `narrate:` line | `narrate:` line |
| Figures (`plot` … `scene`) | drawn | printed as records | printed as records |
| Call graph (`__graph__`/`__trace__`) | ✓ | — | — |
| Runaway-loop guard | step budget | — | — |
| Nonzero exit on a failed assert | ✓ | ✓ (atexit) | ✓ (the host's) |

Every feature goes into both the interpreter and the C backend, and the example suite checks that it does. The direct wasm backend reaches wasm without C (a sum-of-squares program is a -byte module) and is checked beside them as a third engine (`tools/three-engines.sh`).

### The CLI

| Command | What it does |
|---|---|
| `viv f.viv -o f.c` | Compile to portable C99. Then `cc f.c -O2 -lm`. |
| `viv f.viv -o f.wasm` | Compile straight to a wasm module, no C step. |
| `viv f.viv --emit tokens\|ast\|ir` | Dump an intermediate representation. The internals page (https://viv-lang.org/internals) is these dumps. |
| `viv run [--release] f.viv` | Interpret. The compiler's checks run first, so `viv run` refuses what the compile refuses, with the same messages: a call to an unknown function in a branch the run never takes, or a variable that holds a number in one place and an array in another. `--release` runs only the method each contract selects instead of all of them. |
| `viv run --quiet f.viv` | Drop the `__graph__`/`__trace__`/`__steps__`/`__story__` lines. |
| `viv run --json f.viv``viv f.viv --json` | Write the errors to stderr as one line of JSON instead of text: `{"diagnostics":[...]}`, each with `severity`, `message`, `file`, `line`, `col`, `end_line`, `end_col`, `label`, `source` (the line) and `notes`. A run with no errors writes an empty list. For tools and models that read the errors and try again. |
| `viv --use lang/viv/std/contracts.viv f.viv` | Put another file (a prelude) in front of this one before parsing. A program names the modules it uses itself, with import lines. |
| `import nn` &middot; `import hf.scf as scf` &middot; `from nn import rmsnorm, softmax_rows` | Python's forms with Python's bindings: `nn.rmsnorm(X)`, `scf.run(...)`, `rmsnorm(X)`; `import a.b` binds `a`; a `from` list may be parenthesized over several lines. Module `a.b` is `a/b.viv` beside the file, then `lib/a/b.viv`, then `a/b.viv` under the repository root; two different files that match are an error. A module is spliced in once, its tests left out, and its names are its own: a program sees what its import lines bind and not what its imports import, and a name two modules define is two names. Its phrases come with it. `import *` and relative imports are refused, and so is `#! use`, with the import to write. Resolved before anything parses: no run-time cost. |
| `private def f` &middot; `private N = 3` | Keep a top-level def, record, enum or constant in its module: importing it, or reading it from another module, is a compile error. |
| `import post "slug" as name` | A post used as a module. |
| `#! if-used basis.sto3g` &hellip; `#! end-if` | Keep the lines between only when the program imports that module, anywhere in its imports. A lookup can serve whichever of several modules a program loads: `basis.named` finds a basis set by name among the tables the program imports. |
| `#! if-none-used basis.a basis.b` &hellip; `#! end-if` | Keep the lines between, imports included, when the program's other imports reach none of the modules named: a default the program overrides by importing one. `basis.named` loads all eight basis tables for a program that imports none. |
| `viv expand f.viv` &middot; `viv imports --write f.viv` &middot; `viv fix --keywords --write f.viv` | The program as the engines read it, imports spliced and names resolved; rewrite a file's `#! use` lines as imports; name the builtin options a file writes by position. |
| `viv doc f.viv [fn]` | Structured annotations as JSON, or one function's doc card. |
| `viv whatis f.viv name` | One name explained: kind, type, where it came from. |
| `viv lsfunc f.viv` | Every `def` with its parameters. |
| `viv lsvar f.viv [fn]` | Every binding, with its glossary line if it has one. |
| `viv src f.viv fn` | One function's source, annotations stripped. |
| `viv deps f.viv` | The static call graph. |
| `viv canon f.viv` | Canonical dependency order. Idempotent. |
| `viv hash f.viv` | Per-definition content addresses: node hash and body hash. |
| `viv features [name]``viv features --usage [paths]` | What Viv can do, by area, each entry with what it does, when to use it, an example and its output; `--usage` counts each feature's uses in a corpus and lists code that does by hand what a feature does. See the feature catalog (#features). |
| `viv diff a.viv b.viv` | Structural diff — edits, renames, unchanged. |

An error at run time reads like one at compile time: the message, the file, line and column, the line itself, and a caret under the part that failed. A compiled program prints the same lines, byte for byte, because the compiler writes each indexing site's location into the C.

The notes answer habits from Python, `**` (write `x ^ y`) and `for i in range(n)` (write `for i to n`), and one from Julia, `2x` (write `2 * x`).

Every one of these commands also runs in the browser, on this site, at /cli (https://viv-lang.org/cli) — the wasm module exports the same library functions `main.rs` wraps.

### The C runtime

Compiled programs carry their runtime with them. One struct,

```c
typedef struct { int rank; int dims[8]; int size; double* data; } viv_arr;
```

plus a set of small `viv_*` functions — array ops, the eigensolvers, the LU and Cholesky factorizations, quadrature, the Park–Miller generator, and the show/plot/narrate printers — emitted **only on demand**. A program that never touches an array gets none of them; the nine-line example on the internals page (https://viv-lang.org/internals) pulls in exactly one. Building it takes only `cc program.c -O2 -lm`.

Every array read and write in the C is checked: the number of indices against the array's rank, then each index against its dimension. A wrong count or an index out of bounds stops a compiled program with the interpreter's error and location instead of reading memory the array does not own. The checks are on by default. Build with `-DVIV_RELEASE` to drop them along with the other debug-build checks, or with `-DVIV_NO_BOUNDS` to drop only them. What they cost, measured with `tools/bounds-cost.sh` (CPU time, best of 7 runs, `cc -O2`, Apple M2 Pro):

| Program | Checked | Unchecked | Cost |
|---|---|---|---|
| `bench/hf_md_eri_x100.viv` (water 6-31G ERIs, built 100 times) | 0.385 s | 0.275 s | +40% |
| `bench/matmul_loops.viv` (400×400, all indexing) | 0.105 s | 0.082 s | +29% |
| `scf/meta_integrals.viv` | 0.503 s | 0.396 s | +27% |
| `scf/meta_scf.viv` | 0.660 s | 0.368 s | +79% |
| `scf/lock_water.viv` | 0.200 s | 0.115 s | +73% |
| `scf/coverage.viv` | 0.260 s | 0.145 s | +80% |
| `bench/hf_scf_bench.viv` (water cc-pVDZ RHF) | 0.042 s | 0.023 s | +84% |

That is 27% to 84% on these programs, so a build you time for speed should pass `-DVIV_RELEASE`; a build you run to check a program should keep them. The WebAssembly modules this site serves (`dft/wasm/`, `engines/wasm/`) are built with `-DVIV_RELEASE`; the same programs run with the checks on in `dft/check.sh` and in the interpreter and in C. In those modules an index out of range is not caught: a write can land in another array of the same program (it cannot leave the module's memory). The interpreter and the module the compiler in your browser makes when you press Run (the direct wasm backend) check every index and stop with the same error as the C. `-DVIV_RELEASE` also trusts the `disjoint` pins of a memory lens instead of testing them; a checked build records which pass of a pinned `parallel for` writes each element and stops if two passes write the same one.

The location in a run-time error names the program as you gave it when the path is relative (`lib/hf/scf.viv:12:5`) and by its file name when the path is absolute (`scf.viv:12:5`), in the interpreter and in C. The C carries those locations, so a program compiled by absolute path gives the same C on every machine and in every checkout.

The builtin tables were probed at build time ().

## Python habits that do not work

| Habit | Python-style input | What Viv says | Write instead |
|---|---|---|---|
| power | `x = 2.0 ** 3.0` | error: expected an expression | `x = 2.0 ^ 3.0` |
| range loop | `for i in range(3) { ⏎ print(i) ⏎ }` | error: expected `to` | `for i to 3 { ⏎ print(i) ⏎ }` |
| colon blocks | `def f(x): ⏎ return x * x` | error: expected `{` | `def f(x) { ⏎ return x * x ⏎ }` |
| elif | `x = 1.0 ⏎ if x > 2.0 { ⏎ print(1.0) ⏎ } elif x > 0.0 { ⏎ print(2.0) ⏎ }` | error: a brace in an expression is a map comprehension, `{k: v for x in xs}` | `x = 1.0 ⏎ if x > 2.0 { ⏎ print(1.0) ⏎ } else if x > 0.0 { ⏎ print(2.0) ⏎ }` |
| imports | `import numpy as np` | error: line 1: `numpy`: no module by that name (looked for ./numpy.viv, lib/numpy.viv, numpy.viv) | `# no numpy to import: every builtin is always in scope; a module of lib/ is imported as Python does it, from nn import rmsnorm` |
| np. prefix | `a = np.zeros(3)` | error: undefined variable `np` | `a = zeros(3)` |
| chained indexing | `A = zeros(2, 2) ⏎ show A[1][1]` | error: array has rank 2 but was indexed with 1 index/indices | `A = zeros(2, 2) ⏎ show A[1, 1]` |
| .T and methods | `A = zeros(2, 2) ⏎ B = A.T` | error: cannot access field `T` of a array | `A = zeros(2, 2) ⏎ B = transpose(A)` |
| list.append | `a = [] ⏎ a.append(1.0)` | error: `append` gives back the longer array, so it is written `a = append(a, x)`, with the same name on both sides | `a = [] ⏎ a = append(a, 1.0)` |
| dict | `d = {} ⏎ d["a"] = 1.0` | error: a brace in an expression is a map comprehension, `{k: v for x in xs}` | `d = map() ⏎ d["a"] = 1.0` |
| a def name to simpson | `def sq(x) { ⏎ return x * x ⏎ } ⏎ show simpson(sq, 0.0, 1.0, 10)` | error: `simpson`'s first argument must be a `fn(x) -> ...` lambda, not the def `sq` (wrap it: `fn(x) -> sq(x)`) | `def sq(x) { ⏎ return x * x ⏎ } ⏎ show simpson(fn(x) -> sq(x), 0.0, 1.0, 10)` |
| lambda | `f = lambda x: x * x` | error: expected `=` after the type | `f = fn(x) -> x * x` |
| True / False / None | `x = True` | error: undefined variable `True` | `x = true` |

- `show name` prints `name: value`; `show` of any other expression prints the bare value, and so does `print`. Numbers print with 6 significant digits.
- Unary minus works anywhere: `-1.0`, `exp(-t)`, `[-1.0, 2.0]`. (Some examples write `0.0 - x`; both forms are the same.)
- The solvers return their results: `x = lin_solve(A, b)`, `L = cholesky(A)`, `vals, vecs = eig_sym(A)`, `U, s, Vt = svd(A)`. (The older forms that fill arrays you allocate, `eig_sym(A, vals, vecs)`, still work.)
- NumPy's array functions have their NumPy names and results: `reshape(a, 2, 3)`, `sum(A, 0)`, `mean(A, 1, true)` (keepdims), `argmax`, `eye`, `outer`, `einsum("ij,jk->ik", A, B)`, and two arrays broadcast. `A.T` is `transpose(A)`, `a[mask]` is `extract(mask, a)`, `a[idx]` is `take(a, idx, 0)`, `np.max(a, axis=0)` is `amax(a, 0)` (`max(a, b)` is elementwise).
- Set a test's tolerance from an error estimate, or run the program once and read the printed error first. In a trial of this pack (five programs written from it) the one failure was a tolerance guessed too tight (Simpson's rule with 8 panels on exp over [0, 1] is off by 2.3e-6, and the test asked for 1e-6).
- A program should pass in both engines: `viv run f.viv` and `viv f.viv -o f.c && cc -O2 f.c -lm && ./a.out` print the same lines. The interpreter accepts a few things the compiler refuses; section 3 lists them.

## The interpreter runs these; the C compiler refuses them

**deriv of a function parameter.** deriv, simpson and trapezoid take a lambda written at the call. Write the function inside the deriv call, or pass the derivative in as a second function parameter. (deriv cannot see inside a call to a def either.)

```viv
def newton(f, x0) {
    x = x0
    for k to 20 {
        x = x - f(x) / deriv(f, x)
    }
    return x
}
show newton(fn(x) -> x * x - 2.0, 1.0)
```

`viv f.viv -o f.c` says: `` `deriv`'s first argument must be a `fn(x) -> ...` ``

Instead:

```viv
def newton(x0) {
    x = x0
    for k to 20 {
        x = x - (x * x - 2.0) / deriv(fn(t) -> t * t - 2.0, x)
    }
    return x
}
show newton(1.0)
```

## Builtins, with the engines that run each

### Real math

Every one-argument function here also maps elementwise over an array, so `y = sin(x)` needs no loop.

- `exp(x)`: e^x [interpreter yes, C yes, direct wasm yes]
- `log(x)`: natural log [interpreter yes, C yes, direct wasm yes]
- `ln(x)`: alias for log [interpreter yes, C yes, direct wasm yes]
- `sqrt(x)`: square root [interpreter yes, C yes, direct wasm yes]
- `pow(x, y)`: x to the y [interpreter yes, C yes, direct wasm yes]
- `abs(x)`: magnitude [interpreter yes, C yes, direct wasm yes]
- `fabs(x)`: alias for abs [interpreter yes, C yes, direct wasm yes]
- `sin(x)` [interpreter yes, C yes, direct wasm yes]
- `cos(x)` [interpreter yes, C yes, direct wasm yes]
- `tan(x)` [interpreter yes, C yes, direct wasm yes]
- `asin(x)` [interpreter yes, C yes, direct wasm yes]
- `acos(x)` [interpreter yes, C yes, direct wasm yes]
- `atan(x)` [interpreter yes, C yes, direct wasm yes]
- `atan2(y, x)`: quadrant-correct angle [interpreter yes, C yes, direct wasm yes]
- `sinh(x)` [interpreter yes, C yes, direct wasm yes]
- `cosh(x)` [interpreter yes, C yes, direct wasm yes]
- `tanh(x)` [interpreter yes, C yes, direct wasm yes]
- `erf(x)`: error function [interpreter yes, C yes, direct wasm yes]
- `gamma(x)`: Γ, C’s tgamma [interpreter yes, C yes, direct wasm yes]
- `floor(x)` [interpreter yes, C yes, direct wasm yes]
- `ceil(x)` [interpreter yes, C yes, direct wasm yes]
- `round(x)` [interpreter yes, C yes, direct wasm yes]
- `min(a, b) · min(a, b, c, ...)`: the smallest; three or more arguments are min(min(a, b), c) [interpreter yes, C yes, direct wasm yes]
- `max(a, b) · max(a, b, c, ...)`: the largest, as min [interpreter yes, C yes, direct wasm yes]
- `hypot(x, y)`: sqrt(x^2 + y^2) with no overflow on the way (musl's, as abs of a complex number) [interpreter yes, C yes, direct wasm yes]
- `copysign(x, y)`: |x| with the sign of y; -0.0 is negative [interpreter yes, C yes, direct wasm yes]
- `fma(a, b, c)`: a*b + c rounded once (IEEE 754 fusedMultiplyAdd), the same bits in every engine and every build [interpreter yes, C yes, direct wasm yes]
- `fsum(xs)`: the correctly rounded sum, as Python's math.fsum (Shewchuk's partials) [interpreter yes, C yes, direct wasm yes]
- `mod(a, b)`: C’s fmod: the result takes the sign of a. The operator `a % b` is the floored remainder, with the sign of b [interpreter yes, C yes, direct wasm yes]

### Bits

These take a number (a double), truncate it toward zero, wrap it into [0, 2^32) and hand back the result as an unsigned 32-bit value, so they compose without a cast. A shift of 32 or more gives 0, except that `sar` of a negative value gives all 32 bits set (4294967295, which `i32` reads as -1). These are the operations a logic gate, an adder or a RISC-V instruction needs.

- `band(a, b)`: bitwise and [interpreter yes, C yes, direct wasm yes]
- `bor(a, b)`: bitwise or [interpreter yes, C yes, direct wasm yes]
- `bxor(a, b)`: bitwise exclusive or [interpreter yes, C yes, direct wasm yes]
- `bnot(a)`: all 32 bits flipped [interpreter yes, C yes, direct wasm yes]
- `shl(a, k)`: shift left by k [interpreter yes, C yes, direct wasm yes]
- `shr(a, k)`: shift right by k, zeros in [interpreter yes, C yes, direct wasm yes]
- `sar(a, k)`: shift right by k, sign in [interpreter yes, C yes, direct wasm yes]
- `u32(x)`: the low 32 bits, unsigned [interpreter yes, C yes, direct wasm yes]
- `i32(x)`: the same bits, signed [interpreter yes, C yes, direct wasm yes]

### Ints

Int is a 64-bit integer. Overflow, `//` and `%` are checked and stop the run with the same words in every engine; these builtins are the unchecked forms an algorithm on 64-bit patterns needs (a random generator, a hash), with the pattern read as unsigned where the name says so.

- `int_zeros(n)`: an [Int] array of n zeros [interpreter yes, C yes, direct wasm yes]
- `wrap_add(a, b) · wrap_sub(a, b) · wrap_mul(a, b)`: the result mod 2^64, as a pattern: no overflow check [interpreter yes, C yes, direct wasm yes]
- `umulhi(a, b)`: the high 64 bits of the 128-bit unsigned product of the two patterns [interpreter yes, C yes, direct wasm yes]
- `ult(a, b)`: a < b with both read as unsigned [interpreter yes, C yes, direct wasm yes]
- `ustr(a)`: the pattern as an unsigned decimal string [interpreter yes, C yes, direct wasm yes]

### Binary data

A file's bytes are an [Int], one element per byte, 0 to 255. The readers take that array and an offset and decode the bytes there, little-endian, as Python's struct with `<` does; an offset past the end, or an element that is not a byte, stops the run. In the browser `read_bytes` reads the files the page attaches to the run.

- `read_bytes(path) · read_bytes(path, offset, count)`: a file's bytes as an [Int]; with offset and count, only those (count left out: to the end) [interpreter yes, C yes, direct wasm yes]
- `le_u8(b, i) · le_i8(b, i) · le_u16(b, i) · le_i16(b, i) · le_u32(b, i) · le_i32(b, i)`: the Int that the 1, 2 or 4 bytes at offset i encode, unsigned or two's complement [interpreter yes, C yes, direct wasm yes]
- `le_u64(b, i) · le_i64(b, i)`: the Int that the 8 bytes at i encode; le_u64 stops when the value is above 2^63 - 1 [interpreter yes, C yes, direct wasm yes]
- `le_f16(b, i) · le_bf16(b, i) · le_f32(b, i) · le_f64(b, i)`: the number an IEEE half, a bfloat16, a single or a double at i holds, exactly (subnormals, infinities and NaN included) [interpreter yes, C yes, direct wasm yes]
- `fround(x)`: x rounded to the nearest binary32 (ties to even), as a number: JavaScript's Math.fround, NumPy's float64(float32(x)); a format that computes in single precision rounds each step with it [interpreter yes, C yes, direct wasm yes]

### Arrays

- `zeros(d1 … d4)`: allocate, rank 1 to 4 [interpreter yes, C yes, direct wasm yes]
- `len(a)`: the leading dimension [interpreter yes, C yes, direct wasm yes]
- `numel(a)`: the number of elements, every axis multiplied [interpreter yes, C yes, direct wasm yes]
- `shape(a)`: the lengths of the axes, as a vector [interpreter yes, C yes, direct wasm yes]
- `linspace(lo, hi, n)`: n evenly spaced points [interpreter yes, C yes, direct wasm yes]
- `fill(a, v)`: set every element in place [interpreter yes, C yes, direct wasm yes]
- `sort(a)`: in place, ascending, stable; -0.0 before 0.0, NaN last; no value (copy first for a sorted copy) [interpreter yes, C yes, direct wasm yes]
- `transpose(A)` [interpreter yes, C yes, direct wasm yes]
- `sum(a)`: the sum of every element of a real array, in NumPy's order in every exact engine: the entries in row-major order, from 0.0, in blocks of 8192, each block by NumPy's pairwise sum (below 8 entries left to right; to 128, eight running sums combined as a tree; above, two halves, the first a multiple of 8 long). Fast mode is not bound by this order [interpreter yes, C yes, direct wasm yes]
- `histogram(a, lo, hi, nbins)`: counts of a's entries in nbins equal bins over [lo, hi], the last bin closed on the right, entries outside not counted. Draw it with bars() [interpreter yes, C yes, direct wasm yes]
- `roll(a, k) · roll(A, k, axis)`: circular shift: element i moves to i + k, wrapping. A periodic stencil in one line: roll(u, 1) - 2.0 * u + roll(u, -1) [interpreter yes, C yes, direct wasm yes]
- `tridiag_solve(sub, diag, super, rhs)`: the Thomas algorithm: solves a tridiagonal system in O(n), for the diagonally dominant systems that implicit PDE schemes produce. sub[0] and super[n-1] are unused [interpreter yes, C yes, direct wasm yes]
- `max(a, s) · max(a, b)`: elementwise, against a scalar (either side) or another array; `max(rho, 1.0e-30)` floors a density [interpreter yes, C yes, direct wasm yes]
- `min(a, s) · min(a, b)`: elementwise, as max [interpreter yes, C yes, direct wasm yes]
- `max(a) · min(a)`: with ONE array: its largest / smallest entry. `max(abs(u - exact))` is the max-norm error of a solve [interpreter yes, C yes, direct wasm yes]
- `max(xs, key=f) · min(xs, key=f)`: the first element whose key is largest (smallest), as Python's; xs an array, a list of strings, a string or rows(M); f a def, a builtin or fn(x) -> ... [interpreter yes, C yes, direct wasm yes]
- `sorted(xs, key=f, reverse=r)`: a new array or list in the order of the keys (numbers or strings), stable, as Python's; a NaN key last [interpreter yes, C yes, direct wasm yes]
- `argsort(xs, key=f, reverse=r)`: the positions of xs in the order sorted(xs, key=f) puts them [interpreter yes, C yes, direct wasm yes]
- `all(m)`: true when every entry of a real array is nonzero, as in the mask `a == b` gives; true for an empty array [interpreter yes, C yes, direct wasm yes]
- `any(m)`: true when some entry is nonzero; false for an empty array [interpreter yes, C yes, direct wasm yes]
- `where(m, a, b)`: a where m is nonzero, b where it is 0. a and b are each an array or a number, and the three broadcast to one shape, as NumPy's `where`: `where(x > 0.0, x, 0.0)` [interpreter yes, C yes, direct wasm yes]
- `pow(a, s)`: every entry of a real array to the same scalar power: `pow(rho, 1.0 / 3.0)` [interpreter yes, C yes, direct wasm yes]
- `slab(A, k)`: a COPY of the k-th slab along the leading axis: row k of a matrix, plane k of a rank-3 array. Real or complex [interpreter yes, C yes, direct wasm yes]
- `setslab(A, k, S)`: overwrite that slab with S, in place [interpreter yes, C yes, direct wasm yes]
- `take(A, idx)`: gather: A viewed flat (row-major), idx a real array of flat indices; the result has idx's shape [interpreter yes, C yes, direct wasm yes]
- `put(A, idx, V)`: scatter, in place: the inverse of take [interpreter yes, C yes, direct wasm yes]
- `a[lo:hi:s] · A[i, :] · A[:, j]`: a COPY of a block: lo to hi - 1 in steps of s, any part left out; a plain index drops its axis. Real or complex. See Slices above [interpreter yes, C yes, direct wasm yes]
- `a[lo:hi] = v · A[:, j] = v`: write a number, or an array of the slice's shape, into the block, in place [interpreter yes, C yes, direct wasm yes]
- `copy(a)`: a new array with a's shape and elements (assignment aliases; this does not) [interpreter yes, C yes, direct wasm yes]
- `concat(a, b, ...)`: vectors joined end to end [interpreter yes, C yes, direct wasm yes]
- `vcat(A, B, ...)`: rows stacked; a vector is one row [interpreter yes, C yes, direct wasm yes]
- `hcat(A, B, ...)`: matrices side by side; a vector is one column [interpreter yes, C yes, direct wasm yes]
- `[]`: an empty vector, the start of an array grown by append [interpreter yes, C yes, direct wasm yes]
- `a = append(a, x)`: a with x on the end: a number onto a vector, or a vector as a new last row of a matrix ([] takes a first row). Written rebinding a; amortized O(1). See Growable arrays above [interpreter yes, C yes, direct wasm yes]
- `a = extend(a, b)`: vector a followed by vector b, grown like an append [interpreter yes, C yes, direct wasm yes]

### Shapes and constructors

NumPy's names and results (interp/arrays.rs; each case checked against NumPy in tools/arrays/numpy_cases.py). Every result is a new array. An axis may be negative, as in NumPy: -1 is the last.

- `reshape(a, d1 … d4)`: the same entries in row-major order under a new shape; one dimension may be -1 and is worked out [interpreter yes, C yes, direct wasm yes]
- `ravel(a) · flatten(a)`: every entry as a vector, row-major [interpreter yes, C yes, direct wasm yes]
- `transpose(a) · transpose(a, axes=axes)`: the axes reversed (a vector is copied); with axes, axis d of the result is axis axes[d] of a [interpreter yes, C yes, direct wasm yes]
- `swapaxes(a, i, j)`: two axes exchanged [interpreter yes, C yes, direct wasm yes]
- `eye(n) · eye(n, m) · eye(n, m, k)`: ones on diagonal k (k > 0 above the main one) [interpreter yes, C yes, direct wasm yes]
- `diag(v) · diag(v, k=k) · diag(A) · diag(A, k=k)`: a vector onto diagonal k of a square matrix, or diagonal k of a matrix as a vector [interpreter yes, C yes, direct wasm yes]
- `ones(d1 … d4) · zeros_like(a) · ones_like(a)`: filled arrays [interpreter yes, C yes, direct wasm yes]
- `arange(stop) · arange(start, stop) · arange(start, stop, step)`: NumPy's arange: ceil((stop - start) / step) entries, entry i start + i * ((start + step) - start) [interpreter yes, C yes, direct wasm yes]
- `linspace(lo, hi, n, endpoint=endpoint)`: NumPy's linspace: i * step + lo, with the end point included or left out; n may be 0 or 1 [interpreter yes, C yes, direct wasm yes]
- `rand(d1 … d4) · randn() · randn(d1 … d4)`: draws from rand()'s generator in row-major order; randn is Box-Muller from two draws, lib/random's formula [interpreter yes, C yes, direct wasm yes]

### Reductions, products and broadcasting

Arithmetic and comparisons between two arrays (not complex ones) broadcast as NumPy's do: the last axes line up, and each pair of lengths is equal or one of them is 1, so `A - mean(A, axis=0)` subtracts the column means and `A * transpose([w])` scales rows. `sum`, `mean`, `var`, `std` along an axis and `trace` add in NumPy's order (pairwise, in blocks of 8192), so they agree with NumPy to the bit; so does `sum(a)` of a whole array. `keepdims=true` keeps the reduced axis with length 1 (NumPy's keepdims).

- `sum(a, ordered=true)`: the sum of every element left to right from 0.0 (row-major), the order of a Fortran or C loop, in every engine; for a port that must match a reference's last digits. `ordered=false` is `sum(a)`. `dot` and `prod` always go left to right and take `ordered=true` to say so [interpreter yes, C yes, direct wasm yes]
- `sum(a, axis=axis) · sum(a, axis=axis, keepdims=keepdims)`: along one axis, in NumPy's order [interpreter yes, C yes, direct wasm yes]
- `mean(a) · mean(a, axis=axis) · mean(a, axis=axis, keepdims=keepdims)`: NumPy's pairwise sum divided by the count [interpreter yes, C yes, direct wasm yes]
- `var(a) · std(a) · var(a, axis=axis) · std(a, axis=axis)`: population variance (divides by n) and its square root, NumPy's order [interpreter yes, C yes, direct wasm yes]
- `prod(a) · prod(a, axis=axis)`: left to right [interpreter yes, C yes, direct wasm yes]
- `amax(a) · amin(a) · amax(a, axis=axis) · amin(a, axis=axis)`: NumPy's max and min: a NaN wins, -0.0 is below 0.0 (max(a) of one array skips NaN, as C's fmax) [interpreter yes, C yes, direct wasm yes]
- `argmax(a) · argmin(a) · argmax(a, axis=axis) · argmin(a, axis=axis)`: the index of the first largest (smallest) entry, flat or along an axis; the first NaN if there is one [interpreter yes, C yes, direct wasm yes]
- `any(a, axis=axis) · all(a, axis=axis)`: masks of 1.0 and 0.0 along an axis [interpreter yes, C yes, direct wasm yes]
- `count_nonzero(a)`: how many entries are not 0 [interpreter yes, C yes, direct wasm yes]
- `outer(u, v)`: u[i] * v[j], each argument flattened [interpreter yes, C yes, direct wasm yes]
- `dot(u, v) · dot(u, v, ordered=true)`: two vectors of one length, summed left to right (A @ B is the matrix product); `ordered=true` states that order and changes nothing [interpreter yes, C yes, direct wasm yes]
- `norm(a)`: sqrt(dot(a, a)) over every entry (the Frobenius norm of a matrix) [interpreter yes, C yes, direct wasm yes]
- `trace(A) · trace(A, offset=k)`: the sum of diagonal k, NumPy's order [interpreter yes, C yes, direct wasm yes]
- `triu(A) · triu(A, k=k) · tril(A) · tril(A, k=k)`: the entries on and above (below) diagonal k; the rest 0 [interpreter yes, C yes, direct wasm yes]
- `einsum("ij,jk->ik", A, B)`: a contraction of one to three arrays written in NumPy's subscripts (a string literal; an implicit output is the letters written once); a number when the output has no index. Summed left to right over the indices in their order of first appearance, so it agrees with NumPy to rounding [interpreter yes, C yes, direct wasm yes]

### Gathers, masks, sorting and scans

A mask is an array of 1.0 and 0.0, as a comparison gives, so an array cannot say whether it is a mask or a list of indices: selecting by a mask is `extract`, and by indices `take`.

- `take(A, idx, axis=axis)`: NumPy's take along an axis: the result has A's axes before it, idx's shape, then A's axes after it. take(A, idx) with two arguments is flat [interpreter yes, C yes, direct wasm yes]
- `extract(mask, a)`: the entries of a where mask is nonzero, as a vector (NumPy's a[mask]) [interpreter yes, C yes, direct wasm yes]
- `putmask(a, mask, v)`: a[mask] = v in place: v a number, or an array of a's shape whose entries go where mask is nonzero [interpreter yes, C yes, direct wasm yes]
- `flatnonzero(a)`: the flat indices of the nonzero entries [interpreter yes, C yes, direct wasm yes]
- `argsort(a)`: indices that sort each row (the last axis), stable, NaN last, -0.0 equal to 0.0: NumPy's argsort(kind="stable") [interpreter yes, C yes, direct wasm yes]
- `cumsum(a) · cumsum(a, axis=axis) · cumprod(a) · cumprod(a, axis=axis)`: running sums and products, left to right; without an axis over every entry, flat [interpreter yes, C yes, direct wasm yes]
- `diff(a) · diff(a, n=n) · diff(a, n=n, axis=axis)`: a[i + 1] - a[i] along the last axis (or axis), n times [interpreter yes, C yes, direct wasm yes]
- `clip(x, lo, hi)`: a number or an array held to [lo, hi]; a NaN entry stays NaN and a NaN bound is no bound, as NumPy's clip [interpreter yes, C yes, direct wasm yes]

### Strings

Python's semantics, by code point; `a.f(b)` is `f(a, b)`. See Strings above.

- `len(s)`: code points (a list: its strings) [interpreter yes, C yes, direct wasm yes]
- `str(x)`: a number as Python's repr, a boolean as true/false [interpreter yes, C yes, direct wasm yes]
- `repr(x)`: a string in quotes, as Python writes it [interpreter yes, C yes, direct wasm yes]
- `float(s)`: CPython's float(str), correctly rounded [interpreter yes, C yes, direct wasm yes]
- `int(s[, base])`: CPython's int(str, base), to the nearest double [interpreter yes, C yes, direct wasm yes]
- `chr(n)`: the character with code point n [interpreter yes, C yes, direct wasm yes]
- `ord(c)`: the code point of one character [interpreter yes, C yes, direct wasm yes]
- `upper(s)`: full Unicode mapping (ß is SS) [interpreter yes, C yes, direct wasm yes]
- `lower(s)`: full mapping, final sigma [interpreter yes, C yes, direct wasm yes]
- `strip(s[, chars])`: whitespace (or chars) off both ends; lstrip, rstrip one end [interpreter yes, C yes, direct wasm yes]
- `lstrip(s[, chars])` [interpreter yes, C yes, direct wasm yes]
- `rstrip(s[, chars])` [interpreter yes, C yes, direct wasm yes]
- `split(s[, sep[, maxsplit]])`: a list of strings; no sep: runs of whitespace [interpreter yes, C yes, direct wasm yes]
- `rsplit(s[, sep[, maxsplit]])`: from the right [interpreter yes, C yes, direct wasm yes]
- `splitlines(s)`: Python's line boundaries [interpreter yes, C yes, direct wasm yes]
- `join(sep, words)`: Python's sep.join(words) [interpreter yes, C yes, direct wasm yes]
- `find(s, sub[, start[, end]])`: first position, or -1; rfind the last [interpreter yes, C yes, direct wasm yes]
- `rfind(s, sub[, start[, end]])` [interpreter yes, C yes, direct wasm yes]
- `index(s, sub)`: find, but a miss stops; also a list's index [interpreter yes, C yes, direct wasm yes]
- `rindex(s, sub)` [interpreter yes, C yes, direct wasm yes]
- `count(s, sub)`: non-overlapping occurrences [interpreter yes, C yes, direct wasm yes]
- `replace(s, old, new[, count])` [interpreter yes, C yes, direct wasm yes]
- `startswith(s, p)` [interpreter yes, C yes, direct wasm yes]
- `endswith(s, p)` [interpreter yes, C yes, direct wasm yes]
- `contains(s, sub)`: the `in` operator: sub in s [interpreter yes, C yes, direct wasm yes]
- `ljust(s, width[, fill])`: rjust and center likewise [interpreter yes, C yes, direct wasm yes]
- `rjust(s, width[, fill])` [interpreter yes, C yes, direct wasm yes]
- `center(s, width[, fill])` [interpreter yes, C yes, direct wasm yes]
- `zfill(s, width)`: zeros after a sign [interpreter yes, C yes, direct wasm yes]
- `isdigit(s)`: also isspace isalpha isalnum isdecimal isnumeric isupper islower isascii isprintable [interpreter yes, C yes, direct wasm yes]
- `isalpha(s) · isalnum(s) · isdecimal(s) · isnumeric(s)`: Python's character classes, by CPython 3.11's Unicode tables: true when s is not empty and every character is in the class [interpreter yes, C yes, direct wasm yes]
- `isspace(s) · isupper(s) · islower(s) · isascii(s) · isprintable(s)`: Python's, character for character (isascii and isprintable are true for an empty string) [interpreter yes, C yes, direct wasm yes]
- `unique(s) · unique(words) · unique(xs)`: each character of a string, or each string of a list, once, in code point order: Python's sorted(set(s)); of numbers, the distinct values ascending, as np.unique [interpreter yes, C yes, direct wasm yes]
- `sorted(words) · sorted(xs)`: a new list, by code point; of numbers, a new array, ascending and stable, NaN last [interpreter yes, C yes, direct wasm yes]
- `reversed(x)`: x[::-1]: an array, a string or a list of strings, back to front [interpreter yes, C yes, direct wasm yes]
- `pop(xs) · pop(xs, i)`: a statement: `y = pop(xs, i)` takes element i (the last when left out) out of the list xs, as Python's [interpreter yes, C yes, direct wasm yes]
- `read_text(path) · read_lines(path) · write_text(path, s)`: a UTF-8 file; read_lines(path) its lines; write_text(path, s) [interpreter yes, C yes, direct wasm yes]

### Databases (SQLite)

SQLite (third_party/sqlite, public domain) is built into viv and compiled beside a C program that uses these (`viv export --with-sqlite`); the browser does not have them yet. A database is a number, its handle. Parameters are `?` in the SQL, bound from a list of strings (as text) or an array of numbers (NaN as NULL), never pasted into it. `db_exec` runs inside a savepoint: SQL that fails changes nothing. A query is one statement giving one column. NULL reads as NaN in `db_col` and is an error in `db_col_int` and `db_col_str` (select COALESCE(x, k)). SQLite's errors are VivError; a misuse (a bad handle, a wrong count of parameters, a NULL in an Int column) is ValueError.

- `db_open(path) · db_close(db)`: open or create the database file; its handle. db_close closes it [interpreter yes, C yes, direct wasm yes]
- `db_exec(db, sql) · db_exec(db, sql, params)`: run SQL (several statements, or one with `?` parameters) inside a savepoint; the rows it changed [interpreter yes, C yes, direct wasm yes]
- `db_col(db, sql, params) · db_col_int(db, sql, params) · db_col_str(db, sql, params)`: the query's one column as an array of numbers (NULL is NaN), an [Int], or a list of strings; params may be left out [interpreter yes, C yes, direct wasm yes]
- `db_stage(db, name, column) · db_insert(db, table)`: set a column ([Num], [Int] or [Str]) aside, then write every staged column as rows in one transaction (NaN as NULL); the rows written [interpreter yes, C yes, direct wasm yes]

### Maps

A map holds numbers or strings under number keys or string keys (one kind of each per map), in the order the keys were first stored. It is one handle: `b = m`, a def's parameter and a record field all change the same map. See Maps above.

- `map()`: a new empty map [interpreter yes, C yes, direct wasm yes]
- `m[k] = v · m[k]`: store v under k (a key already there keeps its place), or read it; a key that is not there stops the run [interpreter yes, C yes, direct wasm yes]
- `m["key"] = v`: a string key: printed in its repr, `{'the': 2}`; keys(m) is then a list of strings [interpreter yes, C yes, direct wasm yes]
- `has(m, k)`: true when k is a key of m [interpreter yes, C yes, direct wasm yes]
- `get(m, k, d)`: m[k] when k is there, else d [interpreter yes, C yes, direct wasm yes]
- `remove(m, k) · remove(m, k, d)`: take k out and give its value; without d a key that is not there stops the run [interpreter yes, C yes, direct wasm yes]
- `keys(m) · values(m)`: a new vector of the keys, or the values, in insertion order: loop over a map with `ks = keys(m)` then `for i to len(ks)` [interpreter yes, C yes, direct wasm yes]
- `copy(m)`: a new map with the same keys, values and order [interpreter yes, C yes, direct wasm yes]
- `clear(m)`: take every key out [interpreter yes, C yes, direct wasm yes]
- `len(m)`: the number of keys [interpreter yes, C yes, direct wasm yes]

### Group-by on Int keys

Written in Viv (lang/viv/std/groupby.viv) and added to a program that calls them, so every engine runs the same code. The keys stay Ints (no float conversion), hashed into an open-addressing table whose slot holds the key beside its group id. A group-by is factorize, then a dense loop: `s[ids[i]] += v[i]`.

- `uniq, ids = factorize(keys)`: the distinct Ints of an [Int] column in first-appearance order, and each row's position in them (`uniq[ids[i]] == keys[i]`), as pandas.factorize(sort=False); any Int is a key [interpreter yes, C yes, direct wasm yes]
- `ka, kb, ids = factorize2(a, b)`: factorize on two key columns at once: the distinct pairs in first-appearance order, as two columns, and each row's pair id; three columns: factorize2 its ids with the third [interpreter yes, C yes, direct wasm yes]
- `match_rows(build, probe)`: a hash join's probe: for each Int of probe, the first row of build with that key, or -1 (pandas.Index.get_indexer) [interpreter yes, C yes, direct wasm yes]
- `lo, hi = key_range(keys)`: the smallest and largest Int of a column, in one pass (0, 0 for no keys) [interpreter yes, C yes, direct wasm yes]
- `uniq, s, c = group_sum(keys, v)`: per-group sum and row count of an [Int] column in one pass, the groups in first-appearance order as factorize's; keys spanning few values (found inside the call) are slotted by key - lo, others hashed with the sums beside the table; group_sum_num sums a number column [interpreter yes, C yes, direct wasm yes]
- `ka, kb, s, c = group_sum2(a, b, v)`: group_sum keyed by two [Int] columns: one pass over pair slots when both span few values, else factorize2 and a loop [interpreter yes, C yes, direct wasm yes]

### Linear algebra

In debug mode each solve checks its own result (a residual, a pivot check or a trace invariant) and stops with an error instead of returning a wrong answer.

- `eig_sym(A, vals, vecs)`: symmetric spectrum by cyclic Jacobi; values ascending, vectors in columns [interpreter yes, C yes, direct wasm yes]
- `spectrum_sym(A, vals)`: the same eigenvalues eig_sym computes; in debug mode it runs every available method and asserts they agree [interpreter yes, C yes, direct wasm yes]
- `lin_solve(A, b, x)`: LU with partial pivoting; checks the residual in debug mode [interpreter yes, C yes, direct wasm yes]
- `lin_solve_multi(A, B, X)`: factor once, many right-hand sides [interpreter yes, C yes, direct wasm yes]
- `cholesky(A, L)`: factors a symmetric positive-definite matrix; checks every pivot [interpreter yes, C yes, direct wasm yes]
- `x = lin_solve(A, b) · X = lin_solve(A, B) · inv(A)`: the same LU as lin_solve and lin_solve_multi, into a new array: the same bits [interpreter yes, C yes, direct wasm yes]
- `L = cholesky(A)`: the same factorization, into a new array [interpreter yes, C yes, direct wasm yes]
- `vals, vecs = eig_sym(A) · w, V = eig_herm(A)`: the eigensolvers, unpacked into new arrays (written as the call that fills them) [interpreter yes, C yes, direct wasm yes]
- `U, s, Vt = svd(A) · svd(A, U, s, Vt)`: A = U diag(s) Vt with k = min(m, n) singular values, descending: one-sided Jacobi, lib/svd.viv's algorithm. The signs of the vectors can differ from LAPACK's [interpreter yes, C yes, direct wasm yes]
- `eig_general(A, re, im)`: balancing + Hessenberg + Francis double-shift QR; complex pairs split across re/im [interpreter yes, C yes, direct wasm yes]
- `eig_near(A, sigma, k, re, im)`: k eigenvalues near a real shift by shift-invert Arnoldi — one LU instead of a full QR [interpreter yes, C yes, direct wasm yes]

### Higher-order and contracts

A contract is named for what it computes, not for the method. In debug mode it runs every method that applies and asserts they agree; `viv run --release` runs the auto-selected one. The methods are ordinary Viv code, in lang/viv/std/contracts.viv.

- `simpson(f, a, b, n)`: composite Simpson [interpreter yes, C yes, direct wasm yes]
- `trapezoid(f, a, b, n)` [interpreter yes, C yes, direct wasm yes]
- `deriv(f, x0)`: the exact derivative, computed by forward-mode autodiff when the program is compiled [interpreter yes, C yes, direct wasm yes]
- `apply(f, args…)`: call a function held in a variable or an array [interpreter yes, C yes, direct wasm yes]
- `release_mode()`: true under --release, else false, so a program can choose what to run in each mode [interpreter yes, C yes, direct wasm yes]
- `clock()`: wall time in seconds from a monotonic clock; only the difference of two readings means anything, and only a program with a `#! timing` line may call it (every engine refuses it elsewhere; three-engines lists a timing program and does not compare it) [interpreter no, C yes, direct wasm yes]

### Derivatives and records

grad and value_and_grad are reverse mode over a whole def (lang/viv/src/grad.rs): the compiler writes the gradient def before any engine runs, so the three engines run the same generated code.

- `grad(f, x, args…)`: the gradient of the number def f returns with respect to its first argument x (a number, an array, or a record of arrays), with the shape of x; the other arguments are constants [interpreter yes, C yes, direct wasm yes]
- `l, g = value_and_grad(f, x, args…)`: the value and the gradient in one forward and one backward pass, unpacked on a line of its own [interpreter yes, C yes, direct wasm yes]
- `for_leaves(f, r1, r2, …, x)`: call def f on each array of one or more records of the same type at once (the records' arrays first, then x), in field order; f may write the arrays in place: an optimizer step over a model [interpreter yes, C yes, direct wasm yes]

### Complex scalars

The `i` suffix works only on a number, so a loop variable named `i` is still a variable. For arrays of them, see the next group.

- `complex(re, im)` [interpreter yes, C yes, direct wasm yes]
- `re(z)`: real part [interpreter yes, C yes, direct wasm yes]
- `im(z)`: imaginary part [interpreter yes, C yes, direct wasm yes]
- `conj(z)` [interpreter yes, C yes, direct wasm yes]

### Complex arrays and the Fourier transform

A complex array is packed — interleaved (re, im) doubles, which is exactly how C99 lays out `double complex`, so both engines walk the same memory and agree to the byte. `fft` is radix-2 at a power of two, mixed radix (with dedicated radix-3 and radix-5 butterflies) when the primes of the length are all ≤ 7, and Bluestein everywhere else, so no length is a special case and 2·3·5-smooth ones are fast. The twiddle factors and factorization for each length are computed once and reused. Whole-array arithmetic works on complex arrays: `a + b`, `a * b`, a real array or a real or complex scalar on either side of `*`, `/` by a scalar, and `@`. Convention: `fft` is unscaled with exp(−2πikn/N); `ifft` carries the 1/N — the same as numpy and FFTW.

- `czeros(d1 … d4)`: a complex array of zeros, rank 1 to 4 [interpreter yes, C yes, direct wasm yes]
- `fft(a)`: forward transform, in place, any length [interpreter yes, C yes, direct wasm yes]
- `ifft(a)`: inverse transform, in place; ifft(fft(x)) returns x to rounding error [interpreter yes, C yes, direct wasm yes]
- `re(a)`: real parts, as a real array of the same shape [interpreter yes, C yes, direct wasm yes]
- `im(a)`: imaginary parts, as a real array [interpreter yes, C yes, direct wasm yes]
- `abs(a)`: magnitudes — what you plot as a spectrum [interpreter yes, C yes, direct wasm yes]
- `conj(a)`: the conjugate, as another complex array [interpreter yes, C yes, direct wasm yes]
- `complex(re, im)`: a complex array from two real ones of the same shape [interpreter yes, C yes, direct wasm yes]
- `fft3(a)`: forward transform of a rank-3 complex array along all three axes, in place; all-zero lines are skipped [interpreter yes, C yes, direct wasm yes]
- `ifft3(a)`: its inverse, carrying 1/N [interpreter yes, C yes, direct wasm yes]
- `eig_herm(A, w, V)`: the Hermitian eigenproblem: ascending eigenvalues into w, eigenvectors into the COLUMNS of V. Householder to real tridiagonal, then implicit QL [interpreter yes, C yes, direct wasm yes]
- `transpose(A)`: plain transpose, NO conjugation: the adjoint is conj(transpose(A)), spelled out where it is meant [interpreter yes, C yes, direct wasm yes]
- `slab(A, k) · setslab · take · put`: slab, setslab, take and put work on complex arrays too; source and target must both be real or both complex [interpreter yes, C yes, direct wasm yes]

### Randomness

Park–Miller MINSTD in pure f64, so the stream is bit-identical between the interpreter and compiled C — a seeded Monte Carlo gives the same answer either way. `#! use lib/random` adds randn, randn_matrix, randint, shuffle, randperm, choice, arange and pi, written in Viv over rand(), each documented as the formula it computes so Python can replay it (tools/random_replay.py); `#! use lib/stats` adds mean and std.

- `srand(seed)`: seed the stream [interpreter yes, C yes, direct wasm yes]
- `rand()`: uniform on [0, 1) [interpreter yes, C yes, direct wasm yes]

### Figures and geometry

Each of these prints a structured record on stdout, and the page draws the chart from it. Compiled C prints the same record, so a figure is part of the program's output, not a call into a plotting library.

- `plot(x, y, "label")`: a line series; consecutive calls sharing x overlay [interpreter yes, C yes, direct wasm yes]
- `scatter(x, y, "label")`: points [interpreter yes, C yes, direct wasm yes]
- `bars(x, y, "label")`: a bar chart: a histogram's counts at its bin centers. A plot() on the same x afterwards draws a curve over the bars [interpreter yes, C yes, direct wasm yes]
- `plot2d(Z, "label")`: a rank-2 array as a viridis heatmap with a colorbar [interpreter yes, C yes, direct wasm yes]
- `plot3d(x, y, z, "label")`: a curve in 3-space, rotatable [interpreter yes, C yes, direct wasm yes]
- `quiver(x, y, u, v, "label")`: a vector field [interpreter yes, C yes, direct wasm yes]
- `scene("title")`: open a 3D scene; wire/verts output belongs to it [interpreter yes, C yes, direct wasm yes]
- `wire(x, y, z, "group")`: a 3D polyline, optionally colormapped by a 4th array [interpreter yes, C yes, direct wasm yes]
- `verts(x, y, z, "group")`: point markers — joints, atoms, nodes [interpreter yes, C yes, direct wasm yes]
- `figure(title, xlabel, ylabel[, xlo, xhi, ylo, yhi[, style]])`: REQUIRED before the first plot. The program sets the axis labels and limits; the page adds none of its own. On a plot2d the window is the picture's physical extent and gives it real axes; style is "contour" (overlay level lines) or "pixels" (no smoothing: the cells are categories or counts) [interpreter yes, C yes, direct wasm yes]
- `xtick("text", x)`: inside a figure: a tick labeled with the text at x (a high-symmetry point, a category). A chart with any xtick shows only those x ticks. [interpreter yes, C yes, direct wasm yes]
- `xscale("log") · yscale("log")`: a logarithmic x or y axis for the open figure: plot the values, not their logs. Ticks at powers of ten; a value of zero or less is a gap, and the chart says how many it left out. "linear" is the default [interpreter yes, C yes, direct wasm yes]
- `aspect("equal")`: one unit of x as long on the screen as one unit of y, for an orbit or a shape in a plane: the plot takes the height that makes a circle round, and widens a range only if that height would be too tall or too flat. "auto" is the default [interpreter yes, C yes, direct wasm yes]
- `fill_between("label", x, lo, hi)`: a band shaded from lo to hi over x: a confidence band, a gap, a tolerance. Grouped like a plot; takes its label's color, so a band labeled like its curve is drawn in the curve's color [interpreter yes, C yes, direct wasm yes]
- `rangebars("label", x, lo, hi)`: a capped bar from lo to hi at each x: error bars, spreads, ranges. It takes the ends, not a center and an error, so mean ± err is rangebars("m", x, m - err, m + err) and an asymmetric interval is no harder; the center is a plot or scatter of its own [interpreter yes, C yes, direct wasm yes]
- `grid("title", rows, cols) · share("x")`: the next rows * cols figures become the panels of one picture, row by row, each drawn as it would be alone. share("x") gives each column one x range and labels it under the lowest panel; share("y") does the same along each row; share("xy") both [interpreter yes, C yes, direct wasm yes]
- `yright("series label", "axis label")`: the series with that label read against a second, linear y axis on the right, with its own ticks and label: two quantities in different units over one x. The legend marks the entry "(right axis)" [interpreter yes, C yes, direct wasm yes]
- `style("label", "dash s2 ring thick")`: how the open figure draws the series with that label: a line (solid, dash, dot, noline), a palette slot (s1 to s6), a marker (circle, ring, square, cross, nomarker) and a width (thin, thick), at most one of each. Any other word is a compile error. With no style, series take the palette in order, solid [interpreter yes, C yes, direct wasm yes]
- `ytick("text", y)`: inside a figure: a y tick labeled with the text; a chart with any shows only those y ticks [interpreter yes, C yes, direct wasm yes]
- `origin("lower") · colormap("diverging")`: heatmap settings for the open figure: row 0 at the bottom instead of the top; a blue-white-red map centered on zero for a signed field [interpreter yes, C yes, direct wasm yes]
- `clim(lo, hi) · levels(values)`: heatmap color limits (beyond them, the end color), and contour lines at exactly these values [interpreter yes, C yes, direct wasm yes]
- `box("id", "style", x0, y0, x1, y1)`: in a figure with the style "diagram": a box by two corners. Style words, at most one of each: a color (s1 to s6, ink, muted), a fill (plain, tint, clear), a line (solid, dash, dot, noline), a width (thin, thick), corners (round, square) [interpreter yes, C yes, direct wasm yes]
- `arrow("id", "style", xs, ys)`: in a diagram: a line through the points with a head at the last. Words: a color, a line (solid, dash, dot), a width, head or nohead, straight or curved [interpreter yes, C yes, direct wasm yes]
- `note("id", "text", "style", x, y)`: in a diagram: text at a point. Words: a color, left, center or right, top, middle or bottom (which part of the text is at the point), small, normal or large (11, 13, 15 units), bold, mono [interpreter yes, C yes, direct wasm yes]
- `checkpoint()`: save the figures so far as one frame; a sequence of frames plays as an animation [interpreter yes, C yes, direct wasm yes]

### Animation

An animation is a Viv program. Objects are created hidden and named by a string id; transitions reveal, morph and move them. Each call prints one record, and the page plays the records with a scrubber.

- `stage("title", xlo, xhi, ylo, yhi)`: opens an animation; the window is its coordinate system [interpreter yes, C yes, direct wasm yes]
- `axes("x label", "y label")`: ticks and a grid on the stage [interpreter yes, C yes, direct wasm yes]
- `curve("id", x, y)`: a named polyline, created hidden. The same id again replaces its data [interpreter yes, C yes, direct wasm yes]
- `dots("id", x, y) · dots("id", x, y, r)`: named point markers; with a radius, discs of that size in data units (atoms) [interpreter yes, C yes, direct wasm yes]
- `circle("id", x, y, r)`: a disc [interpreter yes, C yes, direct wasm yes]
- `rect("id", x, y, w, h) · rect("id", x, y, w, h, angle)`: a rectangle centered at (x, y), turned by angle radians [interpreter yes, C yes, direct wasm yes]
- `segment("id", x1, y1, x2, y2)`: a line between two points: a pole, a rod, a spring [interpreter yes, C yes, direct wasm yes]
- `poly("id", xs, ys)`: a filled polygon [interpreter yes, C yes, direct wasm yes]
- `frame(seconds)`: a KEYFRAME: advance the clock; an object set again later moves there over the time between. In a keyframed animation objects are visible from creation [interpreter yes, C yes, direct wasm yes]
- `label("id", "text", x, y)`: text at a point; $...$ is TeX, with \\ for a backslash [interpreter yes, C yes, direct wasm yes]
- `hline("id", y) · vline("id", x)`: a dashed reference line across the stage [interpreter yes, C yes, direct wasm yes]
- `draw("id", seconds)`: a curve is drawn on, start to end; anything else fades in [interpreter yes, C yes, direct wasm yes]
- `fade_in("id", seconds) · fade_out("id", seconds)`: reveal or hide any object [interpreter yes, C yes, direct wasm yes]
- `morph("id", x, y, seconds)`: move a curve's or dots' points smoothly to new data over the given seconds; the point count may change [interpreter yes, C yes, direct wasm yes]
- `move_to("id", x, y, seconds)`: slide a label [interpreter yes, C yes, direct wasm yes]
- `say("caption") · say("caption", value)`: the caption until the next say; a value is appended, which is how a caption reports a computed number [interpreter yes, C yes, direct wasm yes]
- `wait(seconds)`: hold [interpreter yes, C yes, direct wasm yes]
- `together() … apart()`: the transitions between them start at the same time [interpreter yes, C yes, direct wasm yes]

### Output and testing

A failed assert prints what failed (for assert_near, the two values, their difference and the tolerance), keeps running, and makes the process exit nonzero at the end, so a build can stop on it while you still see every failure in one run.

- `print(x)`: the bare value [interpreter yes, C yes, direct wasm yes]
- `print("text", x, ..., sep=" ")`: text and values on one line, a space (or sep) between them, as Python's print [interpreter yes, C yes, direct wasm yes]
- `fmt(x)`: a print argument: the shortest digits that read back as x, as Python's repr [interpreter yes, C yes, direct wasm yes]
- `fmt(x, 17)`: a print argument: 1 to 17 significant digits, as %.17g [interpreter yes, C yes, direct wasm yes]
- `fmt(x, ".3f")`: a print argument: a Python format, [[fill]align][sign][0][width][,][.precision] then e, f, g, d, % (or E, F, G) [interpreter yes, C yes, direct wasm yes]
- `narrate "text"`: prints a progress line while the program runs, separate from data output [interpreter yes, C yes, direct wasm yes]
- `text(codes)`: writes an array of Unicode code points as one line of text (the counterpart of chars) [interpreter yes, C yes, direct wasm yes]
- `assert_near(a, b, tol)` [interpreter yes, C yes, direct wasm yes]
- `assert_true(cond)` [interpreter yes, C yes, direct wasm yes]
- `assert_gt(a, b)` [interpreter yes, C yes, direct wasm yes]
- `assert_lt(a, b)` [interpreter yes, C yes, direct wasm yes]
- `assert_allclose(a, b, tol)`: elementwise, over arrays [interpreter yes, C yes, direct wasm yes]

## Examples, by topic

### Numbers & arithmetic

Numbers are 64-bit floats, Int is a 64-bit integer, and Bool holds truth values.

**Integer literals are floats.** 10 is read as 10.0, so division never truncates.

```viv
n = 10
show n / 4.0
```

Prints:
```text
2.5
```

**print vs show.** print gives the bare value. show labels a variable with its name, and prints any other expression bare, as print would.

```viv
x = 2.0
print(x)
show x
show x + 1.0
```

Prints:
```text
2
x: 2
3
```

**Scientific notation.**

```viv
tiny = 2.5e-3
show tiny * 4.0e3
```

Prints:
```text
10
```

**Compound assignment.** x += 1 means x = x + 1, and the same holds for -=, *= and /=.

```viv
x = 9.0
x += 1.0
x -= 2.0
x *= 3.0
x /= 4.0
show x
```

Prints:
```text
x: 6
```

**The math builtins.** 29 functions. pow, atan2, hypot, copysign and mod take two arguments, fma takes three, min and max take two or more, and fsum gives the correctly rounded sum of an array. mod is C's fmod.

```viv
show pow(2.0, 10.0)
show atan2(1.0, 1.0)
show mod(7.5, 2.0)
show round(2.5)
```

Prints:
```text
1024
0.785398
1.5
3
```

**erf and gamma.** The engines agree to double precision; gamma(5) is 4! and gamma(1/2) is √π.

```viv
show erf(1.0)
show gamma(5.0)
show gamma(0.5) * gamma(0.5)
```

Prints:
```text
0.842701
24
3.14159
```

**Precedence.** Multiplication binds tighter than addition; parentheses override.

```viv
show 2.0 + 3.0 * 4.0
show (2.0 + 3.0) * 4.0
```

Prints:
```text
14
20
```

**Division by zero and sqrt of a negative.** Following IEEE float rules, these give inf and nan instead of stopping the program.

```viv
show 1.0 / 0.0
show sqrt(0.0 - 1.0)
show 1.0e300 * 1.0e300
```

Prints:
```text
inf
nan
inf
```

**mod with negatives.** mod is C's fmod: the result takes the sign of the first argument.

```viv
show mod(7.0, 3.0)
show mod(0.0 - 7.0, 3.0)
```

Prints:
```text
1
-1
```

**round halves go away from zero.** round follows C: halves round away from zero, so 2.5 gives 3 and -2.5 gives -3, not the nearest even number.

```viv
show round(2.5)
show round(3.5)
show round(0.0 - 2.5)
```

Prints:
```text
3
4
-3
```

**pow(0, 0) is 1.** pow(0, 0) returns 1, as C does.

```viv
show pow(0.0, 0.0)
```

Prints:
```text
1
```

**atan2 knows the quadrant.** atan2(y, x) uses the signs of both arguments to return the angle in the correct quadrant.

```viv
show atan2(1.0, 1.0)
show atan2(1.0, 0.0 - 1.0)
show atan2(0.0 - 1.0, 0.0 - 1.0)
```

Prints:
```text
0.785398
2.35619
-2.35619
```

**abs, min, max.**

```viv
show abs(0.0 - 3.5)
show min(2.0, 0.0 - 5.0)
show max(2.0, 0.0 - 5.0)
```

Prints:
```text
3.5
-5
2
```

**The hyperbolics and inverse trig.**

```viv
show sinh(1.0)
show cosh(1.0) * cosh(1.0) - sinh(1.0) * sinh(1.0)
show asin(1.0)
show acos(0.0)
```

Prints:
```text
1.1752
1
1.5708
1.5708
```

**% and // round toward minus infinity.** `a % b` takes the sign of b and `a // b` is the floor of the quotient, as in Python, so a == (a // b) * b + a % b. mod(a, b) is C's fmod and keeps the sign of a.

```viv
show -7.0 % 3.0
show 7.0 % -3.0
show -7.0 // 2.0
show mod(-7.0, 3.0)
```

Prints:
```text
2
-2
-4
-1
```

**^ is pow.** x ^ y is pow(x, y). It binds tighter than a minus on its left, so -3 ^ 2 is -9, and it groups from the right.

```viv
show 2.0 ^ 10.0
show -3.0 ^ 2.0
show 2.0 ^ 3.0 ^ 2.0
```

Prints:
```text
1024
-9
512
```

### Booleans

true and false are their own type: what a comparison yields, what if and while take, and a value you can store, pass and return.

**A comparison is a Bool.** It prints as true or false, in every engine. It is not the number 1: there is no arithmetic on it.

```viv
show 2.0 > 1.0
show 2.0 == 3.0
```

Prints:
```text
true
false
```

**Store one in a variable.** A flag is a Bool: a comparison gives one, and `and`, `or` and `not` combine them.

```viv
found = false
n = 7.0
small = n < 10.0
show small
both = small and not found
show both
```

Prints:
```text
small: true
both: true
```

**while not done.** not flips a boolean, so while not done keeps looping until done becomes true.

```viv
done = false
count = 0.0
while not done {
    count = count + 1.0
    done = count >= 4.0
}
show count
```

Prints:
```text
count: 4
```

**Functions can return one.** Bool is an annotation like Num and Arr. A predicate is an ordinary function.

```viv
def on_boundary(x: Num, y: Num, tol: Num) -> Bool {
    return x < tol or x > 1.0 - tol or y < tol or y > 1.0 - tol
}
show on_boundary(0.0, 0.5, 1.0e-9)
show on_boundary(0.5, 0.5, 1.0e-9)
```

Prints:
```text
true
false
```

**Bools compare with == and !=.** Two Bools can be compared with == and !=, and combined with and, or and not; no other operator accepts them.

```viv
a = 1.0 < 2.0
b = 3.0 < 2.0
show a == b
show a != b
show not b
```

Prints:
```text
false
true
true
```

**assert_true takes one.** A test can check a condition directly.

```viv
done = 3.0 > 2.0
test "three exceeds two" {
    assert_true(done)
}
```

Prints:
```text
PASS three exceeds two
```

**and and or stop early.** The right side of `and` runs only when the left side is true, and the right side of `or` only when the left side is false, so a bounds test can guard the index beside it.

```viv
a = [1.0, 2.0]
j = 0.0
while j < 2.0 and a[j] > 0.0 {
    j = j + 1.0
}
show j
```

Prints:
```text
j: 2
```

**Comparing arrays gives a mask.** A comparison with an array on either side is elementwise and gives 1.0 where the comparison is true and 0.0 where it is false. all and any turn a mask into one Bool, and where picks from two values with it.

```viv
x = [-1.0, 2.0, 0.5]
show x > 0.0
show x * (x > 0.0)
show all(x > -2.0)
show any(x > 1.0)
show where(x > 0.0, x, 0.0)
```

Prints:
```text
[0, 1, 1]
[0, 2, 0.5]
true
true
[0, 2, 0.5]
```

### Control flow

if/else, while, the counted for loop, and the loop-control statements.

**if / else if / else.**

```viv
x = 2.0
if x < 1.0 {
    show 1.0
} else if x < 3.0 {
    show 2.0
} else {
    show 3.0
}
```

Prints:
```text
2
```

**The counted loop.** for i to n runs i from 0 to n−1.

```viv
total = 0.0
for i to 5 {
    total += i
}
show total
```

Prints:
```text
total: 10
```

**while and break.**

```viv
p = 1.0
while 1.0 < 2.0 {
    p = p * 2.0
    if p > 100.0 { break }
}
show p
```

Prints:
```text
p: 128
```

**continue.** Sum only the loop indices above 2.

```viv
s = 0.0
for i to 6 {
    if i < 2.5 { continue }
    s += i
}
show s
```

Prints:
```text
s: 12
```

**Conditions compose with and / or / not.**

```viv
x = 5.0
if x > 0.0 and not (x > 10.0) {
    show x
}
```

Prints:
```text
x: 5
```

**Nested loops.** break leaves the innermost loop only.

```viv
hits = 0.0
for i to 3 {
    for j to 5 {
        if j > 1.5 { break }
        hits += 1.0
    }
}
show hits
```

Prints:
```text
hits: 6
```

**Computed and fractional bounds.** The bound is any expression; i counts 0, 1, 2, … while i < bound, so a bound of 2.5 runs three times.

```viv
n = 2.0 * 3.0
c1 = 0.0
for i to n { c1 += 1.0 }
show c1
c2 = 0.0
for k to 2.5 { c2 += 1.0 }
show c2
```

Prints:
```text
c1: 6
c2: 3
```

**The loop variable outlives a top-level loop.** At top level the loop variable is a global; after the loop it holds its last value.

```viv
for i to 3 {
    x = i
}
show i
```

Prints:
```text
i: 2
```

**A while that never runs.**

```viv
ran = 0.0
while 1.0 > 2.0 {
    ran = 1.0
}
show ran
```

Prints:
```text
ran: 0
```

**Collatz length.** A while loop with if/else inside; mod(n, 2) tells even from odd.

```viv
def collatz_len(n) {
    steps = 0.0
    while n > 1.5 {
        if mod(n, 2.0) < 0.5 {
            n = n / 2.0
        } else {
            n = 3.0 * n + 1.0
        }
        steps += 1.0
    }
    return steps
}
show collatz_len(27.0)
```

Prints:
```text
111
```

**An if expression.** Only the branch chosen is computed.

```viv
x = -2.0
show "neg" if x < 0.0 else "pos"
print(x, "and", x * x, sep=" ")
```

Prints:
```text
'neg'
-2 and 4
```

### Functions

def, return, docstrings, recursion, and optional type annotations.

**def and return.**

```viv
def hyp(a, b) {
    return sqrt(a * a + b * b)
}
show hyp(3.0, 4.0)
```

Prints:
```text
5
```

**Docstrings.** A bare string as the first statement documents the function and compiles to nothing.

```viv
def double(x) {
    "twice the input"
    return 2.0 * x
}
show double(21.0)
```

Prints:
```text
42
```

**Recursion.** A double holds every integer up to 2^53 exactly, so a counter can stop the recursion with equality; a comparison such as `n < 1.5` also stops it when n is not a whole number.

```viv
def fact(n) {
    if n < 1.5 {
        return 1.0
    }
    return n * fact(n - 1.0)
}
show fact(6.0)
```

Prints:
```text
720
```

**Type annotations.** Optional; the compiler checks them before the program runs.

```viv
def scale(v: Arr, s: Num) -> Arr {
    return v * s
}
show scale([1.0, 2.0], 10.0)
```

Prints:
```text
[10, 20]
```

**A def with no return gives nothing.** A def that ends without a `return` is called as a statement. Its value is nothing: the compiler refuses `show f(x)` or `f(x) + 1.0` for it.

```viv
def report(x) {
    show x + 1.0
}
report(5.0)
```

Prints:
```text
6
```

**Early return.**

```viv
def clamp01(x) {
    if x < 0.0 { return 0.0 }
    if x > 1.0 { return 1.0 }
    return x
}
show clamp01(0.0 - 2.0)
show clamp01(0.4)
show clamp01(9.0)
```

Prints:
```text
0
0.4
1
```

**Call before definition.** Top-level defs are collected before anything runs, so source order does not matter.

```viv
show later(4.0)
def later(x) {
    return x * 10.0
}
```

Prints:
```text
40
```

**Lambdas pass as arguments.** A def can take a function value and call it.

```viv
def apply2(f, x) {
    return f(f(x))
}
show apply2(fn(t) -> t * t, 3.0)
```

Prints:
```text
81
```

**Parameters shadow globals.**

```viv
x = 100.0
def f(x) {
    return x + 1.0
}
show f(5.0)
show x
```

Prints:
```text
6
x: 100
```

**Several results.** return a, b gives two values, and the caller unpacks them on a line of its own; a, b = b, a swaps.

```viv
def divmod2(a, b) {
    q = floor(a / b)
    return q, a - q * b
}
q, r = divmod2(17.0, 5.0)
show q
show r
q, r = r, q
show q
```

Prints:
```text
q: 3
r: 2
q: 2
```

### Arrays

Literals, zeros, indexing, elementwise arithmetic, and passing arrays through functions.

**Elementwise math on an array.** Every 1-argument math builtin maps elementwise over an array. No loop.

```viv
x = linspace(0.0, 3.14159265, 5)
y = sin(x)
show y
```

Prints:
```text
y: [0, 0.707107, 1, 0.707107, 3.58979e-09]
```

**Broadcasting an array and a scalar.** +, -, *, / broadcast between arrays and scalars, in either order.

```viv
x = linspace(0.0, 2.0, 3)
show 2.0 * x + 1.0
show 1.0 - x
show x / 2.0
```

Prints:
```text
[1, 3, 5]
[1, 0, -1]
[0, 0.5, 1]
```

**Resetting an array in place with fill.** fill(a, v) sets every element. Rebinding (a = zeros(n)) makes a new array, and a caller that passed a in never sees the change; fill changes the array in place.

```viv
a = zeros(3)
a[0] = 7.0
fill(a, 0.0)
show a
fill(a, 1.5)
show a
```

Prints:
```text
a: [0, 0, 0]
a: [1.5, 1.5, 1.5]
```

**Literals and indexing.** Comma-separated subscripts, zero-based.

```viv
v = [10.0, 20.0, 30.0]
M = [[1.0, 2.0], [3.0, 4.0]]
show v[1]
show M[1, 0]
```

Prints:
```text
20
3
```

**zeros and indexed stores.** zeros allocates rank 1 through 4.

```viv
Z = zeros(2, 2)
Z[0, 1] = 7.0
show Z
```

Prints:
```text
Z: [[0, 7],
    [0, 0]]
```

**Elementwise arithmetic.** * between arrays is elementwise (the Hadamard product); array times scalar scales.

```viv
a = [1.0, 2.0]
b = [10.0, 20.0]
show a + b
show a * b
show a * 5.0
```

Prints:
```text
[11, 22]
[10, 40]
[5, 10]
```

**len and linspace.** len reads the leading dimension.

```viv
g = linspace(0.0, 1.0, 5)
show g
show len(g)
```

Prints:
```text
g: [0, 0.25, 0.5, 0.75, 1]
5
```

**Functions take and return arrays.** head_tail takes an array and returns a new one.

```viv
def head_tail(v) {
    w = zeros(2)
    w[0] = v[0]
    w[1] = v[len(v) - 1.0]
    return w
}
show head_tail([3.0, 9.0, 27.0])
```

Prints:
```text
[3, 27]
```

**Scalars broadcast.** A scalar combines elementwise with an array from either side, for every arithmetic operator.

```viv
v = [1.0, 2.0, 4.0]
show v + 10.0
show 1.0 + v
show 2.0 / v
show v / 2.0
```

Prints:
```text
[11, 12, 14]
[2, 3, 5]
[2, 1, 0.5]
[0.5, 1, 2]
```

**Matrices combine elementwise too.**

```viv
A = [[1.0, 2.0], [3.0, 4.0]]
B = [[10.0, 20.0], [30.0, 40.0]]
show A + B
show A * B
```

Prints:
```text
[[11, 22],
 [33, 44]]
[[10, 40],
 [90, 160]]
```

**Rank 3 and 4 exist.** Allocate with zeros; index with one subscript per axis. show prints the flat data for rank 3+.

```viv
z = zeros(2, 2, 2)
z[1, 0, 1] = 5.0
show z[1, 0, 1]
```

Prints:
```text
5
```

**linspace endpoints are inclusive.** Both endpoints are in the array, as in NumPy. With n = 1 it is the start alone, and with n = 0 it is empty.

```viv
show linspace(0.0, 10.0, 5)
show linspace(3.0, 9.0, 1)
```

Prints:
```text
[0, 2.5, 5, 7.5, 10]
[3]
```

**dot is the loop you would write.** dot(a, b) is a builtin. The def below is the loop it runs, left to right, and gives the same sum.

```viv
def mydot(a, b) {
    s = 0.0
    for i to len(a) {
        s += a[i] * b[i]
    }
    return s
}
u = [1.0, 2.0, 3.0]
v = [4.0, 5.0, 6.0]
show mydot(u, v)
show dot(u, v)
```

Prints:
```text
32
32
```

**Negative indices count from the end.** As in Python: a[-1] is the last element, a[:-1] all but it, A[-1, :] the last row.

```viv
a = [10.0, 11.0, 12.0]
show a[-1]
show a[:-1]
A = [[1.0, 2.0], [3.0, 4.0]]
show A[-1, :]
```

Prints:
```text
12
[10, 11]
[3, 4]
```

**A comprehension builds an array.** Its for takes any loop head, and its if keeps the elements it holds for.

```viv
xs = [3.0, 1.0, 4.0, 1.0, 5.0]
show [x * x for x in xs if x > 1.0]
show [i * 2.0 for i to 4]
```

Prints:
```text
[9, 16, 25]
[0, 2, 4, 6]
```

**in, and a key for max and sorted.** x in a is NumPy's: some element equals x. A key is a def, a builtin or fn(x) -> ....

```viv
a = [3.0, -7.0, 4.0]
show 4.0 in a
show max(a, key=abs)
show sorted(a, key=fn(v) -> v * v)
```

Prints:
```text
true
-7
[3, 4, -7]
```

**norm is the square root of dot.** norm(v) is a builtin too, the square root of dot(v, v).

```viv
v = [3.0, 4.0]
show norm(v)
show sqrt(dot(v, v))
```

Prints:
```text
5
5
```

### Complex numbers

A complex number type: the `i` suffix on a number makes it imaginary, arithmetic mixes complex and real values, and every engine (the interpreter, C and wasm) computes the complex functions with the same formulas, so they print the same digits.

**Imaginary literals and arithmetic.** `2.5i` is an imaginary literal (a bare identifier `i` is still a variable — only digits take the suffix). +, -, *, / work between complex and real values; reals promote.

```viv
z = 2.0 + 3.0i
w = complex(1.0, -1.0)
show z
show z * w
show z / w
show z + 1.0
```

Prints:
```text
z: 2+3i
5+1i
-0.5+2.5i
3+3i
```

**re, im, conj, abs.** abs of a complex number is its modulus; z times its conjugate is exactly |z| squared.

```viv
z = 3.0 + 4.0i
show re(z)
show im(z)
show conj(z)
show abs(z)
show re(z * conj(z))
```

Prints:
```text
3
4
3-4i
5
25
```

**Euler's identity, computed.** exp, log, sqrt, sin, cos, tan all work on complex arguments. e^{i pi} comes out as -1 with an imaginary part of 1.2e-16, which is sin of the double nearest pi. sqrt takes the half-angle formulas, so the square root of -1 is exactly i.

```viv
pi = 3.141592653589793
show exp(pi * 1.0i)
show abs(exp(1.0i))
show sqrt(-1.0 + 0.0i)
```

Prints:
```text
-1+1.22465e-16i
1
0+1i
```

**A rotation is a complex multiply.** Multiplying by e^{i theta} rotates the plane by theta. A function that takes and returns complex values needs no annotations; the compiler infers the types.

```viv
def rotate(z, theta) {
    return z * exp(complex(0.0, theta))
}
q = rotate(1.0 + 0.0i, 1.5707963267948966)
show im(q)
test "quarter turn sends 1 to i" {
    assert_near(im(q), 1.0, 1.0e-12)
    assert_near(re(q), 0.0, 1.0e-12)
}
```

Prints:
```text
1
PASS quarter turn sends 1 to i
```

### Complex arrays & the FFT

Packed complex arrays, and a Fourier transform that works at every length and gives identical output in every engine.

**A complex array.** czeros allocates; indexing gives a complex number, and storing a real promotes it.

```viv
a = czeros(3)
a[0] = 1.0 + 2.0i
a[1] = 5.0
show a
show len(a)
```

Prints:
```text
a: [1+2i, 5+0i, 0+0i]
3
```

**A literal with a complex in it.** One complex element makes the whole literal a packed complex array; the reals promote.

```viv
z = [1.0 + 2.0i, 3.0, 0.0 - 1.5i]
show z
```

Prints:
```text
z: [1+2i, 3+0i, 0-1.5i]
```

**Real, imaginary, magnitude.** re, im and abs over a complex array give a REAL array of the same shape. abs is what you plot as a spectrum.

```viv
z = [3.0 + 4.0i, 0.0 - 2.0i]
show re(z)
show im(z)
show abs(z)
```

Prints:
```text
[3, 0]
[4, -2]
[5, 2]
```

**A delta transforms to all ones.** A unit impulse at index 0 transforms to all ones, which also confirms that fft does not scale its output.

```viv
d = czeros(8)
d[0] = 1.0
fft(d)
show d
```

Prints:
```text
d: [1+0i, 1+0i, 1+0i, 1+0i, 1+0i, 1+0i, 1+0i, 1+0i]
```

**A pure tone shows up in two bins.** A cosine at frequency 3 on 16 points puts N/2 = 8 into bins 3 and 13. The test checks that each of the other 14 bins is below 1e-12.

```viv
PI = 3.141592653589793
n = 16
t = czeros(n)
for j to n {
    t[j] = cos(2.0 * PI * 3.0 * j / n)
}
fft(t)
m = abs(t)
show m[3]
show m[13]
test "every other bin is rounding noise" {
    for k to n {
        if k != 3.0 and k != 13.0 {
            assert_lt(m[k], 1.0e-12)
        }
    }
}
```

Prints:
```text
8
8
PASS every other bin is rounding noise
```

**ifft inverts fft.** ifft carries the 1/N, so the round trip returns what you started with — the same convention as numpy and FFTW.

```viv
x = [1.0 + 0.5i, 2.0, 0.0 - 1.0i, 4.0]
y = [1.0 + 0.5i, 2.0, 0.0 - 1.0i, 4.0]
fft(y)
ifft(y)
test "the round trip is the identity" {
    for j to 4 {
        assert_lt(abs(y[j] - x[j]), 1.0e-14)
    }
}
```

Prints:
```text
PASS the round trip is the identity
```

**Any length, not just powers of two.** Seven is prime. Lengths that are not a power of two go through Bluestein's algorithm, so they cost N log N too and nothing about the call changes.

```viv
PI = 3.141592653589793
n = 7
x = czeros(n)
for j to n {
    x[j] = cos(2.0 * PI * 2.0 * j / n)
}
fft(x)
m = abs(x)
show m[2]
test "a prime-length tone still lands in its bin" {
    assert_near(m[2], 3.5, 1.0e-12)
    assert_lt(m[1], 1.0e-12)
}
```

Prints:
```text
3.5
PASS a prime-length tone still lands in its bin
```

**A real array is not a complex one.** fft refuses a real array and says what to do instead. The two are different types and never mix silently.

```viv
x = zeros(8)
fft(x)
```

Stops with an error containing: `` `fft` needs a complex array ``

### Linear algebra

The matrix product, transpose, and the symmetric eigensolver.

**The matrix product.** @ always means matmul, never elementwise.

```viv
A = [[1.0, 2.0], [3.0, 4.0]]
B = [[0.0, 1.0], [1.0, 0.0]]
show A @ B
```

Prints:
```text
[[2, 1],
 [4, 3]]
```

**transpose.**

```viv
A = [[1.0, 2.0], [3.0, 4.0]]
show transpose(A)
```

Prints:
```text
[[1, 3],
 [2, 4]]
```

**The symmetric eigensolver.** eig_sym fills its arguments: values ascending, eigenvectors in columns. Trace = sum of eigenvalues.

```viv
A = [[2.0, 1.0], [1.0, 3.0]]
vals = zeros(2)
vecs = zeros(2, 2)
eig_sym(A, vals, vecs)
show vals
show vals[0] + vals[1]
```

Prints:
```text
vals: [1.38197, 3.61803]
5
```

**Matrix-vector products via matmul.** A 2×2 times a length-2 vector.

```viv
A = [[0.0, 1.0], [1.0, 0.0]]
v = [3.0, 7.0]
show A @ v
```

Prints:
```text
[7, 3]
```

**Chained products.**

```viv
R = [[0.0, 0.0 - 1.0], [1.0, 0.0]]
show R @ R @ R @ R
```

Prints:
```text
[[1, 0],
 [0, 1]]
```

**A 3x3 eigenproblem.** The eigenvalues come out in ascending order and sum to the trace, 9.

```viv
A = [[2.0, 1.0, 0.0], [1.0, 3.0, 1.0], [0.0, 1.0, 4.0]]
vals = zeros(3)
vecs = zeros(3, 3)
eig_sym(A, vals, vecs)
show vals
show vals[0] + vals[1] + vals[2]
```

Prints:
```text
vals: [1.26795, 3, 4.73205]
9
```

**Solvers that return their results.** The same solvers as the forms that fill arrays you pass, which still work: the values ascend, and lin_solve's answer gives back b.

```viv
A = [[2.0, 1.0, 0.0], [1.0, 3.0, 1.0], [0.0, 1.0, 4.0]]
vals, vecs = eig_sym(A)
show vals
x = lin_solve(A, [1.0, 2.0, 3.0])
show A @ x
```

Prints:
```text
vals: [1.26795, 3, 4.73205]
[1, 2, 3]
```

**transpose reverses the axes.** Of a matrix, the rows become columns; of a rank-3 array, all three axes turn around; a vector comes back as it was, as in NumPy. transpose(A, axes) puts them in any order.

```viv
T = transpose(reshape(arange(24), 2, 3, 4))
show shape(T)
show transpose([1.0, 2.0])
```

Prints:
```text
[4, 3, 2]
[1, 2]
```

**The dot product is dot, not @.** @ is the matrix product and refuses two vectors; dot(u, v) sums u[i] * v[i] left to right, and outer(u, v) is the other product of two vectors. The cell prints 11 from dot, then stops at the @ line.

```viv
show dot([1.0, 2.0], [3.0, 4.0])
show [1.0, 2.0] @ [3.0, 4.0]
```

Stops with an error containing: `` `@` of two vectors is refused ``

### Lambdas & higher-order builtins

fn values, quadrature, and exact differentiation.

**fn is a value.**

```viv
square = fn(x) -> x * x
show square(7.0)
```

Prints:
```text
49
```

**Quadrature.** Integral of x² over [0,1]; exact answer 1/3.

```viv
show simpson(fn(x) -> x * x, 0.0, 1.0, 40)
```

Prints:
```text
0.333333
```

**Exact derivatives.** deriv uses forward-mode autodiff, so the slope is exact to rounding: d/dx[x sin x] at 2 against the product rule.

```viv
show deriv(fn(x) -> x * sin(x), 2.0)
show sin(2.0) + 2.0 * cos(2.0)
```

Prints:
```text
0.0770038
0.0770038
```

**deriv knows log.** d/dx log(x) at 3 is 1/3.

```viv
show deriv(fn(x) -> log(x), 3.0)
```

Prints:
```text
0.333333
```

**Two-argument lambdas.**

```viv
add = fn(a, b) -> a + b
show add(3.0, 4.0)
```

Prints:
```text
7
```

**simpson accepts an odd panel count.** Simpson's rule needs an even number of panels, so simpson rounds an odd n up to the next even number: 7 panels give exactly the 8-panel answer. The integrand is exp, which Simpson's rule does not integrate exactly, so the two would differ otherwise.

```viv
odd = simpson(fn(x) -> exp(x), 0.0, 1.0, 7)
even = simpson(fn(x) -> exp(x), 0.0, 1.0, 8)
show odd == even
show trapezoid(fn(x) -> x * x, 0.0, 1.0, 1000)
```

Prints:
```text
true
0.333333
```

**deriv composes through the chain rule.** d/dx exp(sin x) at 1 = exp(sin 1) · cos 1.

```viv
show deriv(fn(x) -> exp(sin(x)), 1.0)
show exp(sin(1.0)) * cos(1.0)
```

Prints:
```text
1.25338
1.25338
```

**The functions deriv can differentiate.** Besides + - * and /, deriv knows exp, sqrt, pow, sin, cos, log and tan. The last line uses the quotient rule.

```viv
show deriv(fn(x) -> sqrt(x), 4.0)
show deriv(fn(x) -> pow(x, 3.0), 2.0)
show deriv(fn(x) -> tan(x), 1.0)
show deriv(fn(x) -> sin(x) / x, 2.0)
```

Prints:
```text
0.25
12
3.42552
-0.435398
```

**deriv wants a constant exponent.** pow(x, c) differentiates only for literal c.

```viv
show deriv(fn(x) -> pow(x, x), 2.0)
```

Stops with an error containing: `constant exponent`

**Quadrature error shrinks with panels.** For the integral of sin on [0, pi], exactly 2, Simpson's error falls from 0.0046 with 4 panels to 6.5e-8 with 64.

```viv
pi_v = 3.141592653589793
show abs(simpson(fn(x) -> sin(x), 0.0, pi_v, 4) - 2.0)
show abs(simpson(fn(x) -> sin(x), 0.0, pi_v, 64) - 2.0)
```

Prints:
```text
0.00455975
6.453e-08
```

**An array of functions.** An array holds numbers, or it holds functions. A function array can be written as a literal, indexed, passed to len, and passed to a def. That is enough to keep a list of methods and pick one while the program runs: because the index is not known in advance, the compiled call goes through a function pointer instead of being fixed at compile time.

```viv
def sq(x) {
    return x * x
}
def cube(x) {
    return x * x * x
}
methods = [sq, cube]
for i to len(methods) {
    show apply(methods[i], 3.0)
}
```

Prints:
```text
9
27
```

**Methods take the function they work on.** Each method's first argument is the function it works on, so one program can hand the same problem to several methods and compare their answers. Here spread returns the gap between the largest and smallest answer.

```viv
def target(x) {
    return x * x - 2.0
}
def left(f, a, b) {
    return f(a)
}
def right(f, a, b) {
    return f(b)
}
def spread(ms, f, a, b) {
    lo = apply(ms[0], f, a, b)
    hi = lo
    for i to len(ms) {
        v = apply(ms[i], f, a, b)
        lo = min(lo, v)
        hi = max(hi, v)
    }
    return hi - lo
}
show spread([left, right], target, 0.0, 3.0)
```

Prints:
```text
9
```

**Every function in an array takes the same arguments.** Every function in one array must have the same C signature: the same number of arguments, in the same positions. Mixing them is an error at the line that builds the array. The interpreter refuses it too, so a program cannot pass in the browser and then fail when compiled to C.

```viv
def one(x) {
    return x
}
def two(x, y) {
    return x + y
}
ms = [one, two]
show len(ms)
```

Stops with an error containing: `the same number of arguments`

### Scope & mutation

Where names live, what functions can see, and which values are shared.

**Functions read globals.** A top-level scalar is visible inside every def.

```viv
g = 5.0
def f() {
    return g + 1.0
}
show f()
```

Prints:
```text
6
```

**Assigning inside a def never touches a global scalar.** Assignment in a function body creates a local variable; the global g keeps its value, 5.

```viv
g = 5.0
def f() {
    g = 99.0
    return g
}
show f()
show g
```

Prints:
```text
99
g: 5
```

**A def can change a global array.** A def that stores into an element of a global array changes the global array; no return is needed.

```viv
tape = zeros(3)
def store(v) {
    tape[0] = v
    return 0.0
}
x = store(42.0)
show tape
```

Prints:
```text
tape: [42, 0, 0]
```

**Assignment aliases arrays.** b = a makes b refer to the same array; it does not copy the data. Changing b changes a.

```viv
a = [1.0, 2.0]
b = a
b[0] = 99.0
show a
```

Prints:
```text
a: [99, 2]
```

**Array parameters alias too.** A function mutating its array parameter mutates the caller's array.

```viv
def poke(v) {
    v[0] = 77.0
    return 0.0
}
a = [1.0, 2.0]
x = poke(a)
show a
```

Prints:
```text
a: [77, 2]
```

**Closures see the current global.** Capture is by reference for globals: rebinding y changes what f computes.

```viv
y = 2.0
f = fn(x) -> x * y
y = 10.0
show f(3.0)
```

Prints:
```text
30
```

### Prose in programs

Comments, docstrings, narration, and annotation markers.

**Comments are ignored by the compiler.**

```viv
# this line compiles to nothing
show 1.0  # and so does this tail
```

Prints:
```text
1
```

**narrate reports progress.** narrate prints a progress line while the program runs, kept separate from the data output.

```viv
narrate "building the matrix"
M = zeros(2, 2)
narrate "done"
show M[0, 0]
```

Prints:
```text
narrate: building the matrix
narrate: done
0
```

**Annotation markers are comments.** The compiler ignores them; viv doc extracts them for code browsers.

```viv
#| level-1 note: what the next line does
answer = 42.0
show answer
```

Prints:
```text
answer: 42
```

**The variable glossary marker.** #: name description documents what a variable IS. Inside a def it belongs to that function. viv doc prints the glossary, and the code viewer shows the description when you hover over any use of the name.

```viv
#| One Newton update for sqrt(a).
#: x  the current estimate
#: a  the number whose root we want
def newton_sqrt_step(x, a) {
    #: fx  the residual x^2 - a; zero at the answer
    fx = x * x - a
    return x - fx / (2.0 * x)
}

show newton_sqrt_step(1.5, 2.0)
```

Prints:
```text
1.41667
```

**Narrating phases.** In a long program, narrate the start of each phase so the output shows where the run is.

```viv
narrate "phase 1: building"
M = zeros(2, 2)
M[0, 0] = 1.0
narrate "phase 2: solving"
v = [1.0, 0.0]
w = M @ v
narrate "phase 3: reporting"
show w[0]
```

Prints:
```text
narrate: phase 1: building
narrate: phase 2: solving
narrate: phase 3: reporting
1
```

### Plotting

plot prints a structured record, and the page draws the chart from it.

**plot emits a structured record.** figure declares the chart (title, x label, y label) and plot adds a series. Both print records on stdout; the page draws the chart from them and adds nothing of its own.

```viv
x = linspace(0.0, 3.0, 7)
y = zeros(7)
for i to 7 {
    y[i] = x[i] * x[i]
}
figure("A parabola", "x", "x squared")
plot(x, y, "x squared")
```

Prints:
```text
figure: A parabola ; x ; x squared
plot: x squared ; [0, 0.5, 1, 1.5, 2, 2.5, 3] ; [0, 0.25, 1, 2.25, 4, 6.25, 9]
```

**A matrix as a heatmap with plot2d.** plot2d(Z) or plot2d(Z, "label") renders a rank-2 array as a color map (viridis) with a colorbar. The format of the printed record is fixed, and this example's test fails if it changes.

```viv
Z = zeros(2, 3)
for i to 2 {
    for j to 3 {
        Z[i, j] = i * 3.0 + j
    }
}
figure("A 2 x 3 field", "column", "row")
plot2d(Z, "a tiny field")
```

Prints:
```text
figure: A 2 x 3 field ; column ; row
plot2d: a tiny field ; 2 3 ; [0, 1, 2, 3, 4, 5]
```

**A curve in three dimensions with plot3d.** plot3d(x, y, z, "label") draws a rotatable 3-D curve — drag the rendered figure to spin it. Three equal-length rank-1 arrays.

```viv
n = 3
xs = linspace(0.0, 1.0, n)
ys = zeros(n)
zs = zeros(n)
for i to n {
    ys[i] = xs[i] * xs[i]
    zs[i] = 1.0 - xs[i]
}
figure("An arc in space", "x", "y")
plot3d(xs, ys, zs, "a short arc")
```

Prints:
```text
figure: An arc in space ; x ; y
plot3d: a short arc ; [0, 0.5, 1] ; [0, 0.25, 1] ; [1, 0.5, 0]
```

**A histogram.** histogram counts into equal bins; bars draws the counts at the bin centers. A plot on the same x afterwards draws a curve over the bars, which is how a measured distribution is shown against a theoretical one.

```viv
v = [0.1, 0.4, 0.45, 0.8, 0.9, 0.95, 1.7, 2.2]
h = histogram(v, 0.0, 2.0, 4)
c = linspace(0.25, 1.75, 4)
figure("Counts per bin", "value", "count")
bars(c, h, "counts")
plot(c, 4.0 - 2.0 * c, "a line over them")
```

Prints:
```text
figure: Counts per bin ; value ; count
bars: counts ; [0.25, 0.75, 1.25, 1.75] ; [3, 3, 0, 1]
plot: a line over them ; [0.25, 0.75, 1.25, 1.75] ; [3.5, 2.5, 1.5, 0.5]
```

**Same x merges into one chart.** Consecutive plots sharing an x-array become labeled series on one figure.

```viv
n = 41
x = linspace(0.0, 6.28, n)
s1 = zeros(n)
s2 = zeros(n)
for i to n {
    s1[i] = sin(x[i])
    s2[i] = cos(x[i])
}
figure("Sine and cosine", "x", "value")
plot(x, s1, "sin")
plot(x, s2, "cos")
```

Prints:
```text
figure: Sine and cosine ; x ; value
plot: sin ; [0, 0.157, 0.314, 0.471, 0.628, 0.785, 0.942, 1.099, 1.256, 1.413, 1.57, 1.727, 1.884, 2.041, 2.198, 2.355, 2.512, 2.669, 2.826, 2.983, 3.14, 3.297, 3.454, 3.611, 3.768, 3.925, 4.082, 4.239, 4.396, 4.553, 4.71, 4.867, 5.024, 5.181, 5.338, 5.495, 5.652, 5.809, 5.966, 6.123, 6.28] ; [0, 0.1563558122752478, 0.30886552009893214, 0.4537776270755451, 0.5875275257138919, 0.706825181105366, 0.8087360605531302, 0.8907533184119663, 0.95085946050647, 0.9875759712809227, 0.9999996829318346, 0.987824991900309, 0.9513513762338286, 0.8914760289023131, 0.8096717882771642, 0.7079509086484322, 0.5888155619677953, 0.4551962883957966, 0.31037990967204165, 0.1579286785248665, 0.0015926529164868282, -0.15478254942208766, -0.3073503470745565, -0.4523578147279585, -0.5862379991700271, -0.7056976606684764, -0.8077982814337493, -0.8900283484858137, -0.9503651328813763, -0.9873244456307223, -0.999997146387718, -0.9880715068572294, -0.9518408788156858, -0.8921964781236668, -0.8106054622323361, -0.7090748404422167, -0.5901021046645756, -0.45661379509021477, -0.31189351195256765, -0.15950114418129782, -0.0031853017931379904]
plot: cos ; [0, 0.157, 0.314, 0.471, 0.628, 0.785, 0.942, 1.099, 1.256, 1.413, 1.57, 1.727, 1.884, 2.041, 2.198, 2.355, 2.512, 2.669, 2.826, 2.983, 3.14, 3.297, 3.454, 3.611, 3.768, 3.925, 4.082, 4.239, 4.396, 4.553, 4.71, 4.867, 5.024, 5.181, 5.338, 5.495, 5.652, 5.809, 5.966, 6.123, 6.28] ; [1, 0.9877007947590948, 0.951105719935495, 0.891114956201324, 0.8092041809880323, 0.7073882691671998, 0.5881717303313752, 0.4544871018390625, 0.3096228130570242, 0.1571422952248764, 0.0007963267107332633, -0.15556923017471813, -0.30810803127799236, -0.45306786455514825, -0.5868829485238594, -0.7062616448200052, -0.8082674272693101, -0.8903911157635676, -0.9506125981029155, -0.9874505215449477, -0.9999987317275395, -0.9879485626258074, -0.9515964292456911, -0.8918365362859643, -0.8101388821239993, -0.7085130991922731, -0.5894590202148757, -0.45590518629602694, -0.3111368094638822, -0.15871496167654411, -0.0023889781122815386, 0.153995770516219, 0.3065924679690989, 0.45164747804424354, 0.5855926780613815, 0.7051332290084228, 0.8073286233439505, 0.8896650168087489, 0.9501170649987787, 0.9871977436181957, 0.9999949269133752]
```

**Different x means different charts.** Each chart gets its own figure.

```viv
a_x = linspace(0.0, 1.0, 5)
b_x = linspace(0.0, 2.0, 9)
a_y = zeros(5)
b_y = zeros(9)
for i to 5 { a_y[i] = a_x[i] }
for i to 9 { b_y[i] = b_x[i] * b_x[i] }
figure("A line", "x", "y")
plot(a_x, a_y, "line")
figure("A parabola", "x", "y")
plot(b_x, b_y, "parabola")
```

Prints:
```text
figure: A line ; x ; y
plot: line ; [0, 0.25, 0.5, 0.75, 1] ; [0, 0.25, 0.5, 0.75, 1]
figure: A parabola ; x ; y
plot: parabola ; [0, 0.25, 0.5, 0.75, 1, 1.25, 1.5, 1.75, 2] ; [0, 0.0625, 0.25, 0.5625, 1, 1.5625, 2.25, 3.0625, 4]
```

**figure is required.** A chart's title and axis labels are part of the program. Nothing outside it can add them, so a plot with no figure declared does not compile — in any engine.

```viv
x = linspace(0.0, 1.0, 5)
plot(x, x, "a line")
```

Stops with an error containing: `` needs a `figure(title, xlabel, ylabel)` before it ``

### Animation

An animation is a Viv program: it creates named objects, and transitions reveal, morph and move them. The program prints records; the page plays them back with a scrubber.

**A stage, an object, a transition.** stage opens an animation and gives it a coordinate window. Objects are created hidden; draw reveals a curve by drawing it on. Each call prints one record, which is what this page plays.

```viv
x = linspace(0.0, 1.0, 5)
stage("A parabola, drawn on", 0.0, 1.0, 0.0, 1.0)
axes("x", "y")
curve("f", x, x * x)
draw("f", 1.5)
```

Prints:
```text
stage: A parabola, drawn on ; 0 1 0 1
axes: x ; y
curve: f ; [0, 0.25, 0.5, 0.75, 1] ; [0, 0.0625, 0.25, 0.5625, 1]
draw: f ; 1.5
```

**Morph, with a caption.** morph moves a curve's points smoothly to new data over the given seconds. say sets the caption until the next say; a number given after the text is appended, which is how a caption reports a computed value.

```viv
x = linspace(0.0, 1.0, 3)
stage("A line becomes a parabola", 0.0, 1.0, 0.0, 1.0)
curve("f", x, x)
draw("f", 1.0)
say("now square it; the midpoint falls to", 0.5 * 0.5)
morph("f", x, x * x, 2.0)
```

Prints:
```text
stage: A line becomes a parabola ; 0 1 0 1
curve: f ; [0, 0.5, 1] ; [0, 0.5, 1]
draw: f ; 1
say: now square it; the midpoint falls to ; 0.25
morph: f ; [0, 0.5, 1] ; [0, 0.25, 1] ; 2
```

**Running transitions together.** Transitions run one after another. Between together() and apart() they start at the same moment, and the group lasts as long as its longest member. A label's text may carry TeX between dollars: write two backslashes for one.

```viv
x = linspace(0.0, 1.0, 3)
stage("Together", 0.0, 1.0, 0.0, 1.0)
curve("f", x, x)
hline("half", 0.5)
label("name", "$y = \\tfrac12$", 0.8, 0.56)
together()
draw("f", 1.0)
fade_in("half", 0.5)
fade_in("name", 0.5)
apart()
wait(1.0)
```

Prints:
```text
stage: Together ; 0 1 0 1
curve: f ; [0, 0.5, 1] ; [0, 0.5, 1]
hline: half ; 0.5
label: name ; $y = \tfrac12$ ; 0.8 0.56
together:
draw: f ; 1
fade_in: half ; 0.5
fade_in: name ; 0.5
apart:
wait: 1
```

**A simulation drawn as it runs, with keyframes.** frame(dt) advances the clock. An object set again later moves there over the time between: a ball drawn at two places with a frame between them glides from one to the other. In a keyframed animation objects are visible from creation; shapes (circle, rect, segment, poly) and dots with a radius are sized in the stage's units.

```viv
stage("A bouncing ball", 0.0, 4.0, 0.0, 2.0)
hline("floor", 0.0)
fall = sqrt(2.0 / 4.9)
for k to 9.0 {
    t = 0.2 * k
    # time from the nearest top of a bounce: the ball falls from y = 2 and comes back up
    s = t - 2.0 * fall * floor(t / (2.0 * fall) + 0.5)
    circle("ball", 1.0 + t, 2.0 - 4.9 * s * s, 0.1)
    frame(0.2)
}
```

Prints:
```text
stage: A bouncing ball ; 0 4 0 2
hline: floor ; 0
circle: ball ; 1 2 0.1
frame: 0.2
circle: ball ; 1.2 1.804 0.1
frame: 0.2
circle: ball ; 1.4 1.216 0.1
frame: 0.2
circle: ball ; 1.6 0.236 0.1
frame: 0.2
circle: ball ; 1.8 0.881585 0.1
frame: 0.2
circle: ball ; 2 1.62198 0.1
frame: 0.2
circle: ball ; 2.2 1.97038 0.1
frame: 0.2
circle: ball ; 2.4 1.92677 0.1
frame: 0.2
circle: ball ; 2.6 1.49117 0.1
frame: 0.2
```

### Testing

test blocks and the assert family: a program checks its own results and prints which tests passed and which failed.

**A passing test.** assert_near(a, b, tol) inside a named test block. What prints is the name — one line per test, however many assertions are inside.

```viv
test "pythagoras" {
    assert_near(sqrt(3.0 * 3.0 + 4.0 * 4.0), 5.0, 1.0e-12)
}
```

Prints:
```text
PASS pythagoras
```

**A failing test says what failed.** The FAIL line names the test and shows the failing assertion's values, so you see which check failed and by how much.

```viv
test "wrong on purpose" {
    assert_near(2.0 + 2.0, 5.0, 1.0e-9)
}
```

Prints:
```text
FAIL wrong on purpose -- assert_near: |4 - 5| = 1 > 1e-09
```

**Test a function against a property.** Check a function against a mathematical fact it must satisfy (here, the partial sums of 1/k! converge to e), not only against fixed values.

```viv
def esum(n) {
    s = 0.0
    for i to n {
        s += 1.0 / gamma(i + 1.0)
    }
    return s
}
test "partial sums of 1/k! approach e" {
    assert_near(esum(20.0), exp(1.0), 1.0e-12)
}
```

Prints:
```text
PASS partial sums of 1/k! approach e
```

### Errors, exactly as the compiler says them

The error messages are part of the language. Each one below is what the compiler prints, and the site's tests fail if a message changes.

**A string is not changed in place.** Strings are immutable, as in Python; build a new one instead.

```viv
s = "hi"
s[0] = "H"
```

Stops with an error containing: `'str' object does not support item assignment`

**deriv cannot differentiate a call to a def.** Write the expression inline in the fn instead of calling f.

```viv
def f(x) {
    return x * x
}
show deriv(fn(a) -> f(a), 1.0)
```

Stops with an error containing: `cannot differentiate through`

**A def does not change a record it is given.** A record is a value: change a copy and return it. An array inside it can still be written in place.

```viv
record P { x: Num }
def f(p: P) {
    p.x = 2.0
    return p.x
}
show f(P { x: 1.0 })
```

Stops with an error containing: `` `p` is a parameter, and a def does not change the values it is given ``

**Unknown function.**

```viv
show frobnicate(1.0)
```

Stops with an error containing: `unknown function`

**Arity is checked.**

```viv
def f(a, b) {
    return a + b
}
show f(1.0)
```

Stops with an error containing: `expects 2 argument`

**Out of bounds says where.**

```viv
a = [1.0, 2.0]
show a[5]
```

Stops with an error containing: `index 5 out of bounds for dimension 0 (size 2)`

**An index before the start stops.** A negative index counts from the end, as in Python: -1 is the last element. Two elements have no element -3.

```viv
a = [1.0, 2.0]
show a[-1]
show a[0.0 - 3.0]
```

Stops with an error containing: `index -3 out of bounds for dimension 0 (size 2)`

**Shape mismatch.**

```viv
show [1.0, 2.0] + [1.0, 2.0, 3.0]
```

Stops with an error containing: `array shape mismatch: [2] vs [3]`

**break outside a loop.**

```viv
break
show 1.0
```

Stops with an error containing: `` `break`/`continue` outside a loop ``

**Undefined variable.**

```viv
show nonexistent + 1.0
```

Stops with an error containing: `undefined variable`

**plot wants arrays.**

```viv
figure("Nothing to see", "x", "y")
plot(1.0, 2.0)
```

Stops with an error containing: `` `plot` takes two arrays ``

**Keywords cannot name functions.** record, test, show and the other keywords are reserved.

```viv
def record(v) {
    return v
}
show record(1.0)
```

Stops with an error containing: `expected function name`

**A number is not a condition.** if and while take a Bool. Write the comparison you mean.

```viv
count = 3.0
if count {
    show count
}
```

Stops with an error containing: `condition must be a boolean`

**No arithmetic on a Bool.** true is not 1. Every engine refuses this.

```viv
x = true + 1.0
```

Stops with an error containing: `boolean`

**An array is not a condition.** A comparison of two arrays is a mask, not a Bool; the note says to write all(a == b).

```viv
a = [1.0, 2.0]
if a == a {
    show 1.0
}
```

Stops with an error containing: `condition must be a boolean, and this is an array`

**% by zero stops the program.** `/` by zero gives inf, as IEEE floats do, but `%` and `//` by zero have no answer and stop.

```viv
show 7.0 % 0.0
```

Stops with an error containing: `` division by zero in `a % b` ``

### Gotchas

Float facts and Viv-specific surprises, each shown by a program that runs.

**Float equality is a trap.** 0.1 + 0.2 is not 0.3 in binary floating point; compare with a tolerance instead. (Integers are another matter: a double holds them exactly up to 2^53, so `n == 1.0` on a counter is an exact test.)

```viv
show (0.1 + 0.2) == 0.3
show abs((0.1 + 0.2) - 0.3) < 1.0e-9
```

Prints:
```text
false
true
```

**The printer rounds to 6 significant digits.** show uses %g formatting; the stored value is full double precision even when the print looks tidy.

```viv
show 0.1 + 0.2
show (0.1 + 0.2) - 0.3
```

Prints:
```text
0.3
5.55112e-17
```

**eig_sym does not check symmetry.** Feed it a non-symmetric matrix and it computes something — for [[1,5],[0,2]] the true eigenvalues are 1 and 2, and this is not them. Symmetrize first.

```viv
A = [[1.0, 5.0], [0.0, 2.0]]
v = zeros(2)
w = zeros(2, 2)
eig_sym(A, v, w)
show v
```

Prints:
```text
v: [-3.52494, 6.52494]
```

**Copying an array means allocating one.** Because assignment aliases, copying the data takes copy(a), which allocates a new array with the same entries.

```viv
a = [1.0, 2.0]
b = copy(a)
b[0] = 99.0
show a
show b
```

Prints:
```text
a: [1, 2]
b: [99, 2]
```

**Adding 1 to 1e16 loses the 1.** In exact arithmetic both lines would print 1; here both print 0. A double carries about 16 significant digits, so the 1 is rounded away when it meets 1e16, whichever way the sum is grouped.

```viv
big = 1.0e16
show (big + 1.0) - big
show big + (1.0 - big)
```

Prints:
```text
0
0
```

**Catastrophic cancellation.** At x = 1e12, sqrt(x+1) - sqrt(x) subtracts two nearly equal numbers and keeps only 5 correct digits; the rewrite 1 / (sqrt(x+1) + sqrt(x)) has no subtraction and gets them all.

```viv
x = 1.0e12
naive = sqrt(x + 1.0) - sqrt(x)
rewritten = 1.0 / (sqrt(x + 1.0) + sqrt(x))
show naive
show rewritten
```

Prints:
```text
naive: 5.00004e-07
rewritten: 5e-07
```

### Stencils and PDE steps

A finite-difference stencil is arithmetic on shifted copies of an array, and an implicit time step is a tridiagonal solve. These are the pieces the PDE solvers are built from.

**A periodic second difference.** roll(u, k) moves element i to i + k, wrapping around, so roll(u, 1)[i] is u[i - 1]. The second difference of i squared is exactly 2; the two ends show the wrap.

```viv
u = [0.0, 1.0, 4.0, 9.0, 16.0]
d2 = roll(u, 1.0) - 2.0 * u + roll(u, -1.0)
show d2
```

Prints:
```text
d2: [17, 2, 2, 2, -23]
```

**roll along one axis of a grid.** A third argument picks the axis: 0 shifts rows, 1 shifts along each row. A five-point Laplacian is four of these.

```viv
A = [[1.0, 2.0, 3.0], [4.0, 5.0, 6.0]]
B = roll(A, 1.0, 1.0)
C = roll(A, 1.0, 0.0)
show B[0, 0]
show C[0, 0]
```

Prints:
```text
3
4
```

**A tridiagonal solve.** tridiag_solve(sub, diag, super, rhs) is the Thomas algorithm, O(n). Here the matrix is the 1-D Laplacian (-1, 2, -1), and the right-hand side is chosen so the answer is all ones. sub[0] and super[n-1] are never read.

```viv
n = 4
sub = zeros(n) - 1.0
diag = zeros(n) + 2.0
sup = zeros(n) - 1.0
x = tridiag_solve(sub, diag, sup, [1.0, 0.0, 0.0, 1.0])
show x
```

Prints:
```text
x: [1, 1, 1, 1]
```

**One implicit heat step.** Backward Euler for u_t = u_xx with dt / h^2 = 1: solve (I - r D2) u_new = u_old. The spike spreads symmetrically, and no time step is too large for it.

```viv
n = 5
r = 1.0
u = [0.0, 0.0, 1.0, 0.0, 0.0]
off = zeros(n) - r
u = tridiag_solve(off, zeros(n) + (1.0 + 2.0 * r), off, u)
show u
```

Prints:
```text
u: [0.0555556, 0.166667, 0.444444, 0.166667, 0.0555556]
```

**The max-norm error.** max and min of ONE array reduce it to its largest or smallest entry. With two arguments they stay elementwise.

```viv
x = linspace(0.0, 1.0, 5)
err = max(abs(x * x - x))
show err
show min(x * x - x)
```

Prints:
```text
err: 0.25
-0.25
```

**An array to a power.** pow(a, s) raises every entry to the same scalar power.

```viv
c = pow([1.0, 8.0, 27.0], 1.0 / 3.0)
show c[1]
show sum(pow([1.0, 2.0, 3.0], 2.0))
```

Prints:
```text
2
14
```

**Negating an array.** Unary minus works on arrays, real or complex.

```viv
a = [1.0, -2.0, 0.5]
b = -a
show b
```

Prints:
```text
b: [-1, 2, -0.5]
```

### Small algorithms

Classic numerical algorithms, each a few lines long.

**Newton's square root.** Repeat x = (x + a/x) / 2; each step roughly doubles the number of correct digits.

```viv
def newton_sqrt(a) {
    x = a
    for i to 20 {
        x = 0.5 * (x + a / x)
    }
    return x
}
show newton_sqrt(2.0)
show sqrt(2.0)
```

Prints:
```text
1.41421
1.41421
```

**Bisection.** Root of cos(x) on [1, 2]. Each step halves the interval, so 40 halvings leave it 2^-40 ≈ 9e-13 wide and give pi/2 to about twelve digits; show prints the first six.

```viv
def bisect_cos(lo, hi) {
    for i to 40 {
        mid = 0.5 * (lo + hi)
        if cos(lo) * cos(mid) < 0.0 {
            hi = mid
        } else {
            lo = mid
        }
    }
    return 0.5 * (lo + hi)
}
show bisect_cos(1.0, 2.0)
```

Prints:
```text
1.5708
```

**gcd via mod.** Euclid's algorithm. The floats hold whole numbers exactly, so mod is exact.

```viv
def gcd(a, b) {
    if b < 0.5 {
        return a
    }
    return gcd(b, mod(a, b))
}
show gcd(1071.0, 462.0)
```

Prints:
```text
21
```

**Horner evaluation.** A polynomial as an array of coefficients, highest power first.

```viv
def horner(coeffs, x) {
    acc = 0.0
    for i to len(coeffs) {
        acc = acc * x + coeffs[i]
    }
    return acc
}
show horner([2.0, 0.0 - 3.0, 5.0], 2.0)
```

Prints:
```text
7
```

**The power method.** The dominant eigenvalue of a symmetric 2x2, found by repeated matrix-vector products; compare eig_sym.

```viv
def normalize(v) {
    n = sqrt(v[0] * v[0] + v[1] * v[1])
    return v / n
}
def rayleigh(A, v) {
    w = A @ v
    return v[0] * w[0] + v[1] * w[1]
}
A = [[2.0, 1.0], [1.0, 3.0]]
v = [1.0, 0.0]
for i to 30 {
    v = normalize(A @ v)
}
show rayleigh(A, v)
```

Prints:
```text
3.61803
```

**Compensated summation.** Kahan's method carries each addition's rounding error into the next one. Adding 0.1 a million times this way gives exactly 100000.

```viv
def kahan_sum(n, val) {
    s = 0.0
    comp = 0.0
    for i to n {
        y = val - comp
        t = s + y
        comp = (t - s) - y
        s = t
    }
    return s
}
show kahan_sum(1000000.0, 0.1) - 100000.0
```

Prints:
```text
0
```

**Fibonacci two ways.** Naive recursion, whose number of calls grows exponentially with n, against a loop that takes n steps.

```viv
def fib_rec(n) {
    if n < 1.5 {
        return n
    }
    return fib_rec(n - 1.0) + fib_rec(n - 2.0)
}
def fib_loop(n) {
    a = 0.0
    b = 1.0
    for i to n {
        t = a + b
        a = b
        b = t
    }
    return a
}
show fib_rec(15.0)
show fib_loop(15.0)
```

Prints:
```text
610
610
```

## Posts

A post is one `.viv` file whose first lines are a fenced header: a line that
is exactly `###`, then `key: value` fields, a blank line, the prose notes, and
`###` again. The lexer skips the block. The site build runs every post and
refuses it unless:

- the header has `post:` (title), `date:` (YYYY-MM-DD) and `tags:` (comma separated); `summary:` is at most 600 characters;
- at least one `test "..." { ... }` block runs, and none fails (a `kind: problem` post is the exception: it must have a failing test);
- a demo (the default kind) or an `kind: article` draws a figure: `figure(...)` then `plot(...)` or another series, or `stage(...)` for an animation;
- it runs in the browser in 10 s or less, and the figure records are at most 160 KB;
- a program written by a model says so with `agent: <model name>`. Such a post starts in the AI lane, and only other people's reviews move it up.

Other kinds: `kind: ladder` (built from `rung "title" { ... }` blocks), `kind: problem`
(stubs and tests that fail), `kind: solution` with `solves: <problem slug>` (the
problem's tests, word for word, all passing). `from: <slug>` with `edge:` (fork,
fix, extend, rewrite, adversarial-review, validate) derives a post from another,
and the build checks the claim by running both.

```viv
###
post: The trapezoid rule's error falls as h squared
date: 2026-09-24
tags: quadrature, calculus
summary: Halving the step divides the trapezoid rule's error on sin by four.
agent: claude-opus-5.5

The integral of sin from 0 to pi is exactly 2. The trapezoid rule's error is
proportional to h^2, so doubling the number of panels should cut it by four.
###

def trap(a, b, n) {
    h = (b - a) / n
    s = 0.5 * (sin(a) + sin(b))
    for k to n - 1 {
        s = s + sin(a + (k + 1) * h)
    }
    return s * h
}

ns = [4.0, 8.0, 16.0, 32.0, 64.0]
errs = zeros(5)
for k to 5 {
    errs[k] = abs(trap(0.0, 3.141592653589793, ns[k]) - 2.0)
}
figure("trapezoid error on sin", "panels n", "|error|")
plot(ns, errs, "error")
ratio = errs[3] / errs[4]
show ratio

test "halving h divides the error by four" {
    assert_near(ratio, 4.0, 1.0e-2)
}
```