Devlog #94 – Wave 73 WebGL/GLSL Special: Mandelbulb, Menger Sponge, Wormhole 3D, Supershape & GPGPU Particles

GPU Technologies Used in This Wave

New Simulations

Mandelbulb: Raymarching a 3D Fractal

The Mandelbulb is the most successful 3D analogue of the Mandelbrot set. It uses "triplex" algebra — a way of multiplying 3D vectors using spherical coordinate exponentiation — to generalise the 2D iteration z → z^n + c into three dimensions.

Triplex Power and Distance Estimator

The iteration at power n maps a point (r, θ, φ) in spherical coordinates to:

r_new    = r^n
theta_new = n * theta
phi_new   = n * phi
x_new    = r_new * sin(theta_new) * cos(phi_new) + cx
y_new    = r_new * cos(theta_new) + cy
z_new    = r_new * sin(theta_new) * sin(phi_new) + cz

The distance estimator uses the derivative magnitude dr, accumulated alongside r during iteration: dr = n·r^(n-1)·dr + 1. The lower bound on distance to the set is de = 0.5 · |log(r)| · r / dr. This safe step size lets the raymarcher advance quickly through empty space and slow down near the surface — sphere tracing rather than fixed-step marching.

Shading: Normals, Shadows, and Ambient Occlusion

Normals are estimated by finite differences: evaluate the SDE at six nearby points (±ε on each axis) and take the gradient. Soft shadows are approximated by marching a secondary ray toward the light and tracking the minimum ratio of step size to distance along the way. Ambient occlusion samples five distances along the surface normal and computes how much the surface curves back on itself — concave crevices are darker.

• Power n slider from 2 (classic Mandelbulb) to 12 (increasingly spiky)
• 96 raymarching steps per pixel; early-out when distance < 0.001
• 4 colour palettes: Inferno (fiery), Viridis (scientific), Mono (clean), Iris (purple-cyan)
• Mouse drag to orbit; scroll to zoom into surface details

Menger Sponge: Iterated SDF Folding

The Menger sponge is defined by recursive subtraction: take a cube, divide it into 3×3×3 = 27 sub-cubes, remove the 7 central cross pieces (face centres and the body centre), and repeat on each remaining sub-cube. After infinite iterations, the result has Hausdorff dimension log(20)/log(3) ≈ 2.7268 — less than 3D but more than 2D.

GLSL SDF via Iterated Folding

The signed distance function for the sponge is computed entirely inside the fragment shader using an elegant iterative fold technique:

for (int i = 0; i < iterations; i++) {
    p = abs(p);                    // fold to positive octant
    if (p.x < p.y) p.xy = p.yx;  // sort axes
    if (p.x < p.z) p.xz = p.zx;
    if (p.y < p.z) p.yz = p.zy;
    p = p * 3.0 - 2.0;            // scale and offset
    if (p.z < -1.0) p.z += 2.0;  // fold back
}
float d = (length(max(abs(p) - 1.0, 0.0))) / pow(3.0, float(iterations));

Each outer loop doubles the effective iteration depth. The resulting SDF can be sphere-traced at interactive frame rates for up to 6 iterations of the sponge.

Three Rendering Modes and Mathematical Properties

• Solid: standard raymarching with Phong shading and ambient occlusion — reveals the intricate surface texture
• Wireframe: render edges of the SDF isosurface slightly offset — exposes the cubic self-similarity
• X-ray: accumulate density along each ray proportional to near-miss distance — shows internal structure like a CT scan

The mathematical properties panel displays live values: total surface area (6·(5/9)^n·20^n → 0 as n → ∞ — the surface shrinks at each step), and the Hausdorff dimension formula log(N_retained)/log(scale_factor) = log(20)/log(3).

Wormhole 3D: Spacetime Embedding Shader

A Morris-Thorne traversable wormhole is a hypothetical solution to Einstein's field equations that connects two distant regions of spacetime via a "throat." It cannot be created without exotic matter (negative energy density), but it makes a spectacular visualisation and a genuine lesson in differential geometry.

Flamm's Paraboloid Embedding

The spatial geometry of a wormhole at a moment in time is described by the Flamm paraboloid. In the equatorial plane, the proper radial distance ℓ maps to Euclidean embedding radius r(ℓ) = √(b₀² + ℓ²) and height z(ℓ) = b₀ · arcsinh(ℓ/b₀), where b₀ is the throat radius. This gives two flared funnels connected at r = b₀.

GLSL Shader: Two-Universe Sky and Throat

The fragment shader determines which side of the throat a ray ends up on by integrating the geodesic path through the embedding. The "sky" of each universe is a discrete set of coloured stars — the two universes have different colour temperatures. As b₀ decreases (tighter throat), more of each universe's sky is visible through the aperture. The Doppler-like blueshift gradient near the throat gives a visual sense of the gravitational depth. Three modes:

• Embedding: third-person view of the full paraboloid surface rotating in space
• FPV (first-person): fly through the throat — watch one universe give way to the other
• Cross-section: slice through the equatorial plane showing the two connected sheets

