Animated groups · reports · clockwork orbits

How to find a time-symmetric animation

How the looping spiral-lattice animation with clockwork symmetry r33₁3₂ was actually found: the five steps of the search that produced it, re-run for this page and reproducing the published file bit for bit. One equation throughout, the Brusselator.

The animation#

Three different things go into this film, and the page keeps them apart: a symmetry, which says how the picture repeats in space and time; an equation, which says how the chemicals change from one moment to the next; and a method, which finds a film that obeys both.

Its symmetry#

In orbifold notation the symmetry group is r33₁3₂. Each 3 is a kind of centre about which the picture can be turned by 120°: plain 3 is γ in the picture above, 3₁ is α and 3₂ is β. The subscript counts the thirds of a period that an anticlockwise turn about that centre costs. So turning the picture about γ changes nothing; about α it shows the film as it will be a third of a period later; about β, two thirds later. The prefix r says that sliding the picture by a third of a cell costs a third of a period as well.

Written out, with $q = (u, v)$ the two chemicals and each of the nine symmetries $g$ of the simulation box carrying its cost $\tau_g \in \{0, \tfrac13, \tfrac23\}$:

Only three of the nine cost nothing: doing nothing, and the two turns about γ. Those are all a single frame shows. The other six relate a frame to the one a third or two thirds of a period later.

The correspondence catalogue lists this group under ³3³3¹3: the same orbifold, with each centre marked by how many different moments its turns reach (³3 for α and β, ¹3 for γ). Its number there, 227, is only its place in Xu and Wu's list of the 275 spacetime groups.

The equation#

One equation, the Brusselator: two chemicals, $u$ and $v$, react and diffuse at the same rate.

Nothing in the equation knows about the symmetry. What has to be found: the film repeats with period $T$, so a single frame $q_0$ and the number $T$ determine all of it; everything else comes from simulating forward. The search is for one frame, 2 × 36 × 36 = 2 592 numbers, and one period. Nobody chooses $T$; it comes out as 6.390159918584649.

Why simulating does not work#

The obvious approach is to start somewhere, simulate, and wait. It fails for two reasons.

• Most of the symmetry is invisible to a simulator. It can impose the three free symmetries at every step, but the other six relate a frame to a later frame, so no single frame can satisfy or violate them.
• The animation is unstable. Even started exactly on it, a simulation drifts away: a disturbance of the right shape grows by a factor of 1.16 every period, and the simulation ends with every point oscillating in unison.

So instead of simulating towards the answer, the method writes the requirement as an equation for the frame and the period and solves it by Newton's method, which never asks whether a solution is stable. Newton's method needs a good first guess, and most of the work goes into that guess.

The method in five steps#

Step 1 · The rest state and how it oscillates

Where nothing varies in space, the chemicals can sit still at $u^* = a = 1$, $v^* = b/a = 2.3$, where production and consumption balance. Linearise the reaction there:

The positive real part means the rest state is unstable: any small disturbance starts to oscillate, and grows into a steady oscillation. (That happens once $b > 1 + a^2 = 2$; at $b = 2.3$ we are just past it.) The imaginary part gives the first guess for the period, $T_0 = 2\pi/0.9887 = 6.3551$, within 1 % of the final answer. The eigenvector, $\mathbf e = (1,\ −1.15 + 0.9887\,i)$, says how $v$ moves relative to $u$ in that oscillation; step 3 uses it.

Step 2 · A starting pattern from the Brusselator's slow envelope

Close to the onset of oscillation, every point of the Brusselator swings almost like a sine wave at the frequency of step 1. What varies across space and time is only how strongly each point swings, and in which phase. Put those two into one complex number per point, $A(x, t)$: the envelope of the oscillation. Near onset the Brusselator reduces to an equation for the envelope alone,

whose constants come from the Brusselator: $\alpha = 0$ because the two chemicals diffuse equally, and $\beta = 0.11$, close to the Brusselator's own value $1/9$ at $a = 1$. The job set $D = 2$ and measured the box in envelope units as $40\sqrt{2\mu} = 30.98$, with $\mu = b - 1 - a^2 = 0.3$ the distance past onset.

