fixed issue with rand causing mem error

6 files changed, 104 insertions, 44 deletions

diff --git a/extern/softbody/src/admmpd_lattice.cpp b/extern/softbody/src/admmpd_lattice.cpp
index 84b65260ded..d97ab9c7bc8 100644
--- a/extern/softbody/src/admmpd_lattice.cpp
+++ b/extern/softbody/src/admmpd_lattice.cpp
@@ -9,12 +9,13 @@
 #include <unordered_map>
 #include <set>
 #include "BLI_task.h" // threading
+#include "BLI_assert.h"
 
 namespace admmpd {
 using namespace Eigen;
 
 typedef struct KeepTetThreadData {
-	AABBTree<double,3> *tree; // of embedded faces
+	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;
@@ -39,13 +40,22 @@ static void parallel_keep_tet(
 	// 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> traverser(
+	TetIntersectsMeshTraverse<double> tet_hits_mesh(
 		tet_pts, td->pts, td->faces);
-	bool hit = td->tree->traverse(traverser);
+	bool hit = td->tree->traverse(tet_hits_mesh);
 	if (!hit)
 	{
-		if (traverser.output.ray_hit_count % 2 == 1)
+		// 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; }
@@ -181,23 +191,28 @@ bool Lattice::generate(
 			bool keep = keep_tet[tet_idx];
 			if (keep)
 			{
-				refd_verts.emplace((*it)[0]);
-				refd_verts.emplace((*it)[1]);
-				refd_verts.emplace((*it)[2]);
-				refd_verts.emplace((*it)[3]);
+				const Vector4i &t = *it;
+				refd_verts.emplace(t[0]);
+				refd_verts.emplace(t[1]);
+				refd_verts.emplace(t[2]);
+				refd_verts.emplace(t[3]);
 				++it;
 			}
 			else { it = tets.erase(it); }
 		}
 
+		printf("Lattice: removed %d tets\n", tet_idx-(int)tets.size());
+
 	} // end remove unnecessary tets
 
 	// Copy data into matrices and remove unreferenced
 	{
 		// Computing a mapping of vertices from old to new
+		// Delete any unreferenced vertices
 		std::unordered_map<int,int> vtx_old_to_new;
 		int ntv0 = verts.size(); // original num verts
 		int ntv1 = refd_verts.size(); // reduced num verts
+		BLI_assert(ntv1 <= ntv0);
 		X->resize(ntv1,3);
 		int vtx_count = 0;
 		for (int i=0; i<ntv0; ++i)
@@ -247,7 +262,7 @@ static void parallel_point_in_tet(
 {
 	FindTetThreadData *td = (FindTetThreadData*)userdata;
 	Vector3d pt = td->pts->row(i);
-	PointInTetTraverse<double> traverser(pt, td->tet_x, td->tets);
+	PointInTetMeshTraverse<double> traverser(pt, td->tet_x, td->tets);
 	bool success = td->tree->traverse(traverser);
 	int tet_idx = traverser.output.prim;
 	if (success && tet_idx >= 0)
@@ -288,6 +303,7 @@ bool Lattice::compute_vtx_tet_mapping(
 	// barycoords for each embedded vertex
 	std::vector<AlignedBox<double,3> > tet_aabbs;
 	tet_aabbs.resize(nt);
+	Vector3d veta = Vector3d::Ones()*1e-12;
 	for (int i=0; i<nt; ++i)
 	{
 		tet_aabbs[i].setEmpty();
@@ -296,6 +312,8 @@ bool Lattice::compute_vtx_tet_mapping(
 		{
 			tet_aabbs[i].extend(x_->row(tet[j]).transpose());
 		}
+		tet_aabbs[i].extend(tet_aabbs[i].min()-veta);
+		tet_aabbs[i].extend(tet_aabbs[i].max()+veta);
 	}
 
 	AABBTree<double,3> tree;
diff --git a/extern/softbody/src/admmpd_solver.cpp b/extern/softbody/src/admmpd_solver.cpp
index 706eb9dce11..03680425427 100644
--- a/extern/softbody/src/admmpd_solver.cpp
+++ b/extern/softbody/src/admmpd_solver.cpp
@@ -11,6 +11,7 @@
 
 #include <stdio.h>
 #include <iostream>
+#include <unordered_map>
 
 #include "BLI_task.h" // threading
 #include "BLI_assert.h"
@@ -163,12 +164,12 @@ void Solver::update_constraints(
     std::vector<Eigen::Triplet<double> > trips_z;
 
 	// TODO collision detection
-	FloorCollider().jacobian(
-		&data->x,
-		&trips_x,
-		&trips_y,
-		&trips_z,
-		&l_coeffs);
+//	FloorCollider().jacobian(
+//		&data->x,
+//		&trips_x,
+//		&trips_y,
+//		&trips_z,
+//		&l_coeffs);
 
 	// Check number of constraints.
 	// If no constraints, clear Jacobian.
@@ -418,6 +419,16 @@ void Solver::compute_masses(
 		data->m[tet[2]] += tet_mass / 4.f;
 		data->m[tet[3]] += tet_mass / 4.f;
 	}
+	// Verify masses
+	int nx = data->m.rows();
+	for (int i=0; i<nx; ++i)
+	{
+		if (data->m[i] <= 0.0)
+		{
+			printf("**Solver::compute_masses Error: unreferenced vertex\n");
+			data->m[i]=1;
+		}
+	}
 } // end compute masses
 
 void Solver::append_energies(
diff --git a/extern/softbody/src/admmpd_solver.h b/extern/softbody/src/admmpd_solver.h
index b10105eba84..e61ddfeb82e 100644
--- a/extern/softbody/src/admmpd_solver.h
+++ b/extern/softbody/src/admmpd_solver.h
@@ -12,23 +12,23 @@
 namespace admmpd {
 
 struct Options {
-    double timestep_s;
+    double timestep_s; // TODO: Figure out delta time from blender api
     int max_admm_iters;
     int max_cg_iters;
     double mult_k; // stiffness multiplier for constraints
     double min_res; // min residual for CG solver
     double density_kgm3; // unit-volume density
-    double youngs; // Young's modulus
-    double poisson; // Poisson ratio
+    double youngs; // Young's modulus // TODO variable per-tet
+    double poisson; // Poisson ratio // TODO variable per-tet
     Eigen::Vector3d grav;
     Options() :
-        timestep_s(1.0/100.0), // TODO: Figure out delta time from blender api
+        timestep_s(1.0/100.0),
         max_admm_iters(20),
         max_cg_iters(10),
         mult_k(1.0),
         min_res(1e-4),
         density_kgm3(1100),
-        youngs(10000000),
+        youngs(100000000),
         poisson(0.399),
         grav(0,0,-9.8)
         {}
diff --git a/intern/softbody/admmpd_api.cpp b/intern/softbody/admmpd_api.cpp
index 1895b973fef..cdcede4ca72 100644
--- a/intern/softbody/admmpd_api.cpp
+++ b/intern/softbody/admmpd_api.cpp
@@ -47,9 +47,7 @@ void admmpd_dealloc(ADMMPDInterfaceData *iface)
   if (iface==NULL)
     return;
 
-  iface->totverts = 0;
-  iface->mesh_totverts = 0;
-  iface->mesh_totfaces = 0;
+  iface->totverts = 0; // output vertices
 
   if (iface->data)
   {
@@ -144,7 +142,6 @@ static int admmpd_init_with_lattice(
 
 int admmpd_init(ADMMPDInterfaceData *iface, float *in_verts, unsigned int *in_faces)
 {
-
   if (iface==NULL)
     return 0;
   if (in_verts==NULL || in_faces==NULL)
@@ -152,6 +149,9 @@ int admmpd_init(ADMMPDInterfaceData *iface, float *in_verts, unsigned int *in_fa
   if (iface->mesh_totverts<=0 || iface->mesh_totfaces<=0)
     return 0;
 
+  // Delete any existing data
+  admmpd_dealloc(iface);
+
   // Generate solver data
   iface->data = new ADMMPDInternalData();
   iface->data->options = new admmpd::Options();
diff --git a/intern/softbody/admmpd_api.h b/intern/softbody/admmpd_api.h
index b90a89360bd..e08e527b0e1 100644
--- a/intern/softbody/admmpd_api.h
+++ b/intern/softbody/admmpd_api.h
@@ -57,7 +57,10 @@ void admmpd_copy_from_bodypoint(ADMMPDInterfaceData*, const BodyPoint *pts);
 // Copies internal vertex position and velocity data
 // to BodyPoints (from SoftBody) AND surface mesh vertices.
 // If pts or vertexCos is null, its skipped
-void admmpd_copy_to_bodypoint_and_object(ADMMPDInterfaceData*, BodyPoint *pts, float (*vertexCos)[3]);
+void admmpd_copy_to_bodypoint_and_object(
+    ADMMPDInterfaceData*,
+    BodyPoint *pts,
+    float (*vertexCos)[3]);
 
 // Copies out_verts and out_verts to internal data
 // Performs solve over the time step
diff --git a/source/blender/blenkernel/intern/softbody.c b/source/blender/blenkernel/intern/softbody.c
index 5d94e9014f3..d97b76ba1e3 100644
--- a/source/blender/blenkernel/intern/softbody.c
+++ b/source/blender/blenkernel/intern/softbody.c
@@ -3580,7 +3580,7 @@ static void init_admmpd_interface(Object *ob, float (*vertexCos)[3])
     MEM_freeN(sb->admmpd);
     sb->admmpd = NULL;
   }
-  sb->admmpd = (ADMMPDInterfaceData*)MEM_callocN(sizeof(ADMMPDInterfaceData), "SoftBody_ADMMPD");
+  sb->admmpd = MEM_callocN(sizeof(ADMMPDInterfaceData), "SoftBody_ADMMPD");
 
   // Resize data
   int totfaces = poly_to_tri_count(me->totpoly, me->totloop);
@@ -3616,7 +3616,7 @@ static void init_admmpd_interface(Object *ob, float (*vertexCos)[3])
   looptri = NULL;
 
   // Initalize solver and free tmp data
-  sb->admmpd->init_mode = 1;
+  sb->admmpd->init_mode = 1; // 0=tetmesh, 1=lattice
   admmpd_init(sb->admmpd, in_verts, in_faces);
   MEM_freeN(in_verts);
   MEM_freeN(in_faces);
@@ -3630,9 +3630,12 @@ static void init_admmpd_interface(Object *ob, float (*vertexCos)[3])
     }
     sb->totpoint = sb->admmpd->totverts;
     sb->totspring = 0;
-    sb->bpoint = MEM_callocN(sb->totpoint * sizeof(BodyPoint), "bodypoint");
+    sb->bpoint = (BodyPoint*)MEM_callocN(sb->totpoint*sizeof(BodyPoint), "bodypoint");
   }
 
+  // admmpd_init will have created tets with extra vertices besides
+  // the ones on the surface (vertexCos). We'll store those vertices
+  // in SoftBody::bpoint.
   admmpd_copy_to_bodypoint_and_object(sb->admmpd,sb->bpoint,NULL);
 }
 
@@ -3660,25 +3663,17 @@ void sbObjectStep_admmpd(
   float timescale = 0.f;
   BKE_ptcache_id_from_softbody(&pid, ob, sb);
   BKE_ptcache_id_time(&pid, scene, framenr, &startframe, &endframe, &timescale);
-  
-  // BKE_ptcache_validate:
-  //  sets flag that simulation is valid (|=)
-  //  and sets simframe to framenr
-  // Calling BKE_ptcache_id_reset calls
-  //  free_softbody_intern
 
-  // Reset cache
-  if (sb->admmpd == NULL ||
-    sb->admmpd->mesh_totverts != numVerts)
-  {
-    cache->flag |= PTCACHE_OUTDATED;
-    BKE_ptcache_id_reset(scene, &pid, PTCACHE_RESET_OUTDATED);
-    BKE_ptcache_validate(cache, 0);
-    cache->last_exact = 0;
-    cache->flag &= ~PTCACHE_REDO_NEEDED;
+  // check for changes in mesh, should only happen in case the mesh
+  // structure changes during an animation.
+  // Hopefully we're short-circuiting the boolean here.
+  if (sb->admmpd != NULL &&
+      sb->admmpd->mesh_totverts != numVerts) {
+    BKE_ptcache_invalidate(cache);
+    return;
   }
 
-  // Clamp simulation frame
+  // clamp frame ranges
   if (framenr < startframe) {
     BKE_ptcache_invalidate(cache);
     return;
@@ -3687,6 +3682,22 @@ void sbObjectStep_admmpd(
     framenr = endframe;
   }
 
+  // BKE_ptcache_validate:
+  //  sets flag that simulation is valid (|=)
+  //  and sets simframe to framenr
+  // Calling BKE_ptcache_id_reset calls
+  //  free_softbody_intern
+
+  // Reset cache
+//  if (sb->admmpd && sb->admmpd->mesh_totverts != numVerts)
+//  {
+//    cache->flag |= PTCACHE_OUTDATED;
+//    BKE_ptcache_id_reset(scene, &pid, PTCACHE_RESET_OUTDATED);
+//    BKE_ptcache_validate(cache, 0);
+//    cache->last_exact = 0;
+//    cache->flag &= ~PTCACHE_REDO_NEEDED;
+//  }
+
   // If it's the first frame we reset the cache and data
   if (framenr == startframe || sb->admmpd == NULL) {
     BKE_ptcache_id_reset(scene, &pid, PTCACHE_RESET_OUTDATED);
@@ -3710,8 +3721,17 @@ void sbObjectStep_admmpd(
       cache_result == PTCACHE_READ_INTERPOLATED ||
       (!can_simulate && cache_result == PTCACHE_READ_OLD)) {
 
+      // Copies bodypoint vertices to ADMM internal verts
       admmpd_copy_from_bodypoint(sb->admmpd,sb->bpoint);
+      // Maps ADMM internal verts to surface verts
       admmpd_copy_to_bodypoint_and_object(sb->admmpd,NULL,vertexCos);
+      if(sb->local==0)
+      {
+        for (int i=0; i<me->totvert; ++i)
+        {
+          mul_m4_v3(ob->imat, vertexCos[i]);
+        }
+      }
 
       BKE_ptcache_validate(cache, framenr);
       if (cache_result == PTCACHE_READ_INTERPOLATED &&
@@ -3724,6 +3744,7 @@ void sbObjectStep_admmpd(
     }
     else if (cache_result == PTCACHE_READ_OLD) {} // pass
     else if (cache->flag & PTCACHE_BAKED) {
+
       // if baked and nothing in cache, do nothing
       if (can_write_cache) {
         BKE_ptcache_invalidate(cache);
@@ -3732,7 +3753,9 @@ void sbObjectStep_admmpd(
     }
 
     if (!can_simulate)
+    {
       return;
+    }
   } // end read from cache
 
   // if on second frame, write cache for first frame
@@ -3744,8 +3767,13 @@ void sbObjectStep_admmpd(
   // Update ADMMPD interface variables from cache
   // and perform a solve.
   admmpd_copy_from_bodypoint(sb->admmpd,sb->bpoint);
+
   admmpd_solve(sb->admmpd);
   admmpd_copy_to_bodypoint_and_object(sb->admmpd,sb->bpoint,vertexCos);
+  for (int i=0; i<me->totvert; ++i)
+  { // TODO move this to admmpd API
+    mul_m4_v3(ob->imat, vertexCos[i]);
+  }
 
   // Write cache
   BKE_ptcache_validate(cache, framenr);