Exploring the Reaction-Diffusion Parameter Space — Four Faces of the Gray-Scott Model

#What Is the Gray-Scott Model?

Leopard spots, zebra stripes, coral branching structures. Nature is full of surprisingly regular patterns. Are these patterns written in genes as "draw a black dot here"? In 1952, Alan Turing presented a completely different hypothesis in his paper "The Chemical Basis of Morphogenesis": patterns emerge spontaneously from nothing more than the diffusion and reaction of two chemical substances.

The Gray-Scott model expresses this Turing theory as two partial differential equations. There are only two characters:

• u (activator): The "raw material" of the reaction. Supplied from outside at a constant rate
• v (inhibitor): The "product" of the reaction. Proliferates by reacting with u and is removed at a constant rate

The equations are:

∂u/∂t = Du·∇²u - u·v² + f(1 - u)
∂v/∂t = Dv·∇²v + u·v² - (f + k)·v

Breaking down each term:

• Du·∇²u / Dv·∇²v: Diffusion terms. u and v spread through space. Du is 1.0, Dv is 0.5 -- u diffuses faster. This asymmetry is the key
• u·v²: Reaction term. When u and v meet, u is consumed and v is produced. The v² means reaction accelerates where v is abundant (autocatalysis)
• f(1 - u): Supply term. f is the "feed rate." u is replenished more strongly where it's scarce
• (f + k)·v: Removal term. k is the "kill rate." v decays naturally

Just four lines of equations. Yet from these four lines, life-like patterns emerge through self-organization.

#Four Preset Patterns

The allure of the Gray-Scott model is that changing just two parameters, f and k, produces entirely different patterns. The ReactionDiffusion component includes four presets.

#Spots (f=0.055, k=0.062) -- Scattered Dots

f is relatively high and k is also high. Raw material supply is abundant, but product removal is fast. As a result, "islands" where v concentration is locally elevated form sporadically. Islands grow to a certain size and stabilize, spreading no further.

Intuitively, the balance of supply and removal is at an equilibrium that "maintains dots but allows nothing more." Think of leopard spots or giraffe patterns. Each spot exists independently, maintaining a fixed distance from others. This is because diffusion creates "territories" for the substance.

On screen, green-to-yellow dots scatter across a dark blue background. The pattern spreads like ripples from the initial seed (a circular patch at the center) until the entire grid is covered with spots.

#Stripes (f=0.035, k=0.065) -- Parallel Stripes

Lower f. When raw material supply is throttled, v can no longer survive as isolated "dots." Instead, adjacent v regions connect to form linear structures. Since k is high, line width stays thin, resulting in stripe patterns.

This transition is counterintuitive. Lowering f by just 0.02 causes the spatial pattern to jump from 0-dimensional (dots) to 1-dimensional (lines). Stripes initially point in random directions but gradually align and approach parallel over time. Think of a zebra's body.

On screen, when switching from the Spots preset, dots begin to melt and connect with each other. After several hundred frames, the screen fills with undulating stripes. The spacing between stripes is nearly uniform, determined by a characteristic wavelength dictated by the diffusion coefficient ratio (Du/Dv = 2.0).

#Coral (f=0.060, k=0.062) -- Coral-like Branching Structures

Raise f further. With abundant raw material supply, v regions grow actively. However, since k is the same as Spots at 0.062, the grown v gradually decomposes from the edges. The result is structures that branch while extending outward from the center. Fractal beauty reminiscent of coral or snowflake crystals.

Spots and Coral share the same k of 0.062. The only difference is f (0.055 vs 0.060). A difference of 0.005 separates "stable dots" from "ever-growing branches." This is the phase boundary of the parameter space, and the most fascinating region of Reaction-Diffusion.

On screen, ripples emanating from the seed branch out like fingers, tracing maze-like structures. Active reactions occur at branch tips, glowing yellow to red. The base is stable blue-green. This color gradient reflects temporal history.

#Mitosis (f=0.0367, k=0.0649) -- Cell Division-like Proliferation

The most biological preset. f is low and k is subtly low. In this region, when spots reach a stable size, they become unstable and split in two. The two split dots each grow and split again. It's mitosis itself.

