[Bf-blender-cvs] [e57d3af] gooseberry: Ported the remaining implicit solver functions for Eigen.
Lukas Tönne
noreply at git.blender.org
Thu Oct 30 15:59:39 CET 2014
Commit: e57d3af847d6f618ab4e03d80b341c7eb9da0db6
Author: Lukas Tönne
Date: Sun Oct 12 16:14:15 2014 +0200
Branches: gooseberry
https://developer.blender.org/rBe57d3af847d6f618ab4e03d80b341c7eb9da0db6
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.
===================================================================
M source/blender/physics/intern/BPH_mass_spring.cpp
M source/blender/physics/intern/implicit.h
M source/blender/physics/intern/implicit_blender.c
M source/blender/physics/intern/implicit_eigen.cpp
===================================================================
diff --git a/source/blender/physics/intern/BPH_mass_spring.cpp b/source/blender/physics/intern/BPH_mass_spring.cpp
index ceb76af..875b7bf 100644
--- a/source/blender/physics/intern/BPH_mass_spring.cpp
+++ b/source/blender/physics/intern/BPH_mass_spring.cpp
@@ -620,7 +620,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);
diff --git a/source/blender/physics/intern/implicit.h b/source/blender/physics/intern/implicit.h
index 4659bb0..287c064 100644
--- a/source/blender/physics/intern/implicit.h
+++ b/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 */
diff --git a/source/blender/physics/intern/implicit_blender.c b/source/blender/physics/intern/implicit_blender.c
index 741aeb0..7421a7a 100644
--- a/source/blender/physics/intern/implicit_blender.c
+++ b/source/blender/physics/intern/implicit_blender.c
@@ -925,7 +925,8 @@ int BPH_mass_spring_solver_numvert(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 ==== */
@@ -1386,7 +1387,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);
@@ -1465,7 +1466,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];
@@ -1487,7 +1488,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);
@@ -1497,11 +1498,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);
@@ -1515,19 +1516,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)
diff --git a/source/blender/physics/intern/implicit_eigen.cpp b/source/blender/physics/intern/implicit_eigen.cpp
index 402ffcb..a6148b6 100644
--- a/source/blender/physics/intern/implicit_eigen.cpp
+++ b/source/blender/physics/intern/implicit_eigen.cpp
@@ -111,7 +111,7 @@ public:
coeffRef(k) = v[k];
}
- fVector &operator = (const ctype &v)
+ fVector& operator = (const ctype &v)
{
for (int k = 0; k < 3; ++k)
coeffRef(k) = v[k];
@@ -142,7 +142,7 @@ public:
coeffRef(l, k) = v[k][l];
}
- fMatrix &operator = (const ctype &v)
+ fMatrix& operator = (const ctype &v)
{
for (int k = 0; k < 3; ++k)
for (int l = 0; l < 3; ++l)
@@ -156,27 +156,63 @@ public:
}
};
-typedef Eigen::VectorXf lVector;
+/* Extension of dense Eigen vectors,
+ * providing 3-float block access for blenlib math functions
+ */
+class lVector : public Eigen::VectorXf {
+public:
+ typedef Eigen::VectorXf base_t;
+
+ lVector()
+ {
+ }
+
+ template <typename T>
+ lVector& operator = (T rhs)
+ {
+ base_t::operator=(rhs);
+ return *this;
+ }
+
+ float* v3(int vertex)
+ {
+ return &coeffRef(3 * vertex);
+ }
+
+ const float* v3(int vertex) const
+ {
+ return &coeffRef(3 * vertex);
+ }
+};
typedef Eigen::Triplet<Scalar> Triplet;
typedef std::vector<Triplet> TripletList;
typedef Eigen::SparseMatrix<Scalar> lMatrix;
+/* Constructor type that provides more convenient handling of Eigen triplets
+ * for efficient construction of sparse 3x3 block matrices.
+ * This should be used for building lMatrix instead of writing to such lMatrix directly (which is very inefficient).
+ * After all elements have been defined using the set() method, the actual matrix can be filled using construct().
+ */
struct lMatrixCtor {
- lMatrixCtor(int numverts) :
- m_numverts(numverts)
+ lMatrixCtor()
+ {
+ }
+
+ void reset()
+ {
+ m_trips.clear();
+ }
+
+ void reserve(int numverts)
{
/* reserve for diagonal entries */
m_trips.reserve(numverts * 9);
}
- int numverts() const { return m_numverts; }
-
- void set(int i, int j, const fMatrix &m)
+ void add(int i, int j, const fMatrix &m)
{
- BLI_assert(i >= 0 && i < m_numverts);
- BLI_assert(j >= 0 && j < m_numverts);
i *= 3;
j *= 3;
for (int k = 0; k < 3; ++k)
@@ -184,15 +220,22 @@ struct lMatrixCtor {
m_trips.push_back(Triplet(i + k, j + l, m.coeff(l, k)));
}
- inline lMatrix construct() const
+ void sub(int i, int j, const fMatrix &m)
+ {
+ i *= 3;
+ j *= 3;
+ for (int k = 0; k < 3; ++k)
+ for (int l = 0; l < 3; ++l)
+ m_trips.push_back(Triplet(i + k, j + l, -m.coeff(l, k)));
+ }
+
+ inline void construct(lMatrix &m)
{
- lMatrix m(m_numverts, m_numverts);
m.setFromTriplets(m_trips.begin(), m_trips.end());
- return m;
+ m_trips.clear();
}
private:
- const int m_numverts;
TripletList m_trips;
};
@@ -247,6 +290,7 @@ BLI_INLINE const float *lVector_v3(const lVector &v, int vertex)
return v.data() + 3 * vertex;
}
+#if 0
BLI_INLINE void triplets_m3(TripletList &tlist, float m[3][3], int i, int j)
{
i *= 3;
@@ -289,6 +333,7 @@ BLI_INLINE void lMatrix_sub_triplets(lMatrix &r, const TripletList &tlist)
t.setFromTriplets(tlist.begin(), tlist.end());
r -= t;
}
+#endif
BLI_INLINE void outerproduct(float r[3][3], const float a[3], const float b[3])
{
@@ -297,11 +342,50 @@ BLI_INLINE void outerproduct(float r[3][3], const float a[3], const float b[3])
mul_v3_v3fl(r[2], a, b[2]);
}
+BLI_INLINE void cross_m3_v3m3(float r[3][3], const float v[3], float m[3][3])
+{
+ cross_v3_v3v3(r[0], v, m[0]);
+ cross_v3_v3v3(r[1], v, m[1]);
+ cross_v3_v3v3(r[2], v, m[2]);
+}
+
+BLI_INLINE void cross_v3_identity(float r[3][3], const float v[3])
+{
+ r[0][0] = 0.0f; r[1][0] = v[2]; r[2][0] = -v[1];
+ r[0][1] = -v[2]; r[1][1] = 0.0f; r[2][1] = v[0];
+ r[0][2] = v[1]; r[1][2] = -v[0]; r[2][2] = 0.0f;
+}
+
+BLI_INLINE void madd_m3_m3fl(float r[3][3], float m[3][3], float f)
+{
+ r[0][0] += m[0][0] * f;
+ r[0][1] += m[0][1] * f;
+ r[0][2] += m[0][2] * f;
+ r[1][0] += m[1][0] * f;
+ r[1][1] += m[1][1] * f;
+ r[1][2] += m[1][2] * f;
+ r[2][0] += m[2][0] * f;
+ r[2][1] += m[2][1] * f;
+ r[2][2] += m[2][2] * f;
+}
+
+BLI_INLINE void madd_m3_m3m3fl(float r[3][3], float a[3][3], float b[3][3], float f)
+{
+ r[0][0] = a[0][0] + b[0][0] * f;
+ r[0][1] = a[0][1] + b[0][1] * f;
@@ Diff output truncated at 10240 characters. @@
More information about the Bf-blender-cvs
mailing list