The envelope has one property the full Brusselator lacks: multiply a solution by a constant phase $e^{i\varphi}$ and you get another solution. So its simplest periodic solutions are rotating waves, $A(x,t) = B(x)\, e^{i\Omega t}$: a fixed pattern whose phase turns steadily. For such a solution, waiting a fraction $\tau$ of the period just multiplies the whole field by $e^{2\pi i \tau}$, so the symmetry of the film becomes

It is linear, so any field can be made to satisfy it by averaging over the nine symmetries, $(PB)(x) = \tfrac19 \sum_g e^{2\pi i \tau_g} B(g x)$. And because the envelope equation respects both the turns and the phase rotation, a field that satisfies the condition keeps satisfying it as it is simulated. The job did exactly that:

1. Start from random complex noise (random seed 1) and apply $P$.
2. Simulate the envelope 50 time units at a time, applying $P$ again after each stretch (only to remove round-off), and test whether the field has become a rotating wave, that is, whether $A_t = i\Omega A$ for a single number $\Omega$. The job allowed 600 time units.
3. Polish: solve $F(B) = i\Omega B$ to round-off with Newton's method, where $F$ is the right-hand side of the envelope equation.

t = 0 · noise, averaged by $P$

t = 2

t = 5

t = 8

t = 12

t = 50 · a rotating wave

The dark dots were predicted before any simulation. A turn about α leaves the point α where it is, so there the condition reads $B(\alpha) = e^{-2\pi i/3} B(\alpha)$, which forces $B(\alpha) = 0$; the same holds at β. They are the points α and β of the film at the top, where the chemicals never move.

By $t = 50$ the envelope is a rotating wave to 1.9e−8; the finite time step keeps it from getting closer than 3.4e−9, short of the job's test of 1e−9, so the run goes the full 600 time units (12 000 steps) and one polishing step finishes at 1.5e−14. The envelope is an approximation of the Brusselator, good near onset: the seed needs the right symmetry and roughly the right shape, and step 4 does the rest.

The sign. The minus in $e^{-2\pi i \tau_g}$ assumes the phase turns forwards. The job did not work that out in advance; it tried $e^{+2\pi i \tau_g}$ first. That also settles into a rotating wave, but with its clock running backwards: ⅓ and ⅔ swap, so α and β trade places. Checked against the catalogue entry it fails, with an error of 1.73 against a limit of 2e−7; the minus sign passes at 1.4e−15.

Step 3 · Turn it into a Brusselator frame

The envelope pattern becomes a frame for the two chemicals through the eigenvector of step 1:

Taking the real part gives a snapshot of the oscillation at one instant; multiplying by $\mathbf e$ makes $v$ lag $u$ the way the oscillation requires; the field is scaled to unit root-mean-square and then by the amplitude $A = 0.25$. (Had Newton failed from 0.25, the job would have tried 0.50; it did not need to.) The period guess is $T_0 = 6.3551$ from step 1. This seed has the right symmetry, inherited from $B$, but not the right shape or size: it misses the equation of step 4 by 6.9e−2.

Step 4 · Solve for one frame and one number

The equation. Take a candidate frame $q_0$ and period $T$. Simulate the Brusselator for a third of a period. Then apply a symmetry $g$ that costs a third of a period, say an anticlockwise turn about α. If $(q_0, T)$ is the answer, you are back where you started:

Here $\Phi_s$ means "simulate forward for time $s$" (RK4, time step 0.02, no symmetry imposed), and $S_g$ moves every value to where $g$ sends it. A third of a period is enough: doing it three times gives the whole period, and the ⅔ relations follow. Two more ingredients make the problem square:

• The three free symmetries are conditions on one frame, so they are imposed on $q_0$ by averaging. That leaves 868 free numbers.
• One more equation, a phase condition. Every frame of the film solves $F = 0$, so the solver must be told which frame to return: the total change from the starting guess must be perpendicular to the direction in which the film moves there.