• b₀ (throat radius) slider: narrow pinhole to wide tunnel
• Throat pulse animation: periodic oscillation in b₀
• Star density and colour temperature independent for each universe

Supershape: Gielis Superformula Geometry

Johan Gielis' 2003 superformula is a single parametric equation that generates an extraordinary variety of shapes — from platonic solids to flowers to alien geometries — by varying six parameters. Applied as a spherical product (once for latitude θ, once for longitude φ), it creates 3D surfaces with remarkable richness.

The Superformula and Spherical Product

The polar radius in 2D is:

r(φ) = (|cos(m·φ/4) / a|^n2 + |sin(m·φ/4) / b|^n3)^(-1/n1)

For the 3D spherical product, evaluate r₁(θ) and r₂(φ) independently. The 3D surface point is:

x = r1(θ) · cos(θ) · r2(φ) · cos(φ)
y = r1(θ) · sin(θ) · r2(φ) · cos(φ)
z = r2(φ) · sin(φ)

The mesh is a 300×300 grid of UV points, each evaluated independently and uploaded to Three.js BufferGeometry. Parameter changes rebuild the geometry in a single JavaScript pass at ~20ms.

Eight Presets and Material System

• Sphere: m=0, n1=n2=n3=1, a=b=1 — the baseline
• Star: m=5, n1=1, n2=n3=2 — pentagonal spiky form
• Twisted Donut: asymmetric n parameters with non-integer m — self-intersecting torus variant
• Crown: high m, low n1 — serrated equatorial band
• Pillow: low m, high n — inflated square cross-section
• Diamond: m=4, all n equal — octahedral approximation
• Flower: m=7, low n1, high n2/n3 — petal radiating form
• Mushroom: asymmetric θ/φ parameters — cap and stem form

• Matte and metallic PBR material modes
• Wireframe overlay toggle
• Save PNG button captures the current view at full canvas resolution
• All six superformula parameters exposed as live sliders

GPGPU Particles: 250,000 Particles on the GPU

General-Purpose GPU computing (GPGPU) repurposes the graphics pipeline as a massively parallel computing engine. Each "pixel" of a floating-point render texture stores the state of one particle. The fragment shader that "draws" this texture actually computes the new particle state — and by reading from one texture while writing to another (ping-pong), we can simulate hundreds of thousands of particles at interactive frame rates entirely on the GPU.

Ping-Pong Architecture

The simulation maintains two pairs of RGBA FloatType WebGLRenderTargets: one for positions (x, y, z, age) and one for velocities (vx, vy, vz, mass). Each frame:

1. Bind "compute" shader with current position/velocity textures as uniforms
2. Render a full-screen quad into the next position/velocity textures
3. Swap current and next textures (ping ↔ pong)
4. Render particles as GL_POINTS, sampling position from the current texture in the vertex shader

No data ever leaves the GPU between frames — the CPU only sets initial conditions and tweaks uniforms.

Four Particle Behaviours

• Curl noise: velocity field derived from the curl of a 3D Simplex noise gradient — divergence-free, so particles never clump or disperse; produces smooth swirling tendrils
• N-body attractors: 4 fixed attractor points with F = m/r²·n̂ plus quadratic drag; particles orbit in chaotic figure-eights
• Boids approximation: separation, alignment, and cohesion rules encoded as distance-weighted functions of the local density field sampled from the position texture
• Lorenz attractor: each particle independently integrates the Lorenz ODE (ẋ=σ(y-x), ẏ=rx-y-xz, ż=xy-bz) from a slightly different initial condition — the butterfly attractor rendered in particle density

• N slider: 1,000 to 250,000 particles (512×512 texture)
• WebGL2 capability check: graceful fallback message on unsupported devices
• Point size and opacity respond to particle age for born/death fading
• FPS counter displays GPU compute cost across N levels

Engineering Notes: WebGL2 and Floating Point Textures

Getting 250,000 particles working correctly required solving several WebGL2-specific challenges:

Float Texture Precision

WebGL2 supports OES_texture_float (32-bit RGBA float textures) as a core feature, removing the extension dance required in WebGL1. However, linear filtering of float textures still requires OES_texture_float_linear — the GPGPU particle sim uses nearest-neighbour sampling since each particle maps to exactly one texel.

MRT and Simultaneous Position/Velocity Update

Updating position and velocity simultaneously without using the previous frame's velocity would require reading and writing the same texture — impossible in WebGL. The solution is separate position and velocity textures, each in their own ping-pong pair. The compute shader reads both from the "current" set and writes both to the "next" set, keeping the full state consistent.

Particle Indexing

Each particle's UV coordinate in the texture is computed from its index: uv = vec2(float(index % WIDTH), float(index / WIDTH)) / float(WIDTH). The vertex shader reads this UV from the position texture and uses the result as gl_Position. With 250,000 particles in a 500×500 texture, every texel is a particle.

Up Next

That's Wave 73 — all 5 simulations, all GPU-native. With 638 simulations now live across all categories, we're pausing to consolidate: fixing mobile layouts, adding missing keyboard shortcuts, and expanding the search index before starting Wave 74. Check the blog archive for all previous devlogs.