Skip to the plate
FORMA PUBLIC DOMAIN GENERATIVE ATLAS / ED. 0.28
Plate 153, Bowl of Integers: a still of the lattices plate as the atlas renders it, in the lattices accent.

PL. 153  ·  LATTICES / LATTICES

Bowl of Integers

Frederick Soddy, 1937

OPEN THE LIVE PLATE ▸

DEFINITION

(b₁+b₂+b₃+b₄+b₅)² = 3(b₁²+b₂²+b₃²+b₄²+b₅²)
b₅′ = b₁+b₂+b₃+b₄ − b₅, and the same on (b·x, b·y, b·z)
root (−1, 2, 2, 3, 3); every curvature is ≡ 0 or 2 (mod 3)
solid drawn: ⋃ { |x − cᵢ| ≤ rᵢ − ρ } over rᵢ > ρ, cut by x·â ≤ s
ray: |o + t·û − cᵢ| = rᵢ − ρ,  n̂ = (x − cᵢ)/(rᵢ − ρ)

NOTES

Soddy gave this packing its name in Nature in January 1937, the year after the poem: the bowl of integers, a unit sphere filled with spheres whose curvatures are all whole numbers, and the hexlet, a ring of six that closes. The construction is the one PL. 152 draws in section and this plate draws whole. Write a sphere as the four numbers curvature, and curvature times each coordinate of the centre, and the fifth sphere tangent to four of a mutually tangent five is the sum of all five minus twice the one it replaces, applied to all four numbers at once. Start from the curvatures minus one, two, two, three and three — a unit sphere holding two of radius one half stacked on an axis and two of radius one third on the waist — and every sphere the recurrence makes has a whole-number curvature forever, because in space the reflection carries a coefficient of one where in the plane it carries two. The six spheres of curvature three are the hexlet: a closed ring at sixty degree spacing round the waist, each of radius one third, each centred two thirds of the way out, each touching both of the halves and its two neighbours and the bowl. That is what the plate opens on, cut through its own middle, and it was counted rather than assumed — exactly six spheres of curvature three exist in the whole packing, at thirty, ninety and a hundred and fifty degrees either side, and nowhere else. What has to be said before anything else about the picture is what was done to it to make a picture possible at all. A sphere packing in three dimensions leaves nothing over: the spheres fill the bowl completely, the part they miss has volume zero, and a lit render of the object as it stands is a smooth ball with a faint crust and nothing to see inside. So every sphere here is drawn shrunk by one number the reader controls, and any sphere smaller than that number is not drawn. That is an exact operation and not an approximation: for a set of balls that touch at most at single points, shrinking the union is the same as shrinking each ball, so what is drawn is precisely the balls of radius r minus rho with the ones at or below rho gone. It is the only way to see into a packing, the gaps between the spheres are the width of it, and the plate would be dishonest if it let a reader think the real packing has room in it. The geometry is exact everywhere else. The spheres are enumerated from the recurrence, so there is no distance estimate here and nothing is approximated: a ray meets a sphere where a quadratic says it does, and the normal at that point is the radius direction, which is why the surfaces come out clean where an estimated field would leave them fused. The rays are cast through a grid of cells that holds the spheres each cell meets, so a ray tests a few dozen of them instead of all of them, and the grid was checked against testing every sphere on every ray rather than trusted. There is no shader on this plate, and the reason is worth stating because every other lit solid in this wave carries one. A fragment shader receives numbers, one per constant, and a list of eleven hundred spheres has no channel into it. The way round that is well known and was tried here in full: the Apollonian group is generated by inversions in the five spheres orthogonal to four of the roots at a time, so a point can be folded into a fundamental domain and the distance read off the root spheres divided by the scale the folding accumulated, which needs no list. Those five mirrors were derived, and they reproduce the curvature-centre reflection exactly. The estimate they give is a heuristic and not a bound, and measured against the enumerated packing it is not good enough: the surfaces it finds sit measurably in front of the real ones and the shrunk spheres come out fused into lumps rather than standing apart. The numbers are under Provenance. A wrong picture that renders cleanly is this atlas worst failure, so the shader is refused by name rather than shipped with a caveat. Where the mathematics stops and the picture begins: the packing, the shrinking, the cut plane, the intersections and the normals are the plate. The light is not — two directions fixed relative to the eye with one slider turning the key round, and where a lamp stands is staging. Lambert on the normal, a halfway-vector specular and a hemisphere of short rays for the ambient term are published models rather than invented ones, which is why they are named under Provenance. The optics of the page run after all of it and outside this function entirely. The camera does not move and it cannot: casting this many rays is far more than one frame holds, so the picture develops instead — every paint traces a batch into a buffer, coarsest first, the first paint laying down a whole eighth-resolution picture and later ones refining it. A camera move would throw all of that away, so the view is constants the reader turns rather than an animation, and turning one restarts a development that takes about a twentieth of a second.