Why do they split? f is low, so growth is slow. But k is also low, so v persists longer. When a dot grows large enough, u at the center becomes depleted (supply can't keep up). The v concentration at the center drops, the dot "pinches," and divides in two.

The screen animation is hypnotic. Slowly swelling dots suddenly become two. The two dots slowly swell again and split again. They increase exponentially until the grid is filled with uniformly sized dots.

#Implementation Details

ReactionDiffusion.tsx is 296 lines. When faithfully translating equations to code, you can clearly see where the line count goes.

#Grid and Double Buffering

const GRID = 200
const DU = 1.0
const DV = 0.5
const DT = 1.0

function createGrid(): Float32Array {
  return new Float32Array(GRID * GRID)
}

200x200 = 40,000 cells. Four Float32Arrays are used: u, v, and their respective next buffers (for writing computation results). Double buffering is a standard technique for separating "read source" and "write destination." If computed on a single buffer, updates to upper-left cells would affect calculations for lower-right cells, creating update-order-dependent distortion.

Six simulation steps are run per frame:

const STEPS_PER_FRAME = 6

for (let s = 0; s < STEPS_PER_FRAME; s++) {
  simulate(u, v, nu, nv, f, k)
  // Swap buffers
  const tmpU = u
  const tmpV = v
  u = nu
  v = nv
  nu = tmpU
  nv = tmpV
}

STEPS_PER_FRAME = 6 is not a stability requirement from the CFL condition (Courant-Friedrichs-Lewy condition) but a visual speed adjustment. With DT = 1.0 and Du = 1.0, the CFL condition for the 2D diffusion equation is dt <= dx^2 / (4*Du), which gives dt <= 0.25 as the stability limit with grid spacing dx = 1. In practice, however, DT = 1.0 remains stable. This is because the Gray-Scott reaction terms work to counteract diffusion, making the stability limit more relaxed than for the pure diffusion equation.

#5-Point Laplacian Stencil

Discrete Laplacian approximation of nabla^2 u. Uses the center cell and its four orthogonal neighbors:

const lapU = u[y * n + xm] + u[y * n + xp]
           + u[ym * n + x] + u[yp * n + x]
           - 4 * uC

A 9-point stencil (including diagonal cells) would improve accuracy but more than doubles computation. 200x200 x 6 steps x 60 FPS = 144 million additions/subtractions per second. If the 5-point stencil provides sufficient accuracy, there's no need for extra computation.

#Periodic Boundary Conditions

const ym = y === 0 ? n - 1 : y - 1
const yp = y === n - 1 ? 0 : y + 1
const xm = x === 0 ? n - 1 : x - 1
const xp = x === n - 1 ? 0 : x + 1

The grid's top and bottom, left and right edges are connected. Mathematically, the simulation runs on a torus (the surface of a donut). This avoids boundary artifacts while keeping patterns seamlessly continuous.

#Color LUT

A lookup table that precomputes the conversion from v concentration to color:

const COLOR_LUT = new Uint8Array(256 * 4)
// dark blue -> cyan -> green -> yellow -> red

Five color stops (dark blue, cyan, green, yellow, red) are linearly interpolated to create a 256-level x RGBA table of 1024 bytes. Since 40,000 pixels need color computation every frame, computing trigonometric functions or exponents in real-time is out of the question. A table lookup takes just one memory access.

When converting v values to 0-255 indices, the value is multiplied by 512 and bitwise-OR'd for integer conversion:

const ci = Math.min(255, Math.max(0, (vVal * 512) | 0))

With a 512x multiplier, v = 0.5 gives an index of 256 (clamped to 255). In actual patterns, v concentration falls within the 0-0.5 range, so the full color space is utilized.

#Seeding with Clicks

The most intuitive experience with the ReactionDiffusion component is clicking (or touching) the canvas to plant pattern "seeds."

const SEED_RADIUS = 5

function seedAt(u: Float32Array, v: Float32Array, gx: number, gy: number) {
  for (let dy = -SEED_RADIUS; dy <= SEED_RADIUS; dy++) {
    for (let dx = -SEED_RADIUS; dx <= SEED_RADIUS; dx++) {
      if (dx * dx + dy * dy > SEED_RADIUS * SEED_RADIUS) continue
      const x = gx + dx
      const y = gy + dy
      if (x >= 0 && x < GRID && y >= 0 && y < GRID) {
        const idx = y * GRID + x
        u[idx] = 0.5 + Math.random() * 0.02
        v[idx] = 0.25 + Math.random() * 0.02
      }
    }
  }
}

A circular patch with radius 5 cells is injected at the pointer position. Adding small random fluctuations to u=0.5 and v=0.25 breaks symmetry, causing different pattern growth from each injection point.

The essence of this feature is "the sensation of strolling through the 'phase diagram' of parameter space." Select f/k with a preset, click to generate new nuclei, and watch patterns spread. Click with Spots selected and dots spread; with Stripes, stripes race across the screen. The fact that the same action yields completely different results captures the essence of Reaction-Diffusion.

Pointer events are handled with onPointerDown/Move/Up, and the conversion from canvas coordinates to grid coordinates maintains the aspect ratio:

const dim = Math.min(rect.width, rect.height)
const ox = (rect.width - dim) / 2
const oy = (rect.height - dim) / 2
pointerPosRef.current = {
  x: ((e.clientX - rect.left - ox) / dim) * GRID,
  y: ((e.clientY - rect.top - oy) / dim) * GRID,
}

The square simulation area is centered on the canvas, with coordinate conversion accounting for margins.

#Parameter Sensitivity

The most astonishing property of Reaction-Diffusion is that tiny parameter changes produce qualitatively different patterns.

Compare Spots and Coral. k is the same at 0.062. Only f differs: 0.055 vs 0.060. This difference of 0.005 is a very narrow distance in parameter space. Yet the results are entirely different topologies: "stable dots" vs "growing branches."

Raise Stripes' f=0.035 slightly to f=0.036, and some stripes begin to break apart into dots. At f=0.037, stripes and dots coexist. At f=0.038, dots dominate. Transitions happen in increments of just 0.001.

This sensitivity is fundamentally different from chaos. The chaos of a double pendulum is the phenomenon where "tiny differences in initial conditions grow exponentially over time." The parameter sensitivity of Reaction-Diffusion is the phenomenon where "tiny differences in parameters select different attractors." The former is unpredictability; the latter is bifurcation.

What's important is that starting from the same parameters and initial conditions always reproduces the same pattern. Deterministic and reproducible. Yet shifting f by just 0.001 qualitatively changes the reproduced pattern. This property of being "reproducible yet sensitive" is what makes Reaction-Diffusion mathematically fascinating.

Viewing the f/k two-dimensional parameter space as a map, there exist "Spots Country," "Stripes Country," "Coral Country," and "Mitosis Country," with surprisingly narrow borders between them. The four preset points are stable representative points, like "capitals" of each country. Choosing parameters near the borders can produce chimera states where multiple patterns coexist.

#Turing's Prescience

In his 1952 paper "The Chemical Basis of Morphogenesis," Turing claimed that "animal patterns form spontaneously through diffusion and reaction of chemical substances." This claim was prescient for three reasons.

First, the double helix structure of DNA hadn't even been discovered yet (Watson and Crick's paper came in 1953). In an era when no one understood how genes control morphology, he declared "patterns can arise from chemical reactions alone."

