Ported the remaining implicit solver functions for Eigen.

Also added a couple of utility wrapper functions for Eigen types to make
interfacing with plain float arrays and blenlib math easier.

source/blender/physics/intern/BPH_mass_spring.cpp

@@ -519,7 +519,10 @@ static void cloth_calc_force(ClothModifierData *clmd, float UNUSED(frame), ListB
 		/* scale gravity force */
 		mul_v3_v3fl(gravity, clmd->scene->physics_settings.gravity, 0.001f * clmd->sim_parms->effector_weights->global_gravity);
 	}
-	BPH_mass_spring_force_gravity(data, gravity);
+	vert = cloth->verts;
+	for (i = 0; i < cloth->numverts; i++, vert++) {
+		BPH_mass_spring_force_gravity(data, i, vert->mass, gravity);
+	}
 #endif
 
 	cloth_calc_volume_force(clmd);

source/blender/physics/intern/implicit.h

@@ -59,6 +59,7 @@ extern "C" {
 //#define IMPLICIT_ENABLE_EIGEN_DEBUG
 
 struct Implicit_Data;
+struct ImplicitSolverInput;
 struct SimDebugData;
 
 typedef struct ImplicitSolverResult {
@@ -118,8 +119,6 @@ void BPH_mass_spring_set_velocity(struct Implicit_Data *data, int index, const f
 void BPH_mass_spring_get_motion_state(struct Implicit_Data *data, int index, float x[3], float v[3]);
 void BPH_mass_spring_set_vertex_mass(struct Implicit_Data *data, int index, float mass);
 
-int BPH_mass_spring_add_block(struct Implicit_Data *data, int v1, int v2);
-
 void BPH_mass_spring_clear_constraints(struct Implicit_Data *data);
 void BPH_mass_spring_add_constraint_ndof0(struct Implicit_Data *data, int index, const float dV[3]);
 void BPH_mass_spring_add_constraint_ndof1(struct Implicit_Data *data, int index, const float c1[3], const float c2[3], const float dV[3]);
@@ -131,9 +130,9 @@ void BPH_mass_spring_apply_result(struct Implicit_Data *data);
 /* Clear the force vector at the beginning of the time step */
 void BPH_mass_spring_clear_forces(struct Implicit_Data *data);
 /* Fictitious forces introduced by moving coordinate systems */
-void BPH_mass_spring_force_reference_frame(struct Implicit_Data *data, int index, const float acceleration[3], const float omega[3], const float domega_dt[3]);
+void BPH_mass_spring_force_reference_frame(struct Implicit_Data *data, int index, const float acceleration[3], const float omega[3], const float domega_dt[3], float mass);
 /* Simple uniform gravity force */
-void BPH_mass_spring_force_gravity(struct Implicit_Data *data, const float g[3]);
+void BPH_mass_spring_force_gravity(struct Implicit_Data *data, int index, float mass, const float g[3]);
 /* Global drag force (velocity damping) */
 void BPH_mass_spring_force_drag(struct Implicit_Data *data, float drag);
 /* Custom external force */

source/blender/physics/intern/implicit_blender.c

@@ -754,7 +754,8 @@ void BPH_mass_spring_solver_free(Implicit_Data *id)
 
 void BPH_mass_spring_solver_debug_data(Implicit_Data *id, struct SimDebugData *debug_data)
 {
-	id->debug_data = debug_data;
+	if (id)
+		id->debug_data = debug_data;
 }
 
 /* ==== Transformation from/to root reference frames ==== */
@@ -1215,7 +1216,7 @@ void BPH_mass_spring_set_vertex_mass(Implicit_Data *data, int index, float mass)
 	mul_m3_fl(data->M[index].m, mass);
 }
 
-int BPH_mass_spring_add_block(Implicit_Data *data, int v1, int v2)
+static int BPH_mass_spring_add_block(Implicit_Data *data, int v1, int v2)
 {
 	int s = data->M[0].vcount + data->num_blocks; /* index from array start */
 	BLI_assert(s < data->M[0].vcount + data->M[0].scount);
@@ -1294,7 +1295,7 @@ void BPH_mass_spring_clear_forces(Implicit_Data *data)
 	data->num_blocks = 0;
 }
 
-void BPH_mass_spring_force_reference_frame(Implicit_Data *data, int index, const float acceleration[3], const float omega[3], const float domega_dt[3])
+void BPH_mass_spring_force_reference_frame(Implicit_Data *data, int index, const float acceleration[3], const float omega[3], const float domega_dt[3], float mass)
 {
 #ifdef CLOTH_ROOT_FRAME
 	float acc[3], w[3], dwdt[3];
@@ -1316,7 +1317,7 @@ void BPH_mass_spring_force_reference_frame(Implicit_Data *data, int index, const
 	sub_v3_v3(f, coriolis);
 	sub_v3_v3(f, centrifugal);
 	
-	mul_m3_v3(data->M[index].m, f); /* F = m * a */
+	mul_v3_fl(f, mass); /* F = m * a */
 	
 	cross_v3_identity(deuler, dwdt);
 	cross_v3_identity(dcoriolis, w);
@@ -1326,11 +1327,11 @@ void BPH_mass_spring_force_reference_frame(Implicit_Data *data, int index, const
 	
 	add_m3_m3m3(dfdx, deuler, dcentrifugal);
 	negate_m3(dfdx);
-	mul_m3_m3m3(dfdx, data->M[index].m, dfdx);
+	mul_m3_fl(dfdx, mass);
 	
 	copy_m3_m3(dfdv, dcoriolis);
 	negate_m3(dfdv);
-	mul_m3_m3m3(dfdv, data->M[index].m, dfdv);
+	mul_m3_fl(dfdv, mass);
 	
 	add_v3_v3(data->F[index], f);
 	add_m3_m3m3(data->dFdX[index].m, data->dFdX[index].m, dfdx);
@@ -1344,19 +1345,14 @@ void BPH_mass_spring_force_reference_frame(Implicit_Data *data, int index, const
 #endif
 }
 
-void BPH_mass_spring_force_gravity(Implicit_Data *data, const float g[3])
+void BPH_mass_spring_force_gravity(Implicit_Data *data, int index, float mass, const float g[3])
 {
-	int i, numverts = data->M[0].vcount;
-	/* multiply F with mass matrix
-	 * force = mass * acceleration (in this case: gravity)
-	 */
-	for (i = 0; i < numverts; i++) {
-		float f[3];
-		world_to_root_v3(data, i, f, g);
-		mul_m3_v3(data->M[i].m, f);
-		
-		add_v3_v3(data->F[i], f);
-	}
+	/* force = mass * acceleration (in this case: gravity) */
+	float f[3];
+	world_to_root_v3(data, index, f, g);
+	mul_v3_fl(f, mass);
+	
+	add_v3_v3(data->F[index], f);
 }
 
 void BPH_mass_spring_force_drag(Implicit_Data *data, float drag)