[Bf-blender-cvs] [6191b3d2db2] soc-2020-soft-body: octree for lattice gen and smaller dt with CD only at init solve

over0219 noreply at git.blender.org
Wed Jul 15 21:58:52 CEST 2020


Commit: 6191b3d2db2d8788e6d8f12cc6e433d247377387
Author: over0219
Date:   Wed Jul 15 14:58:46 2020 -0500
Branches: soc-2020-soft-body
https://developer.blender.org/rB6191b3d2db2d8788e6d8f12cc6e433d247377387

octree for lattice gen and smaller dt with CD only at init solve

===================================================================

M	extern/softbody/src/admmpd_bvh.cpp
M	extern/softbody/src/admmpd_bvh.h
M	extern/softbody/src/admmpd_bvh_traverse.h
M	extern/softbody/src/admmpd_collision.cpp
M	extern/softbody/src/admmpd_embeddedmesh.cpp
M	extern/softbody/src/admmpd_solver.cpp
M	extern/softbody/src/admmpd_solver.h
M	extern/softbody/src/admmpd_types.h

===================================================================

diff --git a/extern/softbody/src/admmpd_bvh.cpp b/extern/softbody/src/admmpd_bvh.cpp
index fffc8c33a74..6ea7a676d27 100644
--- a/extern/softbody/src/admmpd_bvh.cpp
+++ b/extern/softbody/src/admmpd_bvh.cpp
@@ -209,6 +209,133 @@ bool AABBTree<T,DIM>::traverse_children(
 
 } // end traverse children
 
+template<typename T, int DIM>
+void Octree<T,DIM>::clear()
+{
+    m_root = std::make_shared<Node>();
+}
+
+template<typename T, int DIM>
+typename Octree<T,DIM>::AABB Octree<T,DIM>::bounds() const
+{
+    if (m_root)
+        return m_root->bounds();
+    return AABB();
+}
+
+template<typename T, int DIM>
+void Octree<T,DIM>::init(const MatrixXT *V, const Eigen::MatrixXi *F, int stopdepth)
+{
+    BLI_assert(V != nullptr);
+    BLI_assert(F != nullptr);
+    BLI_assert(F->cols()==3);
+    m_root = std::make_shared<Node>();
+
+    if (DIM !=3 || F->cols()!=3)
+    {
+        return;
+    }
+
+    int nf = F->rows();
+    AABB global_box;
+    std::vector<AABB> boxes(nf);
+    std::vector<int> queue(nf);
+    for (int i=0; i<nf; ++i)
+    {
+        Eigen::RowVector3i f = F->row(i);
+		queue[i]=i;
+		boxes[i].extend(V->row(f[0]).transpose());
+		boxes[i].extend(V->row(f[1]).transpose());
+		boxes[i].extend(V->row(f[2]).transpose());
+		global_box.extend(boxes[i]);
+    }
+
+	T halfwidth = global_box.sizes().maxCoeff()*0.5;
+    m_root.reset(
+        create_children(global_box.center(),halfwidth,stopdepth,V,F,queue,boxes)
+    );
+}
+
+template<typename T, int DIM>
+typename Octree<T,DIM>::Node* Octree<T,DIM>::create_children(
+    const VecType &center, T halfwidth, int stopdepth,
+    const MatrixXT *V, const Eigen::MatrixXi *F,
+    const std::vector<int> &queue,
+    const std::vector<AABB> &boxes)
+{
+    BLI_assert((int)queue.size()>0);
+    BLI_assert((int)prim_boxes.size()>0);
+    BLI_assert(F != nullptr);
+    BLI_assert(V != nullptr);
+    BLI_assert(F->cols()==3);
+    BLI_assert(V->cols()==3);
+    if (stopdepth >= 0)
+    {
+        Node *node = new Node();
+        node->center = center;
+        node->halfwidth = halfwidth;
+        node->prims.clear();
+        AABB box = node->bounds();
+        // Set list of intersected prims
+        int nq = queue.size();
+        for (int i=0; i<nq; ++i)
+        {
+            int p_idx = queue[i];
+            if (box.intersects(boxes[p_idx]))
+                node->prims.emplace_back(p_idx);
+        }
+        // Create children only if prims intersect
+        int np = node->prims.size();
+        for (int i=0; i<nchild && np>0; ++i)
+        {
+
+            T step = node->halfwidth * 0.5;
+            VecType offset = VecType::Zero();
+            offset[0] = ((i & 1) ? step : -step);
+            offset[1] = ((i & 2) ? step : -step);
+            if (DIM==3) offset[2] = ((i & 4) ? step : -step);
+
+            node->children[i] = create_children(
+                node->center+offset, step, stopdepth-1,
+                V, F, node->prims, boxes);
+        }
+        return node;
+    }
+    return nullptr;
+}
+
+template<typename T, int DIM>
+Octree<T,DIM>::Node::Node() :
+    center(VecType::Zero()),
+    halfwidth(0)
+{
+	for (int i=0; i<nchild; ++i)
+		children[i] = nullptr;
+}
+
+template<typename T, int DIM>
+Octree<T,DIM>::Node::~Node()
+{
+    for (int i=0; i<nchild; ++i)
+	    if (children[i] != nullptr)
+			delete children[i];
+}
+
+template<typename T, int DIM>
+bool Octree<T,DIM>::Node::is_leaf() const
+{
+    return children[0] == nullptr;
+}
+
+template<typename T, int DIM>
+typename Octree<T,DIM>::AABB Octree<T,DIM>::Node::bounds() const
+{
+    AABB box;
+    box.extend(center + VecType::Ones()*halfwidth);
+    box.extend(center - VecType::Ones()*halfwidth);
+    return box;
+}
+
 // Compile types
 template class admmpd::AABBTree<double,2>;
 template class admmpd::AABBTree<double,3>;
