FDGL: Use inverse cubic charge law and tide force to prevent graph explosion.

git-svn-id: http://svn.drobilla.net/lad/trunk/ganv@5190 a436a847-0d15-0410-975c-d299462d15a1

2 files changed, 56 insertions, 17 deletions

diff --git a/src/Canvas.cpp b/src/Canvas.cpp
index 20180e3..65f2b16 100644
--- a/src/Canvas.cpp
+++ b/src/Canvas.cpp
@@ -734,7 +734,7 @@ GanvCanvasImpl::layout_iteration()
 	static const double SPRING_K = 14.0;
 
 	// A light directional force to push sources to the top left
-	static const double DIR_MAGNITUDE = -300.0;
+	static const double DIR_MAGNITUDE = -600.0;
 	Vector              dir           = { 0.0, 0.0 };
 	switch (_gcanvas->direction) {
 	case GANV_DIRECTION_RIGHT: dir.x = DIR_MAGNITUDE; break;
@@ -795,6 +795,15 @@ GanvCanvasImpl::layout_iteration()
 			continue;
 		}
 
+		/* Add tide force which pulls all objects as if the layout is happening
+		   on a flowing river surface.  This prevents disconnected components
+		   from being ejected, since at some point the tide force will be
+		   greater than distant repelling charges. */
+		const Vector mouth = { -10000.0, -10000.0 };
+		node->impl->force = vec_add(
+			node->impl->force,
+			tide_force(mouth, reg.pos, 4000000000000.0));
+
 		FOREACH_ITEM(_items, j) {
 			if (i == j || (!GANV_IS_MODULE(*i) && !GANV_IS_CIRCLE(*i))) {
 				continue;
@@ -818,8 +827,8 @@ GanvCanvasImpl::layout_iteration()
 
 		GanvNode* const node = GANV_NODE(*i);
 
-		static const float dur  = 0.15;  // Time duration
-		static const float damp = 0.7;  // Velocity damping (momentum loss)
+		static const float dur  = 0.1; // Time duration
+		static const float damp = 0.5; // Velocity damping
 
 		const bool has_edges = (node->impl->has_in_edges ||
 		                        node->impl->has_out_edges);
@@ -827,10 +836,19 @@ GanvCanvasImpl::layout_iteration()
 			node->impl->vel.x = 0.0;
 			node->impl->vel.y = 0.0;
 		} else {
+			node->impl->vel = vec_mult(node->impl->vel, damp);
 			node->impl->vel = vec_add(node->impl->vel,
 			                          vec_mult(node->impl->force, dur));
-			node->impl->vel = vec_mult(node->impl->vel, damp);
 
+			// Clamp velocity
+			const double        vel_mag = vec_mag(node->impl->vel);
+			static const double MAX_VEL = 4000.0;
+			if (vel_mag > MAX_VEL) {
+				node->impl->vel = vec_mult(
+					vec_mult(node->impl->vel, 1.0 / vel_mag),
+					MAX_VEL);
+			}
+				                           
 			// Update position
 			const Vector dpos = vec_mult(node->impl->vel, dur);
 			if (fabs(dpos.x) >= 1.0 || fabs(dpos.y) >= 1.0) {
diff --git a/src/fdgl.hpp b/src/fdgl.hpp
index 42b3538..084ebed 100644
--- a/src/fdgl.hpp
+++ b/src/fdgl.hpp
@@ -16,9 +16,8 @@
 #include <float.h>
 #include <math.h>
 
-static const double SPRING_K    = 16.0;
-static const double CHARGE_KE   = 40000.0;
-static const double AREA_WEIGHT = 0.5;
+static const double SPRING_K  = 10.0;
+static const double CHARGE_KE = 50000000.0;
 
 struct Region {
 	Vector pos;
@@ -52,6 +51,13 @@ vec_mult(const Vector& a, const Vector& b)
 	return a.x * b.x + a.y * b.y;
 }
 
+/** Magnitude. */
+inline double
+vec_mag(const Vector& vec)
+{
+	return sqrt(vec.x * vec.x + vec.y * vec.y);
+}
+
 /** Reciprocal of magnitude. */
 inline double
 vec_rmag(const Vector& vec)
@@ -64,9 +70,9 @@ inline Vector
 spring_force(const Vector& a, const Vector& b, double length, double k)
 {
 	const Vector vec          = vec_sub(b, a);
-	const double rmag         = vec_rmag(vec);
-	const double displacement = length - (1.0 / rmag);
-	return vec_mult(vec, rmag * k * displacement * 0.5);
+	const double mag          = vec_mag(vec);
+	const double displacement = length - mag;
+	return vec_mult(vec, k * displacement * 0.5 / mag);
 }
 
 /** Spring force with a directional force to align with flow direction. */
@@ -80,14 +86,29 @@ edge_force(const Vector& dir,
 	return vec_add(dir, spring_force(hpos, tpos, length, k));
 }
 
-/** Modified Coulomb's law */
+/** Constant tide force, does not vary with distance. */
+inline Vector
+tide_force(const Vector& a, const Vector& b, double power)
+{
+	static const double G   = 0.0000000000667;
+	const Vector        vec = vec_sub(a, b);
+	const double        mag = vec_mag(vec);
+	return vec_mult(vec, G * power / mag);
+}
+
+/**
+   Repelling charge force.
+
+   Many FDGL algorithms use charge according to Coulomb's law, but here we use
+   an inverse cube (not squared) law so influence drops off more rapidly with
+   distance.  This, in conjunction with a tide, keeps the layout compact.
+*/
 inline Vector
 repel_force(const Region& a, const Region& b)
 {
-	const Vector vec      = vec_mult(vec_sub(a.pos, b.pos), 4.0);
-	const double rmag     = vec_rmag(vec);
-	const Vector a_weight = vec_mult(a.area, AREA_WEIGHT);
-	const Vector b_weight = vec_mult(b.area, AREA_WEIGHT);
-	const Vector force    = vec_mult(vec, rmag * rmag * rmag * CHARGE_KE * 0.5);
-	return vec_mult(force, vec_mult(a_weight, b_weight));
+	const Vector vec   = vec_sub(a.pos, b.pos);
+	const double mag   = vec_mag(vec);
+	const Vector force = vec_mult(
+		vec, (CHARGE_KE * 0.5 / (mag * mag * mag * mag * mag)));
+	return vec_mult(force, vec_mult(a.area, b.area));
 }