- parallelized unified solver using graph coloring

janbender · janbender · commit 7be1617945d1 · 2015-10-20T16:12:06.000+02:00
diff --git a/Changelog.txt b/Changelog.txt
@@ -1,3 +1,4 @@
+	- parallelized unified solver using graph coloring
 	- implemented unified solver for rigid bodies and deformable solids 
 	- added documentation for generic constraints
 	- added generic constraint demo
diff --git a/Demos/BarDemo/main.cpp b/Demos/BarDemo/main.cpp
@@ -138,7 +138,7 @@ void timeStep ()
 		return;
 
 	// Simulation code
-	for (unsigned int i = 0; i < 4; i++)
+	for (unsigned int i = 0; i < 8; i++)
 		sim.step(model);
 
 	for (unsigned int i = 0; i < model.getTetModels().size(); i++)
diff --git a/Demos/ClothDemo/main.cpp b/Demos/ClothDemo/main.cpp
@@ -55,8 +55,8 @@ void TW_CALL getVelocityUpdateMethod(void *value, void *clientData);
 SimulationModel model;
 TimeStepController sim;
 
-const int nRows = 30;
-const int nCols = 30;
+const int nRows = 50;
+const int nCols = 50;
 const float width = 10.0f;
 const float height = 10.0f;
 bool doPause = true;
@@ -161,7 +161,7 @@ void timeStep ()
 		return;
 
 	// Simulation code
-	for (unsigned int i = 0; i < 4; i++)
+	for (unsigned int i = 0; i < 8; i++)
 		sim.step(model);
 
 	for (unsigned int i = 0; i < model.getTriangleModels().size(); i++)
diff --git a/Demos/CouplingDemos/RigidBodyClothCouplingDemo.cpp b/Demos/CouplingDemos/RigidBodyClothCouplingDemo.cpp
@@ -191,7 +191,7 @@ void timeStep ()
 		return;
 
 	// Simulation code
-	for (unsigned int i = 0; i < 4; i++)
+	for (unsigned int i = 0; i < 8; i++)
 		sim.step(model);
 
 	for (unsigned int i = 0; i < model.getTriangleModels().size(); i++)
diff --git a/Demos/FluidDemo/main.cpp b/Demos/FluidDemo/main.cpp
@@ -147,7 +147,7 @@ void timeStep ()
 		return;
 
 	// Simulation code
-	for (unsigned int i = 0; i < 4; i++)
+	for (unsigned int i = 0; i < 8; i++)
 		simulation.step(model);
 }
 
diff --git a/Demos/GenericConstraintsDemos/GenericParticleConstraintsDemo.cpp b/Demos/GenericConstraintsDemos/GenericParticleConstraintsDemo.cpp
@@ -37,8 +37,8 @@ void TW_CALL getBendingStiffness(void *value, void *clientData);
 GenericConstraintsModel model;
 TimeStepController simulation;
 
-const int nRows = 20;
-const int nCols = 20;
+const int nRows = 30;
+const int nCols = 30;
 const float width = 10.0f;
 const float height = 10.0f;
 bool doPause = true;
@@ -126,7 +126,7 @@ void timeStep ()
 		return;
 
 	// Simulation code
-	for (unsigned int i = 0; i < 1; i++)
+	for (unsigned int i = 0; i < 8; i++)
 		simulation.step(model);
 
 	for (unsigned int i = 0; i < model.getTriangleModels().size(); i++)
diff --git a/Demos/RigidBodyDemos/ChainDemo.cpp b/Demos/RigidBodyDemos/ChainDemo.cpp
@@ -130,7 +130,7 @@ void timeStep ()
 		return;
 
 	// Simulation code
-	for (unsigned int i = 0; i < 4; i++)
+	for (unsigned int i = 0; i < 8; i++)
 		sim.step(model);
 }
 
diff --git a/Demos/RigidBodyDemos/JointDemo.cpp b/Demos/RigidBodyDemos/JointDemo.cpp
@@ -131,7 +131,7 @@ void timeStep ()
 		return;
 
 	// Simulation code
-	for (unsigned int i = 0; i < 4; i++)
+	for (unsigned int i = 0; i < 8; i++)
 		sim.step(model);
 }
 