Second, Turing's theory presupposed numerical verification. Analytically solving the partial differential equations is difficult, requiring computer simulation. But computers in 1952 (the Manchester Mark I) couldn't dream of real-time simulation on a 200x200 grid.

Third, Turing's theory was experimentally confirmed more than half a century after the paper. In 2006, Turing patterns were confirmed in Roe deer hair follicle cells, and in 2012, wrinkle formation via reaction-diffusion mechanisms was demonstrated in mouse palates.

In 2026, we can run a 200x200 grid at 60 frames per second in a browser. From Float32Array memory allocation on the Canvas API to color conversion, everything is done in JavaScript. The patterns Turing could only imagine in his head can be generated with a single click, explored while adjusting parameters. This experience itself is a testament to 74 years of technological progress.

#Summary: Lessons from Reaction-Diffusion

The biggest insight from implementing Reaction-Diffusion is that complexity doesn't necessarily arise from complex rules.

The equations are four lines. The parameters are essentially two (f and k). The implementation is 296 lines. What accounts for the difference between equations and implementation? Boundary condition handling, double buffering, color LUT, pointer interaction, canvas rendering via requestAnimationFrame. In other words, the "transformation layer for showing computation to humans" accounts for over 90% of the total.

This ratio is suggestive. Natural laws are concise, but the apparatus to make them observable is complex. Browser-based simulation has the same structure.

Another lesson: exploring parameter space took more time than writing the code. The values f=0.055, k=0.062 were found through dozens of trial-and-error runs as points where "beautiful patterns emerge." Implementing the equations takes an hour. Finding beautiful parameters took half a day. The boundary between mathematics and art lies in parameter selection.

With just two parameters, life-like patterns emerge. Witnessing this fact on screen is an experience that can never be gained from reading equations alone.