PROVENANCE

Origin
F. Soddy, "The Bowl of Integers and the Hexlet", Nature 139 (9 January 1937), pages 77 to 79, doi:10.1038/139077a0 — re-resolved against the live Crossref record for this plate: title, sole author Frederick Soddy, container Nature, volume 139, issue 3506, pages 77 to 79, dated 1937-01-09, all matching. It is the paper that names both the integral packing and the ring of six, and Graham, Lagarias, Mallows, Wilks and Yan call it perhaps the first observation of integrality in these packings. The five-sphere relation itself is the closing stanza of the same author "The Kiss Precise", Nature 137 (20 June 1936), 1021, doi:10.1038/1371021a0, which is what PL. 152 cites; the n-dimensional form is Nature 139 (1937), 62, doi:10.1038/139062a0, attributed to Thorold Gosset throughout the literature although the Crossref record for that note carries no author field at all — the caveat PL. 152 records, repeated here rather than dropped.
The one it is not
PL. 152 draws a plane travelling through this same packing and strokes the circles it cuts, which are almost never tangent to one another. This plate draws the spheres themselves as a lit solid and cuts them open, so the flat faces in the picture ARE that section, filled in. PL. 36 is the two-dimensional gasket, four circles rather than five spheres, with a reflection of coefficient two.
Curvature-centre coordinates
R. L. Graham, J. C. Lagarias, C. L. Mallows, A. R. Wilks and C. H. Yan, "Apollonian circle packings: geometry and group theory III. Higher dimensions", Discrete and Computational Geometry 35 (2006), pages 37 to 72, doi:10.1007/s00454-005-1197-8, re-resolved for this plate against the live record: all five authors, volume 35 issue 1, pages 37 to 72, published online 17 October 2005 and printed in the January 2006 volume. The atlas cites the printed volume, which is the dating caveat PL. 152 records. The four numbers per sphere are Definition 3.2 and the reflection is equation 4.1. Theorem 4.1 of the same paper is why three dimensions is the last one: the group is discrete for n = 2 and n = 3 and not discrete above, so in four dimensions the orbit spheres overlap and it stops being a packing.
The hexlet, counted
Not asserted, counted. The recurrence was run to exhaustion below a radius of 0.0018 — 67,219 spheres — and every sphere of curvature exactly 3 was listed: there are six, all at distance two thirds from the centre, all of radius one third, all in the plane through the centre perpendicular to the axis of the two curvature-2 spheres, at azimuths minus 150, minus 90, minus 30, 30, 90 and 150 degrees. Sixty degree spacing, exactly, and none anywhere else in the packing. The published quintuple relation was checked on the root at the same time — the square of the sum is 81 and three times the sum of the squares is 81, exactly, in integers — and so was the residue claim PL. 152 makes: over the whole enumeration every curvature is congruent to 0 or 2 modulo 3 and none is ever 1.
Two independent enumerations agree
The packing is built here by the curvature-centre recurrence. It was checked against a completely different construction that shares no line of reasoning with it: apply the five inversions of the Apollonian group to the five root spheres AS SPHERES, breadth first, and keep everything reachable in at most N applications. At N = 4 that gives 299 spheres of positive curvature; every one of them is present in the curvature-centre enumeration, worst match error 5.8e-16 over centre and radius together, and the same six curvature-3 members come out at the same six azimuths. The counts by word length are 5, 10, 30, 90, 300, 990, 3330 — one new sphere per generator at the first step, because a generator moves one sphere of a mutually tangent five and fixes the other four.
Why there is no shader, measured rather than assumed
Every other lit solid in this wave carries a GLSL body; this one does not, and the route that would have given it one was built and measured before being refused. A fragment shader gets one float per constant and cannot be handed a sphere list, so the shader-shaped way to draw this object is the iterated inversion of K. Nakamura, "Iterated inversion system: an algorithm for efficiently visualizing Kleinian groups and extending the possibilities of fractal art", Journal of Mathematics and the Arts 15(2), 2021, pages 106 to 136, doi:10.1080/17513472.2021.1943998 (resolved against the live Crossref record for this plate): fold a point into a fundamental domain by repeated inversion, accumulate the conformal scale factor, and divide the distance to the base configuration by it. The five mirrors were derived from the root quintuple rather than taken from anywhere — the sphere orthogonal to four of the five roots, found by solving |c − cⱼ|² = R² + rⱼ² for the other four, which for two of the five is inconsistent and says the dual is a PLANE through the four centres concerned. What comes out is three spheres, centred at (√3/6, 0, 0) with radius √3/6 and at (2√3/3, 0, ∓1) with radius 2√3/3, and two planes through the origin with normals (1/2, ±√3/2, 0). Each of the five inversions reproduces the curvature-centre reflection on all five roots to within 1.1e-16, and the duals are orthogonal to their four roots to 0.0 exactly. One structural fact fell out and is worth recording because it is not the two-dimensional case: in space the group is NOT a free product of five copies of the two-element group. Applied to the root quintuple, the product of any two generators has order three — (S₀S₁)³ returns (−1, 2, 2, 3, 3) exactly, checked in integers — so the five mirrors meet at sixty degrees and the fundamental domain is a truncated simplex in four-dimensional hyperbolic space whose five vertices are all hyperideal, the truncating spheres being the five roots. The fold terminates: over 20,000 uniform points the median point folds 5 times, the 99th percentile 17, and 3 in 20,000 reach a cap of 64. And the estimate it gives fails. Marched at the same views as the picture that ships, against a brute-force ray/sphere trace of the enumerated packing at the same shrink radius — 131,500 rays over five parameter tuples — it agrees on hit against miss for 96.47 per cent of rays, but 28.95 per cent of the rays that hit stop more than 0.01 of the bounding radius IN FRONT of the real surface, mean absolute error 3.14e-2 and worst 1.29. That is not a silhouette wobble, it is the shrunk spheres fusing into lumps: the estimate keeps residue from spheres far below the shrink radius, and adding the shrink back turns the residue into material. Capping the fold by the accumulated scale changes it by nothing measurable (96.47 per cent and 28.95 per cent, the same to three figures). Without the shrink the estimate is much better — over 288 exterior points sampled at controlled clearances it never once exceeded the true distance — but then the object is a filled ball with no picture in it. So the shader is refused, the JS path is the only path, and the numbers are here rather than a sentence saying it was hard.
The rendering, and every part of it named
Rays are cast, not marched, because the geometry is exact: a ray meets a sphere where |o + t·û − c| = r has a real root, which is a quadratic, and the normal at the hit is the radius direction. No distance estimate and no epsilon, so nothing here inherits sphere tracing and Hart 1996 is deliberately NOT cited — it is the paper the refused route would have rested on. Shading: Lambert, Photometria, 1760, pre-DOI and public domain, of which the whole content here is max(0, n dot l); the specular is the halfway vector of J. F. Blinn, "Models of light reflection for computer synthesized pictures", Proceedings of SIGGRAPH 77, pages 192 to 198, doi:10.1145/563858.563893 (Crossref carries the title, the author and the pages and leaves volume and issue empty, which is why it is cited as the proceedings rather than in the Computer Graphics 11(2) form). The ambient term is S. Zhukov, A. Iones and G. Kronin, "An ambient light illumination model", Eurographics Rendering Techniques 98, Springer 1998, pages 45 to 55, doi:10.1007/978-3-7091-6453-2_5: W(P) is the mean over the hemisphere of a monotone function of the distance to the nearest occluder, zero at zero and one at infinity, and what runs here is four cosine-weighted directions on a golden-angle set, each cast to a radius R with the simplest member of the family they give. It is deliberately not the five-tap distance-field trick that circulates as ambient occlusion in demo code, which has no published derivation. The shadow is a single ray to the key and is hard rather than soft, because with exact geometry there is no unbounding sphere to take a cone from. All three identifiers were resolved against the live Crossref record for this plate.
The shrink, and the acceleration, both checked
Shrinking is exact for this object and that is a property of packings rather than a convenience: the balls meet at most at single points, so the erosion of their union is the union of their erosions -- a ball of radius rho placed anywhere near a tangency sticks out into the interstice and so lies in neither ball, hence in neither erosion nor the erosion of the union. What is drawn is therefore precisely the balls of radius r minus rho with everything at or below rho absent, and nothing about that step is approximate. The grid that carries the rays is checked rather than trusted, and the check is the strong one: 160,000 rays over four parameter tuples spanning both ends of the shrink dial, both named view axes, three lenses and two spins, each traced BOTH through the cell walk and against every sphere in the packing with no acceleration at all. Zero disagreements -- same sphere, same curvature, and the hit distances equal to the last bit, worst absolute difference exactly 0.0. At the default constants a ray tests about 19 spheres of the 1,182 present and visits about 12 cells. Sphere counts and build cost by the shrink radius, measured through the atlas own test harness in node rather than in a browser: 129 spheres at 0.08 building in 1.4 ms, 293 at 0.05, 1,182 at 0.03, 2,956 at 0.02, 6,247 at 0.015 in 5.6 ms. The build is spread across paints at 2,600 quintuples each rather than run in one, which is the Gray-Scott lesson: a plate should not lurch when it scrolls into view.
The colour, and why it stays in the bright lobe
The order palette is brightest near t = 0 and t = 1 and essentially black at t = 0.5, and a lit solid may look like the FIELD case that is allowed to cross the trough. It is not, for the reason PL. 147 records: the material that would land near black is the material a reader would read as a hole, and on a packing a black region is what a missing sphere looks like. What ships is a rising monotone segment of the bright lobe, 0.72 to 0.97, mapped from the logarithm of the curvature of the sphere the ray landed on — an integer, known exactly, because the sphere is enumerated rather than estimated. At the default shrink radius the two spheres of curvature 2 sit at the dark end of that segment, the six of curvature 3 one step up from them, and the smallest drawn at the bright end, so the hexlet reads as a band and not as an accident of the light.
The framing, which is definitional rather than measured
Every other raymarched plate in this wave derives its camera distance from a silhouette radius it has to compute. This is the one plate whose bounding radius is not a measurement: the packing is contained in the sphere of curvature minus one, which is the unit sphere exactly, by construction. The camera sits at 1.12 of it on the frame half-height, which leaves the object filling about nine tenths of the height with air around it. The view axis runs from the hexlet own axis, where the ring of six is exactly six-fold and countable, to the plane of the ring seen edge on. The cut plane is perpendicular to the view axis, PL. 147 arrangement, and it is INTERSECTED with the solid rather than applied to the ray, so the flat faces are real surfaces of a real object and the shadow and the ambient term see what the eye sees.
Standing
Public domain — a theorem of 1936 and a paper of 1937, with nineteenth-century antecedents. No patent ever applied to a sphere packing, and none applies to any of the three rendering models named above.
What one paint can afford
One paint cannot afford a frame of this, so a paint is a batch, and the batch IS the eighth-resolution pass rather than a number that happens to land near it -- so the first paint is a whole coarse picture by construction. The grid is capped at 420 cells across, which is 318 by 256 at a card (the card is narrower than the cap, so it renders one cell per CSS pixel) and 420 by 280 at the drawer. Measured through the atlas own test harness in node at the default constants: a card is 81,408 rays in 64 batches, warm median 0.59 ms and 95th percentile 0.92, 57 ms in total, then 0.01 ms a frame for ever because a developed picture is a blit; the drawer is 117,600 rays in 64 batches at a median of 0.70 ms and 51 ms in total; a 240-square thumbnail is 57,600 rays in 65 batches and 32 ms; the 2048-wide export is 117,600 rays at a median of 0.68 ms. The first paint is dearer than the rest and that is the packing being built rather than the rays being cast: 1.0 to 1.7 ms depending on the shrink radius, spread over one to three paints. The batch count is the constraint that decides the cap, because renderStill develops an exposure over 241 draws and a cap whose full pass needed more than that would ship a half-traced thumbnail and a half-traced export. 64 to 67 batches leaves it a factor of 3.6. There is no shader frame to report beside the sibling solids, for the reason two blocks above.
Constants, measured
Nothing is locked, and that had to be measured rather than assumed, because grain looks exactly like a cost dial. It is not one. The sphere count runs over a factor of forty across its declared range -- 6,247 at 0.015 against 129 at 0.08 -- and the batch cost barely notices, because a ray tests what the cell it is in holds and not what the packing holds: warm median 0.58 ms a batch at 0.015, 0.56 at 0.02, 0.68 at 0.03, 0.54 at 0.045, 0.45 at 0.06 and 0.39 at 0.08, all at a card frame. What it does price is the build, 1.4 ms to 5.6 ms, and that is spread across paints. Nor is the dial dead at either end: it is the difference between a hundred and twenty-nine spheres and six thousand, and every value in the range renders the same packing at a different cutoff rather than a different figure. A dial that is neither dear nor dead is one regenerate can have -- PL. 147 own re-measured reasoning, applied here from nothing. The other five are the camera and the light and all five are free. A camera whose distance is fixed by a bounding radius that is exactly one cannot be jittered into looking at nothing, and the key light rides a cone at 42 degrees about the eye so that lightaz turns it round the viewer rather than away from the object: a surface facing the camera is never darker than the cosine of 42 degrees whatever the slider says, which is what makes that constant safe to jitter where an azimuth measured in world coordinates would not be. That was found by measuring rather than by taste -- the first rig had a fixed elevation with the azimuth turning only the horizontal part, and at 90 degrees off axis it took the mean of the frame from 66 to 19 in 0..255. Swept as regenerate actually jitters, all six at once, 24 seeds at two sizes over 140 frames: 0 throws, 0 plates below the contrast floor, pixelRange worst 196 and median 212 against the gate floor of 14. Swept over the WHOLE declared box instead, 40 uniform tuples at two sizes: still 0 and 0, with a worst of 201. Every figure in this entry was taken in node through the atlas own test harness, which is not a browser and does not pretend to be one.
Source
doi:10.1038/139077a0

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. 153 · BOWL OF INTEGERS — Frederick Soddy, 1937
//   (b₁+b₂+b₃+b₄+b₅)² = 3(b₁²+b₂²+b₃²+b₄²+b₅²)
//   b₅′ = b₁+b₂+b₃+b₄ − b₅, and the same on (b·x, b·y, b·z)
//   root (−1, 2, 2, 3, 3); every curvature is ≡ 0 or 2 (mod 3)
//   solid drawn: ⋃ { |x − cᵢ| ≤ rᵢ − ρ } over rᵢ > ρ, cut by x·â ≤ s
//   ray: |o + t·û − cᵢ| = rᵢ − ρ,  n̂ = (x − cᵢ)/(rᵢ − ρ)
// 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=bowl

