Belousov–Zhabotinsky reaction

Reaction-diffusion model on cellular automata

Alasdair Turner (2009) A Simple Model of the Belousov-Zhabotinsky Reaction from First Principles
implementation note: discovery.ucl.ac.uk/id/eprint/17241/1/17241.pdf

en.wikipedia.org/wik…zhabotinsky_reaction

Log in to post a comment.

const W = 300;

const N1 = 0;
const N2 = 10;
const N3 = 20;
const N4 = 40;
const N5 = 60;
const N6 = 80;
const N7 = 100;
const N8 = 200;
const N9 = 400;

const N = [N1, N2, N3, N4, N5, N6, N7, N8, N9];

let done = 0;

let D = 60;


Canvas.setpenopacity(-0.75);

const t = new Turtle();

for (let x = 0; x < 3; x++) {
    for (let y = 0; y < 3; y++) {
        let cx = -100 + 200 / 3 * (x + 0.5);
        let cy = -100 + 200 / 3 * (y + 0.5);
        
        for (let r = D/2; r <= D/2 + 1; r += 0.1) {
            t.jump(cx - r, cy);
            t.seth(-90);
            t.circle(r);
        }
    }
}



let a = [], a_buf = [];
let b = [], b_buf = [];
let c = [], c_buf = [];
for (let x = 0; x < W; x++) {
    a[x] = []; a_buf[x] = [];
    b[x] = []; b_buf[x] = [];
    c[x] = []; c_buf[x] = [];
    for (let y = 0; y < W; y++) {
        a[x][y] = Math.random(); a_buf[x][y] = 0.0;
        b[x][y] = Math.random(); b_buf[x][y] = 0.0;
        c[x][y] = Math.random(); c_buf[x][y] = 0.0;
    }
}


function walk(frame) {
    if (frame == N[done]) {
        let cx = -100 + 200 / 3 * (done % 3 + 0.5);
        let cy = -100 + 200 / 3 * (Math.floor(done / 3) + 0.5);;
        
        for (let x = 0; x < W; x++) {
            for (let y = 0; y < W; y++) {
                if ((x - W/2)**2 + (y - W/2)**2 < (W/2)**2) {
                    let ac = a[x][y];
                
                    t.jump(cx + (x - W/2) * D / W, cy + (y - W/2) * D / W);
                    t.seth(360 * Math.random());
                    t.forward(ac * 0.2);
                }
            }
        }
        
        done++;
    }
    
    for (let x = 0; x < W; x++) {
        for (let y = 0; y < W; y++) {
            let an = 0.0, bn = 0.0, cn = 0.0;
            
            for (let dx = -1; dx <= +1; dx++) {
                for (let dy = -1; dy <= +1; dy++) {
                    an += a[(x + dx + W) % W][(y + dy + W) % W];
                    bn += b[(x + dx + W) % W][(y + dy + W) % W];
                    cn += c[(x + dx + W) % W][(y + dy + W) % W];
                }
            }
            
            an /= 9;
            bn /= 9;
            cn /= 9;
            
            a_buf[x][y] = an + an * (bn - cn);
            b_buf[x][y] = bn + bn * (cn - an);
            c_buf[x][y] = cn + cn * (an - bn);
        }
    }

    for (let x = 0; x < W; x++) {
        for (let y = 0; y < W; y++) {
            a[x][y] = Math.min(Math.max(a_buf[x][y], 0.0), 1.0);
            b[x][y] = Math.min(Math.max(b_buf[x][y], 0.0), 1.0);
            c[x][y] = Math.min(Math.max(c_buf[x][y], 0.0), 1.0);
        }
    }

    return done < N.length;
}