@@ -218,5 +345,9 @@ template class admmpd::TraverserFromFunctionPtrs<double,2>;
 template class admmpd::TraverserFromFunctionPtrs<double,3>;
 template class admmpd::TraverserFromFunctionPtrs<float,2>;
 template class admmpd::TraverserFromFunctionPtrs<float,3>;
+template class admmpd::Octree<double,2>;
+template class admmpd::Octree<double,3>;
+template class admmpd::Octree<float,2>;
+template class admmpd::Octree<float,3>;
 
 } // namespace admmpd
\ No newline at end of file
diff --git a/extern/softbody/src/admmpd_bvh.h b/extern/softbody/src/admmpd_bvh.h
index 19f47545435..20b810ac7be 100644
--- a/extern/softbody/src/admmpd_bvh.h
+++ b/extern/softbody/src/admmpd_bvh.h
@@ -42,8 +42,6 @@ public:
 		std::function<void(const AABB&, bool&, const AABB&, bool&, bool&)> t,
 		std::function<bool(const AABB&, int)> s) const;
 
-protected:
-
 	struct Node
 	{
 		AABB aabb;
@@ -60,6 +58,8 @@ protected:
 		}
 	};
 
+protected:
+
 	std::shared_ptr<Node> root;
 
 	void create_children(
@@ -77,6 +77,58 @@ protected:
 
 }; // class AABBtree
 
+
+// Octree is actually a quadtree if DIM=2
+template<typename T, int DIM>
+class Octree
+{
+protected:
+	typedef Eigen::AlignedBox<T,DIM> AABB;
+	typedef Eigen::Matrix<T,DIM,1> VecType;
+	typedef Eigen::Matrix<T,Eigen::Dynamic,Eigen::Dynamic> MatrixXT;
+	static const int nchild = std::pow(2,DIM);
+public:
+
+	// Removes all BVH data
+	void clear();
+
+	// Creates the Octree up to depth with smart subdivision
+	// (only create children if it contains prims) and does not
+	// create a cell if it is outside the mesh.
+	// ** Assumes a closed mesh and only defined for 3D
+	void init(const MatrixXT *V, const Eigen::MatrixXi *F, int stopdepth);
+
+	// Returns bounding box of the root node
+	AABB bounds() const;
+
+	struct Node
+	{
+		VecType center;
+		T halfwidth;
+		Node *children[4*DIM];
+		std::vector<int> prims; // includes childen
+		bool is_leaf() const;
+		AABB bounds() const;
+		Node();
+		~Node();
+	};
+
+	// Return ptr to the root node
+	// Becomes invalidated after init()
+	std::shared_ptr<Node> root() { return m_root; }
+
+protected:
+
+	std::shared_ptr<Node> m_root;
+
+	Node* create_children(
+		const VecType &center, T halfwidth, int stopdepth,
+		const MatrixXT *V, const Eigen::MatrixXi *F,
+		const std::vector<int> &queue,
+		const std::vector<AABB> &boxes);
+
+}; // class Octree
+
 } // namespace admmpd
 
 #endif // ADMMPD_BVH_H_
