1 files changed, 37 insertions, 11 deletions

diff --git a/src/shader/compute.glsl b/src/shader/compute.glsl
index e322a08..d8714b6 100644
--- a/src/shader/compute.glsl
+++ b/src/shader/compute.glsl
@@ -6,8 +6,33 @@ layout (std430, binding=1) buffer bufferA{ float data[]; };
 
 uniform vec2 world;
 
-float rand(vec2 co){
-	return fract(sin(dot(co.xy, vec2(12.9898,78.233))) * 43758.5453);
+// vector field definition
+
+vec2 f(vec2 v) {
+	return vec2(
+		cos(v.x*sin(v.y)),
+		sin(v.x-v.y)
+	);
+}
+
+// ODE solver
+
+vec2 explicitEuler(float h, vec2 v) {
+	return v + h * f(v);
+}
+
+// pseudo random numbers for particle placement
+
+float rand(vec2 v){
+	return fract(sin(dot(v, vec2(12.9898,78.233))) * 43758.5453);
+}
+
+float mapUnitToWorldX(float s) {
+	return -(world.x/2.) + s * world.x;
+}
+
+float mapUnitToWorldY(float s) {
+	return -(world.y/2.) + s * world.y;
 }
 
 bool insideWorld(vec2 v) {
@@ -18,17 +43,18 @@ bool insideWorld(vec2 v) {
 }
 
 void main() {
-	uint idx = 3*gl_GlobalInvocationID.x;
-	vec2 v = vec2(data[idx+0], data[idx+1]);
+	const uint i = 3*gl_GlobalInvocationID.x;
+	const vec2 v = vec2(data[i+0], data[i+1]);
+	const vec2 w = explicitEuler(0.01, v);
 
-	if ( data[idx+2] < 5. && insideWorld(v) ) {
-		data[idx+0] += 0.01 * cos(v.x*sin(v.y));
-		data[idx+1] += 0.01 * sin(v.x-v.y);
-		data[idx+2] += 0.01;
+	if ( data[i+2] < 5. && insideWorld(v) ) {
+		data[i+0]  = w.x;
+		data[i+1]  = w.y;
+		data[i+2] += 0.01;
 	} else {
-		data[idx+0] = -(world.x/2.) + rand(vec2(data[idx+1], data[idx+0])) * world.x;
-		data[idx+1] = -(world.y/2.) + rand(vec2(data[idx+0], data[idx+1])) * world.y;
-		data[idx+2] = uint(rand(v) * 5.);
+		data[i+0] = mapUnitToWorldX(rand(v));
+		data[i+1] = mapUnitToWorldY(rand(w));
+		data[i+2] = rand(v+w) * 5.;
 	}
 }
 )";