[Bf-blender-cvs] [9ceb2981560] soc-2020-soft-body: still pretty buggy. The solver should be working, but lattice is too lousy to see. The interface between blender and solver needs work too, doesn't seem to reset as often as I think it should be?
over0219
noreply at git.blender.org
Wed Jun 3 06:26:45 CEST 2020
Commit: 9ceb298156044e616bcea97b3c82f1c416ec4385
Author: over0219
Date: Tue Jun 2 23:26:39 2020 -0500
Branches: soc-2020-soft-body
https://developer.blender.org/rB9ceb298156044e616bcea97b3c82f1c416ec4385
still pretty buggy. The solver should be working, but lattice is too lousy to see. The interface between blender and solver needs work too, doesn't seem to reset as often as I think it should be?
===================================================================
M extern/softbody/CMakeLists.txt
A extern/softbody/src/admmpd_collision.cpp
A extern/softbody/src/admmpd_collision.h
M extern/softbody/src/admmpd_energy.cpp
M extern/softbody/src/admmpd_lattice.cpp
M extern/softbody/src/admmpd_math.cpp
M extern/softbody/src/admmpd_math.h
M extern/softbody/src/admmpd_solver.cpp
M extern/softbody/src/admmpd_solver.h
M intern/softbody/admmpd_api.cpp
===================================================================
diff --git a/extern/softbody/CMakeLists.txt b/extern/softbody/CMakeLists.txt
index 08a36882f28..8f6ec97d7dd 100644
--- a/extern/softbody/CMakeLists.txt
+++ b/extern/softbody/CMakeLists.txt
@@ -35,6 +35,8 @@ set(SRC
src/admmpd_lattice.cpp
src/admmpd_solver.h
src/admmpd_solver.cpp
+ src/admmpd_collision.h
+ src/admmpd_collision.cpp
)
set(LIB
diff --git a/extern/softbody/src/admmpd_collision.cpp b/extern/softbody/src/admmpd_collision.cpp
new file mode 100644
index 00000000000..2bce56fa9f0
--- /dev/null
+++ b/extern/softbody/src/admmpd_collision.cpp
@@ -0,0 +1,33 @@
+
+
+#include "admmpd_collision.h"
+
+namespace admmpd {
+
+void FloorCollider::detect(const Eigen::MatrixXd *x)
+{
+ (void)(x);
+ // Can just do detection in jacobian I guess.
+}
+
+void FloorCollider::jacobian(
+ const Eigen::MatrixXd *x,
+ std::vector<Eigen::Triplet<double> > *trips_x,
+ std::vector<Eigen::Triplet<double> > *trips_y,
+ std::vector<Eigen::Triplet<double> > *trips_z,
+ std::vector<double> *l)
+{
+ int nx = x->rows();
+ for (int i=0; i<nx; ++i)
+ {
+ Eigen::Vector3d p = x->row(i);
+ if (p[2]>0)
+ continue;
+
+ int c_idx = l->size();
+ trips_z->emplace_back(c_idx,i,1.0);
+ l->emplace_back(0);
+ }
+} // end floor collider Jacobian
+
+} // namespace admmpd
\ No newline at end of file
diff --git a/extern/softbody/src/admmpd_collision.h b/extern/softbody/src/admmpd_collision.h
new file mode 100644
index 00000000000..274be91be6d
--- /dev/null
+++ b/extern/softbody/src/admmpd_collision.h
@@ -0,0 +1,42 @@
+
+
+
+#ifndef _ADMMPD_COLLISION_H
+#define _ADMMPD_COLLISION_H
+
+#include <Eigen/Sparse>
+#include <vector>
+
+namespace admmpd {
+
+// I'll update this class/structure another day.
+// For now let's get something in place to do floor collisions.
+// Probably will work better to use uber-collision class for
+// all self and obstacle collisions, reducing the amount of
+// for-all vertices loops.
+class Collider {
+public:
+ virtual void detect(
+ const Eigen::MatrixXd *x) = 0;
+ virtual void jacobian(
+ const Eigen::MatrixXd *x,
+ std::vector<Eigen::Triplet<double> > *trips_x,
+ std::vector<Eigen::Triplet<double> > *trips_y,
+ std::vector<Eigen::Triplet<double> > *trips_z,
+ std::vector<double> *l) = 0;
+};
+
+class FloorCollider : public Collider {
+public:
+ void detect(const Eigen::MatrixXd *x);
+ void jacobian(
+ const Eigen::MatrixXd *x,
+ std::vector<Eigen::Triplet<double> > *trips_x,
+ std::vector<Eigen::Triplet<double> > *trips_y,
+ std::vector<Eigen::Triplet<double> > *trips_z,
+ std::vector<double> *l);
+};
+
+} // namespace admmpd
+
+#endif //_ADMMPD_COLLISION_H
diff --git a/extern/softbody/src/admmpd_energy.cpp b/extern/softbody/src/admmpd_energy.cpp
index 978fd565aa4..ea7c68b6dca 100644
--- a/extern/softbody/src/admmpd_energy.cpp
+++ b/extern/softbody/src/admmpd_energy.cpp
@@ -8,7 +8,7 @@ using namespace Eigen;
Lame::Lame() : m_model(0)
{
- set_from_youngs_poisson(10000000,0.299);
+ set_from_youngs_poisson(10000000,0.399);
}
void Lame::set_from_youngs_poisson(double youngs, double poisson)
diff --git a/extern/softbody/src/admmpd_lattice.cpp b/extern/softbody/src/admmpd_lattice.cpp
index 4323f8fce33..be7a5bc63d2 100644
--- a/extern/softbody/src/admmpd_lattice.cpp
+++ b/extern/softbody/src/admmpd_lattice.cpp
@@ -3,6 +3,7 @@
#include "admmpd_lattice.h"
#include "admmpd_math.h"
+#include <iostream>
//#include "vdb.h"
diff --git a/extern/softbody/src/admmpd_math.cpp b/extern/softbody/src/admmpd_math.cpp
index 98b1fcc2520..7c76dfb8b3a 100644
--- a/extern/softbody/src/admmpd_math.cpp
+++ b/extern/softbody/src/admmpd_math.cpp
@@ -6,26 +6,26 @@ namespace admmpd {
namespace barycoords {
Eigen::Matrix<double,4,1> point_tet(
- const Eigen::Matrix<double,3,1> &p,
- const Eigen::Matrix<double,3,1> &a,
- const Eigen::Matrix<double,3,1> &b,
- const Eigen::Matrix<double,3,1> &c,
- const Eigen::Matrix<double,3,1> &d)
+ const Eigen::Vector3d &p,
+ const Eigen::Vector3d &a,
+ const Eigen::Vector3d &b,
+ const Eigen::Vector3d &c,
+ const Eigen::Vector3d &d)
{
auto scalar_triple_product = [](
- const Eigen::Matrix<double,3,1> &u,
- const Eigen::Matrix<double,3,1> &v,
- const Eigen::Matrix<double,3,1> &w )
+ const Eigen::Vector3d &u,
+ const Eigen::Vector3d &v,
+ const Eigen::Vector3d &w )
{
return u.dot(v.cross(w));
};
- Eigen::Matrix<double,3,1> ap = p - a;
- Eigen::Matrix<double,3,1> bp = p - b;
- Eigen::Matrix<double,3,1> ab = b - a;
- Eigen::Matrix<double,3,1> ac = c - a;
- Eigen::Matrix<double,3,1> ad = d - a;
- Eigen::Matrix<double,3,1> bc = c - b;
- Eigen::Matrix<double,3,1> bd = d - b;
+ Eigen::Vector3d ap = p - a;
+ Eigen::Vector3d bp = p - b;
+ Eigen::Vector3d ab = b - a;
+ Eigen::Vector3d ac = c - a;
+ Eigen::Vector3d ad = d - a;
+ Eigen::Vector3d bc = c - b;
+ Eigen::Vector3d bd = d - b;
double va6 = scalar_triple_product(bp, bd, bc);
double vb6 = scalar_triple_product(ap, ac, ad);
double vc6 = scalar_triple_product(ap, ad, ab);
@@ -37,11 +37,11 @@ namespace barycoords {
// Checks that it's on the "correct" side of the normal
// for each face of the tet. Assumes winding points inward.
bool point_in_tet(
- const Eigen::Matrix<double,3,1> &p,
- const Eigen::Matrix<double,3,1> &a,
- const Eigen::Matrix<double,3,1> &b,
- const Eigen::Matrix<double,3,1> &c,
- const Eigen::Matrix<double,3,1> &d)
+ const Eigen::Vector3d &p,
+ const Eigen::Vector3d &a,
+ const Eigen::Vector3d &b,
+ const Eigen::Vector3d &c,
+ const Eigen::Vector3d &d)
{
using namespace Eigen;
auto check_face = [](
diff --git a/extern/softbody/src/admmpd_math.h b/extern/softbody/src/admmpd_math.h
index d777df5ce2e..25ff2d3114f 100644
--- a/extern/softbody/src/admmpd_math.h
+++ b/extern/softbody/src/admmpd_math.h
@@ -9,19 +9,19 @@
namespace admmpd {
namespace barycoords {
-Eigen::Matrix<double,4,1> point_tet(
- const Eigen::Matrix<double,3,1> &p,
- const Eigen::Matrix<double,3,1> &a,
- const Eigen::Matrix<double,3,1> &b,
- const Eigen::Matrix<double,3,1> &c,
- const Eigen::Matrix<double,3,1> &d);
+Eigen::Vector4d point_tet(
+ const Eigen::Vector3d &p,
+ const Eigen::Vector3d &a,
+ const Eigen::Vector3d &b,
+ const Eigen::Vector3d &c,
+ const Eigen::Vector3d &d);
bool point_in_tet(
- const Eigen::Matrix<double,3,1> &p,
- const Eigen::Matrix<double,3,1> &a,
- const Eigen::Matrix<double,3,1> &b,
- const Eigen::Matrix<double,3,1> &c,
- const Eigen::Matrix<double,3,1> &d);
+ const Eigen::Vector3d &p,
+ const Eigen::Vector3d &a,
+ const Eigen::Vector3d &b,
+ const Eigen::Vector3d &c,
+ const Eigen::Vector3d &d);
} // namespace barycoords
diff --git a/extern/softbody/src/admmpd_solver.cpp b/extern/softbody/src/admmpd_solver.cpp
index eb2222ef3bb..052008121b9 100644
--- a/extern/softbody/src/admmpd_solver.cpp
+++ b/extern/softbody/src/admmpd_solver.cpp
@@ -4,11 +4,13 @@
#include "admmpd_solver.h"
#include "admmpd_lattice.h"
#include "admmpd_energy.h"
+#include "admmpd_collision.h"
#include <Eigen/Geometry>
#include <Eigen/Sparse>
#include <stdio.h>
+#include <iostream>
namespace admmpd {
using namespace Eigen;
@@ -29,34 +31,6 @@ bool Solver::init(
return true;
} // end init
-void Solver::init_solve(
- const ADMMPD_Options *options,
- ADMMPD_Data *data)
-{
- double dt = std::max(0.0, options->timestep_s);
- int nx = data->x.rows();
- if (data->M_xbar.rows() != nx)
- data->M_xbar.resize(nx,3);
-
- // velocity and position
- data->x_start = data->x;
- for( int i=0; i<nx; ++i )
- {
- data->v.row(i) += options->grav;
- data->M_xbar.row(i) =
- data->m[i] * data->x.row(i) +
- dt*data->m[i]*data->v.row(i);
- }
-
- // ADMM variables
- data->Dx.noalias() = data->D * data->x;
- data->z = data->Dx;
- data->z_prev = data->z;
- data->u.setZero();
- data->u_prev.setZero();
-
-} // end init solve
-
int Solver::solve(
const ADMMPD_Options *options,
ADMMPD_Data *data)
@@ -67,17 +41,14 @@ int Solver::solve(
init_solve(options,data);
int ne = data->rest_volumes.size();
- Lame lame;
+ Lame lame; // TODO lame params as input
// Begin solver loop
int iters = 0;
- int max_iters = options->max_iters < 0 ? 10000 : options->max_iters;
- for (; iters < max_iters; ++iters)
+ for (; iters < options->max_admm_iters; ++iters)
{
// Local step
data->Dx.noalias() = data->D * data->x;
- data->z_prev.noalias() = data->z;
- data->u_prev.noalias() = data->u;
for(int i=0; i<ne; ++i)
{
EnergyTerm().update(
@@ -92,18 +63,162 @@ int Solver::solve(
}
// Global step
+ update_constraints(options,data);
data->b.noalias() = data->M_xbar + data->DtW2*(data->z-data->u);
- data->x.noalias() = data->ldltA.solve(data->b);
+ solve_conjugate_gradients(options,data);
} // end solver iters
- double dt = std::max(0.0, options->timestep_s);
+ double dt = options->timestep_s;
if (dt > 0.0)
data->v.noalias() = (data->x-data->x_start)*(1.0/dt);
return iters;
} // end solve
+void Solver::init_solve(
+ const ADMMPD_Options *options,
+ ADMMPD_Data *data)
+{
+ int nx = data->x.rows();
+ if (data->M_xbar.rows() != nx)
+ data->M_xbar.resize(nx,3);
+
+ // velocity and position
+ double dt = std::max(0.0, options->timestep_s);
+ data->x_start = data->x;
+ for( int i=0; i<nx; ++i )
+ {
+ data->v.row(i) += options->grav;
+ data->M_xbar.row(i) =
+ data->m[i] * data->x.row(i) +
+ dt*data->m[i]*data->v.row(i);
+ }
+
+ // ADMM variables
+ data->Dx.noalias() = data->D * data->x;
+ data->z = data->Dx;
+ data->u.setZero();
+
+} // end init solve
+
+void Solver::update_constraints(
+ const ADMMPD_Options *options,
+ ADMMPD_Data *data)
+{
+
+ std::vector<double> l_coeffs;
+ std::vector<Eigen::Triplet<double> > trips_x;
+ std::vector<Eigen::Triplet<double> > trips_y;
+ std::vector<Eigen::Triplet<double> > trips_z;
+
+ // TODO collision detection
+ FloorCollider().jacobian(
+ &data->x,
+ &trips_x,
+ &trips_y,
+ &trips_z,
+ &l_coeffs);
+
+ // Check number of constraints.
+ // If no constraints, clear Jacobian.
+ int nx = data->x.rows();
+ int nc = l_coeffs.size();
+ if
@@ Diff output truncated at 10240 characters. @@
More information about the Bf-blender-cvs
mailing list