Newton–Krylov. Newton's method linearises $F$ at the current guess and solves the resulting linear system for a correction. The matrix of that system has one column per unknown, and each column would cost a simulation, so it is never built. A Krylov solver (LGMRES) needs only the matrix times a vector, and that is one extra simulation: $\big(F(z + \varepsilon v) - F(z)\big)/\varepsilon$. That is all "Newton–Krylov" means. The job used SciPy's root(method='krylov'), with the period as $\log T$ so it stays positive.

Root-mean-square residual of $F$ against Newton iteration, on a log scale. Hover over the chart for each iteration's values, or open the table below.

For 35 iterations almost nothing happens: the residual goes from 6.9e−2 to 6.5e−2 and the period wanders between 6.355 and 6.363. The guess is outside the region where Newton's method converges quickly, and the steps are small and cautious. Then it gets close, and 8 iterations take the residual to 5.5e−16; in the last few, the number of correct digits roughly doubles each step, the signature of Newton's method. The whole solve took 43 iterations and 462 simulations of a third of a period, about 52 seconds on an ordinary laptop. The period it lands on is 6.390159918584649.

seed · from step 3

iteration 35 · end of the stall
6.5e−2

iteration 36
5.3e−2

iteration 37
2.2e−2

iteration 39
1.1e−3

iteration 43 · solved
5.5e−16

Notice what never comes up in this step: stability. $F = 0$ holds for an unstable film just as well as for a stable one, which is exactly why this works where simulating did not.

Step 5 · Make the film and check it

• The film. Simulate the Brusselator for one period from the solved frame, with no symmetry imposed, and save the frames (96 in the job; the published file keeps 48).
• The symmetry, checked. All nine relations $q(gx, t + \tau_g T) = q(x, t)$ hold on the film to 4.4e−9, against a limit of 2e−7. The film was simulated with nothing imposed, so this checks the answer rather than assuming it.
• Visible motion. A solution can satisfy the equation and still look static. The audit asks that every symmetry that costs time move the picture visibly at every moment, by at least 5 % of its range; this one moves it by at least 34 %.
• An independent replay. Simulate again from the first frame at two time steps without any averaging; the result must match the film to 1e−5. It matches to 9.3e−10 at the halved step, and closes up after one period to 4.2e−9. The published single-precision file replays to 4.7e−8, the round-off of single precision.
• A finer mesh. Resample to 48 × 48 nodes and solve again: 2 Newton steps, and the period moves by 6.2e−6 of itself, so the 36-node answer is resolved.

The published file is the solved frame polished once more, 48 frames of one period, averaged over the nine symmetries so that symmetric samples are identical, and stored in single precision. Re-running all five steps for this page reproduces it exactly: the same SHA-256 hash, e47c12f9099f…, and a maximum difference of 0.

Where the symmetry comes from#

Does the Brusselator make this symmetric pattern by itself, or does the mesh? Neither. The symmetry is put in on purpose, in two steps, and after that only checked.

The averaging in step 2 is what makes the pattern. Run the same noise through the same envelope equation without it, and the symmetric part, 11 % of the noise at the start, is down to 3.5e−6 of the field by $t = 50$; every point then oscillates in unison. The mesh does nothing beyond letting each symmetry carry mesh points onto mesh points, which is why its side, 36, is a multiple of 3.

Run it yourself#

The search ran on 7 September 2026 as one job of 204 in a batch, g227-bruss-a1b23-L40-s1: 53.75 seconds on a two-core cloud container, about $0.0023, no graphics processor. The script code/original.py runs the same job on your machine, one step at a time, with the repository's own functions, and records everything this page quotes, the control of the last section included:

git clone https://github.com/yaroslavvb/animated-groups-fable
cd animated-groups-fable/docs/reports/clockwork-orbits/code
python3 original.py   # ~2 minutes
python3 figures.py    # the pictures

original.py needs NumPy and SciPy and writes its record into original-run/; figures.py needs Pillow and redraws this page's picture panels from that record.

The functions it calls live in rdlab.py: cgl_relative_equilibrium and cgl_polish (step 2), seed_from_complex (step 3), rpo_polish (step 4) and audit (step 5). The job itself is run_brusselator in the batch driver, and the export is in build_equation_atlas.py. The code README explains what each output means.