diff --git a/extern/softbody/src/admmpd_bvh_traverse.h b/extern/softbody/src/admmpd_bvh_traverse.h
index ecbb85f3277..dce83d07a20 100644
--- a/extern/softbody/src/admmpd_bvh_traverse.h
+++ b/extern/softbody/src/admmpd_bvh_traverse.h
@@ -136,6 +136,7 @@ public:
 	struct Output {
 		std::vector< std::pair<int,T> > hits; // [prim,t]
 		int num_hits() const { return hits.size(); }
+		bool is_inside() const { return hits.size()%2==1; }
 	} output;
 
 	PointInTriangleMeshTraverse(
diff --git a/extern/softbody/src/admmpd_collision.cpp b/extern/softbody/src/admmpd_collision.cpp
index 5e14bc3a756..a05b27287cc 100644
--- a/extern/softbody/src/admmpd_collision.cpp
+++ b/extern/softbody/src/admmpd_collision.cpp
@@ -308,9 +308,10 @@ void EmbeddedMeshCollision::linearize(
 
 	for (int i=0; i<np; ++i)
 	{
-		Vector2i pair_idx = vf_pairs[i];
+		const Vector2i &pair_idx = vf_pairs[i];
 		VFCollisionPair &pair = per_vertex_pairs[pair_idx[0]][pair_idx[1]];
 		int emb_p_idx = pair.p_idx;
+		Vector3d p_pt = mesh->get_mapped_vertex(x,emb_p_idx);
 
 		//
 		// If we collided with an obstacle
@@ -334,8 +335,17 @@ void EmbeddedMeshCollision::linearize(
 				};
 				q_n = ((q_tris[1]-q_tris[0]).cross(q_tris[2]-q_tris[0]));
 				q_n.normalize();
+
+				// Update constraint linearization
+				pair.q_pt = geom::point_on_triangle<double>(p_pt,
+					q_tris[0], q_tris[1], q_tris[2]);
 			}
 
+			// Is the constraint active?
+			bool active = (p_pt-pair.q_pt).dot(q_n) <= 0.0;
+			if (!active)
+				continue;
+
 			// Get the four deforming verts that embed
 			// the surface vertices, and add constraints on those.
 			RowVector4d bary = mesh->emb_barys.row(emb_p_idx);
diff --git a/extern/softbody/src/admmpd_embeddedmesh.cpp b/extern/softbody/src/admmpd_embeddedmesh.cpp
index b16d92b640a..7e677b42455 100644
--- a/extern/softbody/src/admmpd_embeddedmesh.cpp
+++ b/extern/softbody/src/admmpd_embeddedmesh.cpp
@@ -18,55 +18,6 @@
 namespace admmpd {
 using namespace Eigen;
 
-typedef struct KeepTetThreadData {
-	const AABBTree<double,3> *tree; // of embedded faces
-	const MatrixXd *pts; // of embedded verts
-	const MatrixXi *faces; // embedded faces
-	const std::vector<Vector3d> *tet_x;
-	const std::vector<Vector4i> *tets;
-	std::vector<int> *keep_tet; // 0 or 1
-} KeepTetThreadData;
-
-static void parallel_keep_tet(
-	void *__restrict userdata,
-	const int i,
-	const TaskParallelTLS *__restrict UNUSED(tls))
-{
-	KeepTetThreadData *td = (KeepTetThreadData*)userdata;
-	RowVector4i tet = td->tets->at(i);
-	Vector3d tet_pts[4] = {
-		td->tet_x->at(tet[0]),
-		td->tet_x->at(tet[1]),
-		td->tet_x->at(tet[2]),
-		td->tet_x->at(tet[3])
-	};
-
-	// Returns true if the tet intersects the
-	// surface mesh. Even if it doesn't, we want to keep
-	// the tet if it's totally enclosed in the mesh.
-	TetIntersectsMeshTraverse<double> tet_hits_mesh(
-		tet_pts, td->pts, td->faces);
-	bool hit = td->tree->traverse(tet_hits_mesh);
-	if (!hit)
-	{
-		// We only need to check if one vertex of the
-		// tet is inside the mesh. If a subset of
-		// vertices are inside the mesh, then there would
-		// be tri/tri intersection.
-		PointInTriangleMeshTraverse<double> pt_in_mesh(
-			tet_pts[0], td->pts, td->faces);
-		td->tree->traverse(pt_in_mesh);
-		if (pt_in_mesh.output.num_hits() % 2 == 1)
-		{
-			hit = true;
-		}
-	}
-
-	if (hit) { td->keep_tet->at(i) = 1; }
-	else { td->keep_tet->at(i) = 0; }
-
-} // end parallel test if keep tet
-
 // Gen lattice with subdivision
 struct LatticeData {
 	//SDF<double> *emb_sdf;
@@ -114,110 +65,87 @@ static inline void merge_close_vertices(LatticeData *data, double eps=1e-12)
 	}
 }
 
-static inline void subdivide_cube(
-	LatticeData *data,
-	const std::vector<Vector3d> &cv,
-	const std::vector<int> &faces,
-	int level)
+static inline void add_tets_from_box(
+	const Vector3d &min,
+	const Vector3d &max,
+	std::vector<Vector3d> &verts,
+	std::vector<Vector4i> &tets)
 {
-	BLI_assert((int)cv.size()==8);
-	auto add_tets_from_box = [&]()
+	std::vector<Vector3d> v

@@ Diff output truncated at 10240 characters. @@



More information about the Bf-blender-cvs mailing list