diff --git a/Demos/Simulation/SimulationModel.cpp b/Demos/Simulation/SimulationModel.cpp
@@ -20,6 +20,8 @@ SimulationModel::SimulationModel()
 	m_solid_poissonRatio = 0.3f;
 	m_solid_normalizeShear = false;
 	m_solid_normalizeStretch = false;
+
+	m_groupsInitialized = false;
 }
 
 SimulationModel::~SimulationModel(void)
@@ -86,6 +88,11 @@ SimulationModel::ConstraintVector & SimulationModel::getConstraints()
 	return m_constraints;
 }
 
+SimulationModel::ConstraintGroupVector & SimulationModel::getConstraintGroups()
+{
+	return m_constraintGroups;
+}
+
 void SimulationModel::updateConstraints()
 {
 	for (unsigned int i = 0; i < m_constraints.size(); i++)
@@ -98,7 +105,10 @@ bool SimulationModel::addBallJoint(const unsigned int rbIndex1, const unsigned i
 	BallJoint *bj = new BallJoint();
 	const bool res = bj->initConstraint(*this, rbIndex1, rbIndex2, pos);
 	if (res)
+	{
 		m_constraints.push_back(bj);
+		m_groupsInitialized = false;
+	}
 	return res;
 }
 
@@ -107,7 +117,10 @@ bool SimulationModel::addBallOnLineJoint(const unsigned int rbIndex1, const unsi
 	BallOnLineJoint *bj = new BallOnLineJoint();
 	const bool res = bj->initConstraint(*this, rbIndex1, rbIndex2, pos, dir);
 	if (res)
+	{
 		m_constraints.push_back(bj);
+		m_groupsInitialized = false;
+	}
 	return res;
 }
 
@@ -116,7 +129,10 @@ bool SimulationModel::addHingeJoint(const unsigned int rbIndex1, const unsigned
 	HingeJoint *hj = new HingeJoint();
 	const bool res = hj->initConstraint(*this, rbIndex1, rbIndex2, pos, axis);
 	if (res)
+	{
 		m_constraints.push_back(hj);
+		m_groupsInitialized = false;
+	}
 	return res;
 }
 
@@ -125,7 +141,10 @@ bool SimulationModel::addUniversalJoint(const unsigned int rbIndex1, const unsig
 	UniversalJoint *uj = new UniversalJoint();
 	const bool res = uj->initConstraint(*this, rbIndex1, rbIndex2, pos, axis1, axis2);
 	if (res)
+	{
 		m_constraints.push_back(uj);
+		m_groupsInitialized = false;
+	}
 	return res;
 }
 
@@ -135,7 +154,10 @@ bool SimulationModel::addRigidBodyParticleBallJoint(const unsigned int rbIndex,
 	RigidBodyParticleBallJoint *bj = new RigidBodyParticleBallJoint();
 	const bool res = bj->initConstraint(*this, rbIndex, particleIndex);
 	if (res)
+	{
 		m_constraints.push_back(bj);
+		m_groupsInitialized = false;
+	}
 	return res;
 }
 
@@ -144,7 +166,10 @@ bool SimulationModel::addDistanceConstraint(const unsigned int particle1, const
 	DistanceConstraint *c = new DistanceConstraint();
 	const bool res = c->initConstraint(*this, particle1, particle2);
 	if (res)
+	{
 		m_constraints.push_back(c);
+		m_groupsInitialized = false;
+	}
 	return res;
 }
 
@@ -154,7 +179,10 @@ bool SimulationModel::addDihedralConstraint(const unsigned int particle1, const
 	DihedralConstraint *c = new DihedralConstraint();
 	const bool res = c->initConstraint(*this, particle1, particle2, particle3, particle4);
 	if (res)
+	{
 		m_constraints.push_back(c);
+		m_groupsInitialized = false;
+	}
 	return res;
 }
 
@@ -164,7 +192,10 @@ bool SimulationModel::addIsometricBendingConstraint(const unsigned int particle1
 	IsometricBendingConstraint *c = new IsometricBendingConstraint();
 	const bool res = c->initConstraint(*this, particle1, particle2, particle3, particle4);
 	if (res)
+	{
 		m_constraints.push_back(c);
+		m_groupsInitialized = false;
+	}
 	return res;
 }
 
@@ -174,7 +205,10 @@ bool SimulationModel::addFEMTriangleConstraint(const unsigned int particle1, con
 	FEMTriangleConstraint *c = new FEMTriangleConstraint();
 	const bool res = c->initConstraint(*this, particle1, particle2, particle3);
 	if (res)
+	{
 		m_constraints.push_back(c);
+		m_groupsInitialized = false;
+	}
 	return res;
 }
 
@@ -184,7 +218,10 @@ bool SimulationModel::addStrainTriangleConstraint(const unsigned int particle1,
 	StrainTriangleConstraint *c = new StrainTriangleConstraint();
 	const bool res = c->initConstraint(*this, particle1, particle2, particle3);
 	if (res)
+	{
 		m_constraints.push_back(c);
+		m_groupsInitialized = false;
+	}
 	return res;
 }
 
@@ -194,7 +231,10 @@ bool SimulationModel::addVolumeConstraint(const unsigned int particle1, const un
 	VolumeConstraint *c = new VolumeConstraint();
 	const bool res = c->initConstraint(*this, particle1, particle2, particle3, particle4);
 	if (res)
+	{
 		m_constraints.push_back(c);
+		m_groupsInitialized = false;
+	}
 	return res;
 }
 
@@ -204,7 +244,10 @@ bool SimulationModel::addFEMTetConstraint(const unsigned int particle1, const un
 	FEMTetConstraint *c = new FEMTetConstraint();
 	const bool res = c->initConstraint(*this, particle1, particle2, particle3, particle4);
 	if (res)
+	{
 		m_constraints.push_back(c);
+		m_groupsInitialized = false;
+	}
 	return res;
 }
 
@@ -214,7 +257,10 @@ bool SimulationModel::addStrainTetConstraint(const unsigned int particle1, const
 	StrainTetConstraint *c = new StrainTetConstraint();
 	const bool res = c->initConstraint(*this, particle1, particle2, particle3, particle4);
 	if (res)
+	{
 		m_constraints.push_back(c);
+		m_groupsInitialized = false;
+	}
 	return res;
 }
 
@@ -223,7 +269,10 @@ bool SimulationModel::addShapeMatchingConstraint(const unsigned int numberOfPart
 	ShapeMatchingConstraint *c = new ShapeMatchingConstraint(numberOfParticles);
 	const bool res = c->initConstraint(*this, particleIndices, numClusters);
 	if (res)
+	{
 		m_constraints.push_back(c);
+		m_groupsInitialized = false;
+	}
 	return res;
 }
 
@@ -270,3 +319,66 @@ void SimulationModel::addTetModel(
 }
 
 
+
+void SimulationModel::initConstraintGroups()
+{
+	if (m_groupsInitialized)
+		return;
+
+	const unsigned int numConstraints = (unsigned int) m_constraints.size();
+	const unsigned int numParticles = (unsigned int) m_particles.size();
+	const unsigned int numRigidBodies = (unsigned int) m_rigidBodies.size();
+	const unsigned int numBodies = numParticles + numRigidBodies;
+	m_constraintGroups.clear();
+
+	// Maps in which group a particle is or 0 if not yet mapped
+	std::vector<unsigned char*> mapping;
+
+	for (unsigned int i = 0; i < numConstraints; i++)
+	{
+		Constraint *constraint = m_constraints[i];
+
+		bool addToNewGroup = true;
+		for (unsigned int j = 0; j < m_constraintGroups.size(); j++)
+		{
+			bool addToThisGroup = true;
+
+			for (unsigned int k = 0; k < constraint->m_numberOfBodies; k++)
+			{
+				if (mapping[j][constraint->m_bodies[k]] != 0)
+				{
+					addToThisGroup = false;
+					break;
+				}
+			}
+
+			if (addToThisGroup)
+			{
+				m_constraintGroups[j].push_back(i);
+
+				for (unsigned int k = 0; k < constraint->m_numberOfBodies; k++)
+					mapping[j][constraint->m_bodies[k]] = 1;
+
+				addToNewGroup = false;
+				break;
+			}
+		}
+		if (addToNewGroup)
+		{
+			mapping.push_back(new unsigned char[numBodies]);
+			memset(mapping[mapping.size() - 1], 0, sizeof(unsigned char)*numBodies);
+			m_constraintGroups.resize(m_constraintGroups.size() + 1);
+			m_constraintGroups[m_constraintGroups.size()-1].push_back(i);
+			for (unsigned int k = 0; k < constraint->m_numberOfBodies; k++)
+				mapping[m_constraintGroups.size() - 1][constraint->m_bodies[k]] = 1;
+		}
+	}
+
+	for (unsigned int i = 0; i < mapping.size(); i++)
+	{
+		delete[] mapping[i];
+	}
+	mapping.clear();
+
+	m_groupsInitialized = true;
+}
diff --git a/Demos/Simulation/SimulationModel.h b/Demos/Simulation/SimulationModel.h
@@ -23,6 +23,8 @@ namespace PBD
 			typedef std::vector<RigidBody*, Eigen::aligned_allocator<RigidBody*> > RigidBodyVector;
 			typedef std::vector<TriangleModel*, Eigen::aligned_allocator<TriangleModel*> > TriangleModelVector;
 			typedef std::vector<TetModel*, Eigen::aligned_allocator<TetModel*> > TetModelVector;
+			typedef std::vector<unsigned int> ConstraintGroup;
+			typedef std::vector<ConstraintGroup> ConstraintGroupVector;
 
 
 		protected:
@@ -31,6 +33,7 @@ namespace PBD
 			TetModelVector m_tetModels;
 			ParticleData m_particles;
 			ConstraintVector m_constraints;
+			ConstraintGroupVector m_constraintGroups;
 
 			float m_cloth_stiffness;
 			float m_cloth_bendingStiffness;
@@ -57,6 +60,8 @@ namespace PBD
 			TriangleModelVector &getTriangleModels();
 			TetModelVector &getTetModels();
 			ConstraintVector &getConstraints();
+			ConstraintGroupVector &getConstraintGroups();
+			bool m_groupsInitialized;
 
 			void addTriangleModel(
 				const unsigned int nPoints,
@@ -73,6 +78,7 @@ namespace PBD
 				unsigned int* indices);
 
 			void updateConstraints();
+			void initConstraintGroups();
 
 			bool addBallJoint(const unsigned int rbIndex1, const unsigned int rbIndex2, const Eigen::Vector3f &pos);
 			bool addBallOnLineJoint(const unsigned int rbIndex1, const unsigned int rbIndex2, const Eigen::Vector3f &pos, const Eigen::Vector3f &dir);
diff --git a/Demos/Simulation/TimeStepController.cpp b/Demos/Simulation/TimeStepController.cpp
diff --git a/Doxyfile b/Doxyfile
diff --git a/README.md b/README.md

Original file line number	Diff line number	Diff line change
`@@ -1,3 +1,4 @@`
	`1`	`+ - parallelized unified solver using graph coloring`
`1`	`2`	`- implemented unified solver for rigid bodies and deformable solids`
`2`	`3`	`- added documentation for generic constraints`
`3`	`4`	`- added generic constraint demo`
Original file line number	Diff line number	Diff line change
`@@ -147,7 +147,7 @@ void timeStep ()`
`147`	`147`	`return;`
`148`	`148`
`149`	`149`	`// Simulation code`
`150`		`- for (unsigned int i = 0; i < 4; i++)`
	`150`	`+ for (unsigned int i = 0; i < 8; i++)`
`151`	`151`	`simulation.step(model);`
`152`	`152`	`}`
`153`	`153`
Original file line number	Diff line number	Diff line change
`@@ -130,7 +130,7 @@ void timeStep ()`
`130`	`130`	`return;`
`131`	`131`
`132`	`132`	`// Simulation code`
`133`		`- for (unsigned int i = 0; i < 4; i++)`
	`133`	`+ for (unsigned int i = 0; i < 8; i++)`
`134`	`134`	`sim.step(model);`
`135`	`135`	`}`
`136`	`136`
Original file line number	Diff line number	Diff line change
`@@ -131,7 +131,7 @@ void timeStep ()`
`131`	`131`	`return;`
`132`	`132`
`133`	`133`	`// Simulation code`
`134`		`- for (unsigned int i = 0; i < 4; i++)`
	`134`	`+ for (unsigned int i = 0; i < 8; i++)`
`135`	`135`	`sim.step(model);`
`136`	`136`	`}`
`137`	`137`