float p_grain  = 0.03 + chf('grain_tweak');     // shrink radius — smallest sphere drawn · live 0.015 .. 0.08
float p_orient = 0.18 + chf('orient_tweak');    // view axis — hexlet axis to its plane · live 0 .. 1
float p_cut    = 0.5 + chf('cut_tweak');        // section cut along the view axis · live 0 .. 0.6
float p_spin   = 0 + chf('spin_tweak');         // view spin about the hexlet axis ° · live -180 .. 180

// The plate's own colour: FORMA's LATTICES 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.4185 + 0.4549 * cos(6.28318530718 * (t + 0.05)),
    0.1389 + 0.151 * cos(6.28318530718 * (t + 0.1)));
}

// The bowl of integers itself, as geometry: one point per sphere of the
// packing, carrying the radius the plate draws it at and the colour the plate
// gives it. Copy a sphere to the points and the object is there; nothing about
// the picture is baked in.
//
// The recurrence is the one PL. 152 ports one plate over, and the reasoning is
// the same: the Soddy-Gosset relation (b1+..+b5)^2 = 3(b1^2+..+b5^2) read as a
// quadratic in one curvature has two roots that sum to the other four, so the
// OTHER sphere tangent to a given four is b1+b2+b3+b4-b5 with no root taken and
// no sign to choose. In curvature-centre coordinates (b, b*x, b*y, b*z) --
// Graham, Lagarias, Mallows, Wilks and Yan, part III, Definition 3.2 -- the
// identical reflection applies to all four numbers at once (their equation
// 4.1), and the centres come back by dividing by b. The root quintuple is
// -1, 2, 2, 3, 3 and 9^2 = 3 * 27 exactly, so every curvature below it is a
// whole number. VEX has no recursion, so the reflection tree is an explicit
// queue of quintuples walked with a read cursor -- breadth first, for the
// reason soddy.vex records: depth first plus a duplicate test builds a smaller
// packing, because a sphere first reached down a deep branch keeps that deep
// level and the short branch is then thrown away as a duplicate.
//
// TWO THINGS THE PLATE DOES THAT THIS PORT DOES DIFFERENTLY, both on the
// record. First, the shrink. A sphere packing in three dimensions leaves
// nothing over, so the plate draws every sphere shrunk by p_grain and drops
// anything at or below it -- exact for a set of balls that touch at single
// points, and the only way to see into a packing. pscale carries the shrunk
// radius, so a copied sphere is the plate's sphere. The unshrunk radius is on
// the point as forma_r should a reader want the packing as it really is, and
// the integer curvature is on it as forma_b, which is the point of the plate.
// Second, the cut. The plate intersects the solid with a half-space and the
// flat faces in its picture are that section; a point cloud cannot carry a
// half-cut sphere. So this port DROPS the spheres the plane removes outright
// and leaves the ones it grazes whole, and writes the plane itself onto the
// detail as forma_cutn and forma_cutd -- one Boolean node away from the
// plate's own figure.
//
// The bounding sphere is not emitted. Its curvature is -1, so it is the frame
// rather than a member: the plate does not draw it either, and a copied sphere
// of radius 1 would swallow the packing.
//
// waived: lens, lightaz — the camera throw and the key light. Both are pure
// staging on the page (the plate says so in its own note: where a lamp stands
// is not a fact about the packing), and a Houdini scene owns its own camera
// and its own lights. Every constant that decides which spheres exist and
// where they are -- grain, orient, spin, cut -- reaches the body below.

