Rauzy Fractal

DEFINITION

σ: 1→12, 2→13, 3→1
M (letter counts) = [[1,1,1],[1,0,0],[0,1,0]]
plot Σᵢ₌₁ⁿ e_{wᵢ}, projected onto the contracting eigenplane of M

NOTES

Iterate the substitution 1→12, 2→13, 3→1 from a single "1" and it converges to one infinite word — the Tribonacci word, the three-letter sibling of the Fibonacci word two plates over. Walk that word letter by letter, keeping a running tally of how many 1s, 2s and 3s have gone by, and project each tally into the plane spanned by the substitution matrix’s two complex eigenvectors — the ones whose eigenvalues sit inside the unit circle, so the projected walk stays bounded instead of running off to infinity the way the raw tally does. Rauzy did this in 1982 while asking a different question — whether a substitution’s dynamics could be realised as a rotation on a compact group — and the picture that fell out was this one: a single connected tile with a fractal shore, built from three self-similar copies of itself, that tiles the plane by translation. The letter tally’s escaping direction is governed by the third eigenvalue, the Tribonacci constant β ≈ 1.8393 sitting outside the unit circle — a Pisot number, and the reason the other two directions contract at all.

PROVENANCE

Origin

G. Rauzy, "Nombres algébriques et substitutions", Bulletin de la Société Mathématique de France 110, 1982, 147–178 (numdam.org/item/BSMF_1982__110__147_0)

Standing

Public domain — implemented here from the substitution and the incidence-matrix eigenproblem set out in the paper, not from any existing picture of the tile

The tile, not just the picture

Rauzy’s own paper already proves the Tribonacci case tiles the plane by translation and decomposes into three rescaled copies of itself. The stronger claim, that every Pisot substitution has such a tile, was established afterward by P. Arnoux & S. Ito, "Pisot substitutions and Rauzy fractals", Bulletin of the Belgian Mathematical Society — Simon Stevin 8(2), 2001, 181–207. This plate renders the original 1982 case, not the later generalisation.

Checked, not just plotted

The three letters’ long-run frequencies in the word are the Perron eigenvector of M — β, 1, 1/β normalised to sum to 1, i.e. 0.5437 / 0.2956 / 0.1607 — and a direct count over the first 20,000 letters here lands on 0.5437 / 0.2956 / 0.1607 exactly, independent evidence the projection below is the one the eigenproblem actually specifies rather than a shape that merely looks plausible

Constants

points is a cost dial, not a look dial, on fibword’s and gosper’s precedent: more prefix points only ever resolve more of the same fixed tile, never a different or a dead one, so it is locked from regenerate the way their generation counts are

Source

doi:10.24033/bsmf.1957

HOUDINI · VEX

The same published mathematics as a Detail Wrangle body. Paste it into a Wrangle with Run Over set to Detail; every constant is the published value plus a tweak channel, so Create Spare Parameters gives a slider that starts where the paper does.

// FORMA — PL. 97 · RAUZY FRACTAL — Gérard Rauzy, 1982
//   σ: 1→12, 2→13, 3→1
//   M (letter counts) = [[1,1,1],[1,0,0],[0,1,0]]
//   plot Σᵢ₌₁ⁿ e_{wᵢ}, projected onto the contracting eigenplane of M
// Paste into a Detail Wrangle (Run Over: Detail), no inputs needed.
// Written from the published mathematics, not adapted from any code.
// Constants arrive at their published values. Press the node's Create
// Spare Parameters button and every tweak becomes a slider — starting
// at 0, the published figure, and moving in the constant's own units.
// https://forma-gen.com/#plate=rauzy

float p_points  = 5500 + chf('points_tweak');     // prefix points plotted · live 2000 .. 9000
float p_bearing = 205 + chf('bearing_tweak');     // bearing (deg) · live 0 .. 360

// The plate's own colour: FORMA's FRACTALS accent as a cosine ramp,
// brightest near t = 0 and t = 1, near-black around t = 0.5.
vector forma_ramp(float t){
  return set(
    0.46 + 0.5 * cos(6.28318530718 * (t + 0)),
    0.1389 + 0.151 * cos(6.28318530718 * (t + 0.05)),
    0.1912 + 0.2078 * cos(6.28318530718 * (t + 0.1)));
}

// Rauzy's tile from the 1982 paper's own route: iterate 1→12, 2→13, 3→1
// from a single 1 until the Tribonacci word covers the prefix count, then
// project every abelianised prefix through the left eigenvector of the
// incidence matrix for its contracting complex eigenvalue — orthogonal to
// the expanding direction by the eigenpair identity, so the escaping
// component projects away with nothing separate to get wrong. β by Newton
// from the Tribonacci polynomial, α from the deflated quadratic, exactly as
// the plate derives them. The word is grown as an int array — VEX strings
// have no cheap indexing and substr does not exist. One point per prefix,
// the plate's own single bright tone; its faint underlay and HILITE head
// are the reveal's, not the tile's. The plate adds its page phase to the
// bearing; a port has no page, so the dial alone turns the tile.
// Deterministic by construction — nothing here is random.
int N = int(rint(p_points));

// β, the Tribonacci constant, by Newton's method — cubic convergence
float beta = 2.0;
for (int i = 0; i < 30; i++){
    float f  = beta * beta * beta - beta * beta - beta - 1.0;
    float fp = 3.0 * beta * beta - 2.0 * beta - 1.0;
    beta -= f / fp;
}
// α, the contracting complex root of the deflated quadratic x² + (β−1)x + 1/β
float br = beta - 1.0;
float disc = 4.0 / beta - br * br;             // negative: α is complex
float aRe = -br / 2.0, aIm = sqrt(max(0.0, disc)) / 2.0;
// the left eigenvector u for α: u·M = αu, so u = (1, α−1, 1/α)
float u2Re = aRe - 1.0, u2Im = aIm;
float aAbs2 = aRe * aRe + aIm * aIm;
float u3Re = aRe / aAbs2, u3Im = -aIm / aAbs2; // 1/α = conj(α)/|α|²

// the fixed point of the substitution, grown until it covers the prefix
int w[];
push(w, 1);
while (len(w) < N){
    int nx[];
    foreach (int c; w){
        if (c == 1){ push(nx, 1);  push(nx, 2); }
        else if (c == 2){ push(nx, 1);  push(nx, 3); }
        else push(nx, 1);
    }
    w = nx;
}

float ang = p_bearing * M_PI / 180.0;
float co = cos(ang), si = sin(ang);
vector col = forma_ramp(1.02);
int n1 = 0, n2 = 0, n3 = 0;
for (int i = 0; i < N; i++){
    int c = w[i];
    if (c == 1) n1++;
    else if (c == 2) n2++;
    else n3++;
    float zx = float(n1) + float(n2) * u2Re + float(n3) * u3Re;
    float zy = float(n2) * u2Im + float(n3) * u3Im;
    float rx = zx * co - zy * si, ry = zx * si + zy * co;
    // canvas y runs down; negated so the tile leans as the plate shows it
    int pt = addpoint(0, set(rx, -ry, 0.0));
    setpointattrib(0, "Cd", pt, col);
    setpointattrib(0, "Alpha", pt, 0.75);
}

AFTER EFFECTS · DECLINED

Points, not a path: thousands of plotted prefix points of a substitution sequence; the line that visits them is a scribble, not the tile.