Strange Attractors

25 Jun 2021 (updated 22 Jul 2025)

---

§ Overview

In this article, I explain how I implemented these visualizations. I'll defer to Wikipedia for a more detailed treatment on attractors, strange or not.

For these attractors, the motion at any point in space is defined by an ODE dependent only on the position at that point.

This means that if we visualize a strange attractor via particles, then there is no interaction between particles – i.e. embarrassingly parallel.

Using WebGL2, we can take advantage of the parallel nature of the system by simulating hundreds of thousands or millions of particles at once.

To do so, we can create two programs: one to compute particle positions, and one to render the particle to the framebuffer.

§ Halvorsen Attractor

One example of an attractor is Halvorsen's cyclically symmetric attractor.

The particle velocities of this attractor are given by the following ODEs,

$$\begin{align*} \dot{x}&=-\alpha x - 4y - 4z - y^2\\ \dot{y}&=-\alpha y - 4z - 4x - z^2\\ \dot{z}&=-\alpha z - 4x - 4y - x^2 \end{align*}$$

In the above demo, the constant $\alpha$ is set to $1.89.$

§ Compute Shader

In WebGL2 we are able to use transform feedback, this allows us to update vertex attributes in the vertex shader and have them persist in the next rendering pass.

In our case, we can use transform feedback on particle positions; we can treat each particle as a single vertex, then in the vertex shader apply the update rules according to the attractor, the result will be fed back into the shader in another pass.

#version 300 es
precision highp float;

uniform float u_Alpha;
uniform float u_Speed;
uniform sampler2D u_RgbNoise;

in vec3 i_Position;

out vec3 o_Position;

The uniform inputs:

• u_Alpha corresponds to the $\alpha$ parameter inside the Halvorsen equations.
• u_Speed is a scalar that is used for scaling the velocity vector.
• u_RgbNoise is an RGB32F texture with random values for each channel.

Then, the compute_velocity() function is a straight forward implementation of the Halvorsen equations.

vec3 compute_velocity()
{
    vec3 velo = vec3(0.);

    velo.x = -u_Alpha * i_Position.x
             - 4. * i_Position.y
             - 4. * i_Position.z
             - i_Position.y * i_Position.y;

    velo.y = -u_Alpha * i_Position.y
             - 4. * i_Position.z
             - 4. * i_Position.x
             - i_Position.z * i_Position.z;

    velo.z = -u_Alpha * i_Position.z
             - 4. * i_Position.x
             - 4. * i_Position.y
             - i_Position.x * i_Position.x;

    return velo;
}

Inside main(), we compute the new position via the following:

1. Compute position by adding the particle's velocity times u_Speed to its current position.
2. Sample the u_RgbNoise texture to obtain a random Vec3, and add it to position.
   • The long-term behaviour of many attractors is for the particles to converge into distinct lines, by adding some noise onto the position at each timestep we can delay or stop this from occurring.
3. Check if the particle's position exceeds a certain distance from the origin, if it does, teleport it back to a random position as sampled from u_RgbNoise.
   • Some initial conditions may cause particles to shoot out very far, out of the viewable area, so we bring it back with this step.
4. Check if the particle's life span exceeds a certain limit
   • See note on step 2.

In lieu of implementing a hash function on the GPU to use as a PRNG, we generate u_RgbNoise with the CPU and pass it to the GPU. Some less robust hash functions/PRNG methods on the GPU lead to obvious biases in the PRNG output, which can be noticable in some applications. Using the texture enables us to generate more "random randomness" with the built in Math.random() or crypto.getRandomValues().

void main()
{
    vec3 position = i_Position + compute_velocity() * u_Speed;

    // Implementation detail: u_RgbNoise is generated as a 512x512 texture
    ivec2 uv = ivec2(int(i_Position[0]) % 512, int(i_Position[1]) % 512);
    vec3 noise = texelFetch(u_RgbNoise, uv, 0).rgb;

    position = mix(position, noise, float(length(position) > 25.));
    position += noise;

    o_Position = position;
}

Because we do not want to draw anything to the screen while computing the particle positions, we can simply discard all the pixels in the frag shader:

#version 300 es
precision highp float;

void main()
{
    discard;
}