float TAU = 6.28318530718;
int   forma_cap = 24000;      // the plate's own ceiling; no setting in range
                              // reaches it (6,247 at the finest grain)

// The packing, as parallel arrays. The curvature is an int on purpose: it
// comes out an exact whole number at every step, because the recurrence is a
// sum and difference of integers and never divides, and that is what makes it
// a safe bucket key below.
int   cbs[];
float bxs[], bys[], bzs[];
float rt3 = sqrt(3.0);
push(cbs, -1); push(bxs, 0.0); push(bys,  0.0); push(bzs,  0.0);
push(cbs,  2); push(bxs, 0.0); push(bys,  0.0); push(bzs,  1.0);
push(cbs,  2); push(bxs, 0.0); push(bys,  0.0); push(bzs, -1.0);
push(cbs,  3); push(bxs, rt3); push(bys,  1.0); push(bzs,  0.0);
push(cbs,  3); push(bxs, rt3); push(bys, -1.0); push(bzs,  0.0);

// Spheres bucketed by curvature, as head-of-chain plus next-link, so the
// duplicate test costs a short walk instead of a scan of the whole list. A
// sphere is only ever kept when 1/b is greater than p_grain, so b is below
// 1/p_grain and the table sized from the grain cannot be overrun.
int slots = int(ceil(1.0 / p_grain)) + 2;
int bhead[], bnext[];
for (int i = 0; i < slots; i++){
    push(bhead, -1);
}
for (int i = 0; i < 5; i++){
    int k = cbs[i] > 0 ? cbs[i] : 0;   // the bounding sphere parks in slot 0
    push(bnext, bhead[k]);
    bhead[k] = i;
}

