Even with a 20 point Laplacian operator there are still coordinate system artifacts with a circularly symmetric seed.
That is one reason to try a spectral solver.
The main code for the above mentioned simplest Laplacian with a forward Euler solver is:
#define A(U) texture(iChannel0,(U)/iResolution.xy)
void mainImage( out vec4 Q, in vec2 U )
{
// Lookup Field
Q = A(U);
// Mean Field
// Two way: horizontal, vertical
vec4 sum2 = A(U+vec2(0,1))+A(U+vec2(1,0))+A(U-vec2(0,1))+A(U-vec2(1,0));
vec4 mean2 = 1./4.*(sum2);
// Laplacian
vec4 laplacian2 = (mean2 - Q);
// Diffuse each variable differently :
Q += laplacian2 * vec4(1, .4, 1, 1);
// Compute reactions:
Q.x = Q.x * .99 + 0.01 * Q.y;
Q.y = Q.y + .05 * Q.y * (1. - Q.y) - .03 * Q.x - 1e-3;
// Prevent Negative Values (depends on system):
Q = max(Q, 0.);
}
How can this be rewritten to a spectral solver on Shadertoy?