int q0[] = {0};
int q1[] = {1};
int q2[] = {2};
int q3[] = {3};
int q4[] = {4};
int head = 0;

while (head < len(q0) && len(cbs) < forma_cap){
    int i0 = q0[head], i1 = q1[head], i2 = q2[head], i3 = q3[head], i4 = q4[head];
    head += 1;
    int quint[] = array(i0, i1, i2, i3, i4);
    int sb = 0;
    float sx = 0.0, sy = 0.0, sz = 0.0;
    for (int j = 0; j < 5; j++){
        int m = quint[j];
        sb += cbs[m];  sx += bxs[m];  sy += bys[m];  sz += bzs[m];
    }
    for (int j = 0; j < 5; j++){
        int m = quint[j];
        int nb = sb - 2 * cbs[m];
        // Curvature 0 is a plane and negative curvature is an enclosing
        // sphere: neither is a child inside a bounded packing, and the
        // bounding sphere is already in the list. These are the degeneracies
        // of the quadratic, refused here by name. The grain floor is what
        // terminates the walk -- every child has a larger curvature, so the
        // radius falls under it eventually and nothing more is pushed.
        if (nb <= 0) continue;
        if (1.0 / float(nb) <= p_grain) continue;
        float nx = sx - 2.0 * bxs[m];
        float ny = sy - 2.0 * bys[m];
        float nz = sz - 2.0 * bzs[m];
        // Duplicate test. Two distinct spheres of curvature b have centres at
        // least 2/b apart, so in these coordinates they differ by at least 2 --
        // far outside the tolerance below, and far outside anything the
        // arithmetic could drift.
        int dup = 0;
        for (int k = bhead[nb]; k >= 0 && dup == 0; k = bnext[k]){
            if (abs(bxs[k] - nx) + abs(bys[k] - ny) + abs(bzs[k] - nz) < 0.5) dup = 1;
        }
        if (dup == 1) continue;
        int ni = len(cbs);
        push(cbs, nb);  push(bxs, nx);  push(bys, ny);  push(bzs, nz);
        push(bnext, bhead[nb]);
        bhead[nb] = ni;
        push(q0, j == 0 ? ni : i0);
        push(q1, j == 1 ? ni : i1);
        push(q2, j == 2 ? ni : i2);
        push(q3, j == 3 ? ni : i3);
        push(q4, j == 4 ? ni : i4);
    }
}

// The cutting plane, built the way the plate builds its view axis: a rotation
// of the z axis about y through orient * 90 degrees -- the arc from the
// hexlet's own axis to the plane the ring of six lies in -- and then about z
// by spin. cut = 0 takes nothing, 0.5 puts the plane through the centre. The
// bounding radius is exactly 1 by construction, which is why the offset is
// quoted in it directly with nothing measured.
float th = p_orient * 90.0 * M_PI / 180.0;
float sp = p_spin * M_PI / 180.0;
float ct = cos(th), st = sin(th), cs = cos(sp), sn = sin(sp);
vector cutn = set(st * cs, st * sn, ct);
float  cutd = 1.0 - 2.0 * p_cut;
setdetailattrib(0, "forma_cutn", cutn, "set");
setdetailattrib(0, "forma_cutd", cutd, "set");

// The colour is the curvature of the sphere -- an integer, known exactly --
// on a RISING segment of the palette's bright lobe, 0.72 to 0.97. Crossing
// the dark middle of the ramp would put material near black, and on a packing
// black is what a missing sphere looks like, so the picture would state its
// own subject backwards. The two spheres of curvature 2 sit at the dark end,
// the six of curvature 3 one step up, and the smallest drawn at the bright
// end -- so the hexlet reads as a band.
float kspan = log(0.5 / p_grain);

for (int i = 1; i < len(cbs); i++){       // 1: the container is the frame
    float b = float(cbs[i]);
    float r = 1.0 / b;
    vector c = set(bxs[i], bys[i], bzs[i]) / b;
    float rs = r - p_grain;
    if (rs <= 0.0) continue;              // shrunk away
    // Dropped only when the plane removes the whole sphere; a sphere the
    // plane grazes is emitted whole and the detail attributes above say
    // where to cut it.
    if (dot(c, cutn) - r > cutd) continue;
    float band = clamp(log(max(2.0, b) / 2.0) / kspan, 0.0, 1.0);
    int pt = addpoint(0, c);
    setpointattrib(0, "pscale", pt, rs);
    setpointattrib(0, "Cd", pt, forma_ramp(0.72 + 0.25 * band));
    setpointattrib(0, "forma_b", pt, b);
    setpointattrib(0, "forma_r", pt, r);
}