Rigid Body class design

This topic is 3613 days old which is more than the 365 day threshold we allow for new replies. Please post a new topic.

Recommended Posts

Well, I was not sure where to create this topic, as it is related to both general programming( design issues ) and to Math and Physics. I am just expecting that the experience of some of you in designing classes related to Physics Simulation will help me. So, I've designed two RigidBody classes, each using different state variables: the two have position and orientation( Matrix or quaternion ) as state variables, the first one, linear and angular momentums in addition and the second one, linear and angular velocities...I designed the second class because computing forces and torques to have a rigid body, after updating, in a certain state( velocities ) was much more intuitive. When designing both classes, I notice that in fact, designing a function to compute the state derivate and update the auxiliary variables( for the momentum approach, these are the velocities( linear and angular ) and the inverse of the world inertia tensor ), made the code tighter and eased the debugging process. For the state updating, I do it outside the class: to avoid hard coding the integration scheme inside the class, and hence enabling me to use any integration method. I've defined some explicit Runge-Kutta methods, and it was easy to change the integration scheme during a simulation. To achieve that, I wrote a function template inside my RigidBody class, with the integration method as a template parameter. Well, to make the RigidBody class, more general, the idea came to me, to create a template RigidBody class, with a StateVariables class and an AuxiliaryVariables class as template parameters, and to make it inherit from both class templates...I avoid using composition, since not knowing state variables and auxiliary variables in advance, in the general case, forces me to make the RigidBody class inherit from both classes. I simply define inside the RigidBody class, constructors and an update function( template ) as following:
template< class Properties,
class StateVariables,
class AuxiliaryVariables >

class RigidBody: public Properties, public StateVariables, public AuxiliaryVariables
{
public:
RigidBody(){}
RigidBody( const Properties& props,
const StateVariables& state,
const AuxiliaryVariables& auxiliary );

template< class IntegrationScheme,
class ForceGenerator,
class TorqueGenerator,
typename Real >

void update( IntegrationScheme& scheme,
ForceGenerator& force,
TorqueGenerator& torque,
const Real& time,
const Real& deltaTime )
{
scheme.integrate( force, torque, *this, *this, *this, time, deltaTime );
}
////////////////////////////////////////////////////////////////////////////
//functions?
};



template<class Vector,class Matrix> class Auxiliary;
template<typename Real,class Matrix> class Properties;
class ForceGenerator;
class TorqueGenerator;

template<class Vector,class Matrix>
class State
{
private:
Vector position_;
Matrix orientM_;
Vector linMom_;
Vector angMom_;
////////////////////////////////////////////////////////////////////////////////
public:
////////////////////////////////////////////////////////////////////////////////
//constructors...
////////////////////////////////////////////////////////////////////////////////
//accessors
////////////////////////////////////////////////////////////////////////////////
//utility functions
////////////////////////////////////////////////////////////////////////////////
template< class ForceGenerator, class TorqueGenerator, typename Real >
void computeDerivate( const ForceGenerator&,
const TorqueGenerator&,
const Properties&,
Auxiliary&,
const Real&,
State<Vector,Matrix>& stateDerivate )const;
void updateAuxiliary( Properties&,
Auxiliary& )const;
};


template<class Vector,class Matrix>
class Auxiliary
{
public:
//constructors, destructors
////////////////////////////////////////////////////////////////////////////////
//accessors
////////////////////////////////////////////////////////////////////////////////
private:
Vector angVel_,
Vector linVel_,
Matrix invInertia_;
};


template<class Real,class Matrix>
class Properties
{
public:
//constructors,destructors
////////////////////////////////////////////////////////////////////////////////
//accessors
////////////////////////////////////////////////////////////////////////////////
private:
Real invMass_;
Matrix bodyInertia_;
Matrix invBodyInertia_;
};


class IntegrationScheme
{
public:
IntegrationScheme();

template< class ForceGenerator,
class TorqueGenerator,
class Properties,
class State,
class Auxiliary,
typename Real>

void integrate( const ForceGenerator&,
const TorqueGenerator&,
const Properties&,
State&,
Auxiliary&,
const Real&,
const Real& );
};


I just want to have your opinions, on a design point of view... Thanks for any remark. [Edited by - johnstanp on July 23, 2008 1:03:05 PM]

Share on other sites
As you may have guessed, the State, Auxiliary and Properties classes are simple implementation examples, to illustrate the idea...

Share on other sites
I'd rather not reply with 'advice' as I can hardly say I'm an expert in this subject, but I have written a physics thing or two. Instead of giving advice, I'm just going to post my rigid body definitions. After a somewhat cursory overview of your code it seems that we've approached things similarly, but mine has a lot of extra 'stuff' added in. I was reading GamePhysics by Eberly when learning how to write mine.

The PredictNewPositionandOrientation, force_and_torque and Convert functions all relate to solving the equations of motion. Nearly every object (in my code) uses an implementation of RK4 called from 'PredictNewPositionandOrientation(),' but simpler objects can use a simple Euler integration. The force-and-torque virtual method (I see you have force and torque 'generators') is used to calculate the intermediate steps of a diff. eqn. solver. These intermediate steps may also be referred to as 'curvature' functions (as forces and torques are basically derivatives of the equations of motion).

This setup is to make it easier to write implementations of vastly different objects, e.g. Hovercrafts, air planes and cars (so far I've just written cheap ways to simulate air cushion vehicles and cars with suspension).

Similarly with you I keep my RK4 implementation (and other solvers) outside of the rigid body class.

#ifndef	RIGID_BODY_H#define	RIGID_BODY_H#include	"Vector.h"#include	"Matrix.h"#include	"Quaternion.h"struct	StateVector{	Vector3	mPosition;	Vector3	mLinearMomentum;	Vector3	mLinearVelocity;	Vector3	mAngularMomentum;	Vector3	mAngularVelocity;	Quaternion	mQuat_Orientation;	Matrix3x3	mMat_Orientation;	Vector3		mEuler_Orientation;};class	RigidBody{	public:	RigidBody()		{		mMass = 1; mMassInv = 1;	mRadius	=	1.0f;		mStiffness			=	1;		mStaticFriction		=	1;		mDynamicFriction	=	1;	}	~RigidBody()	{	}			virtual	void	PredictNewPosAndOrientation(double	DT) = 0;	virtual void	Force_and_Torque(Vector3 *force,Vector3 *torque,Vector3	& position, Quaternion & quat_orientation,Vector3 & linear_momentum,		Vector3 & angular_momentum,Matrix3x3 & mat_orientation,Vector3 & linear_velocity,Vector3 & angular_velocity) = 0;	virtual void	TouchFunc(RigidBody *pOther, Vector3 & Normal, Vector3 & HitPoint)=0;	//if NULL, then we hit the world	void	SetInertia(float Ixx,float	Iyy,float	Izz)	{		mJBody.mat[0][0]	=	Ixx;		mJBody.mat[1][1]	=	Iyy;		mJBody.mat[2][2]	=	Izz;		if(Ixx != 0)			mJBodyInverse.mat[0][0]	=	1.0f/Ixx;		else			mJBodyInverse.mat[0][0]	=	0.0f;		if(Iyy != 0)			mJBodyInverse.mat[1][1]	=	1.0f/Iyy;		else			mJBodyInverse.mat[1][1]	=	0.0f;		if(Izz != 0)			mJBodyInverse.mat[2][2]	=	1.0f/Izz;		else			mJBodyInverse.mat[2][2]	=	0.0f;	}	//p = linear momentum, l = angular momentum	//R = output matrix orientation	//V = output linear velocity	//W = output angular velocity	//The first three parameters converted to last 3, respectively	void	Convert(Quaternion & Q, Vector3	& P,Vector3	&L,Matrix3x3 & R,Vector3 & V,Vector3 & W);		void	ApplyImpulse(Vector3	&	dMomentum,Vector3	&	ContactPoint);	//Applies LINEAR impulse only	void	ApplyLinearImpulse(Vector3	&	dMomentum)	{		CurrentState.mLinearMomentum	+=	dMomentum;	}	//Applies ANGULAR impulse only	void	ApplyAngularImpulse(Vector3	&	dMomentum,Vector3 &	ContactPoint)	{		CurrentState.mAngularMomentum += CrossProduct(&dMomentum,&ContactPoint);	}	void	AddForce(Vector3	&	Force,Vector3	& Location)	//note that location must be in local coordinates	{		mForces		+=	Force;		mMoments	+=	CrossProduct(&Location,&Force);	}	float	GetMass() const { return	mMass; }	float	GetMassInv() const { return	mMassInv;	}	void	SetMass(float	IN);	/*		Constant physical properties		*/	float		mRadius;	//Radius that encloses all collision hulls of the body	Matrix3x3	mJBody;			//Inertia tensor in body coordinates	Matrix3x3	mJBodyInverse;	//Inverse inertia tensor in body coordinates (never changes)	float		mStiffness;		//Used to calculate coefficient of restitution	float		mStaticFriction;	float		mDynamicFriction;		/**		DYNAMIC PHYSICAL PROPERTIES	**/	StateVector	CurrentState;	/*		Predicted states, filled in by the built-in differential equation solver	*/	StateVector	PredictedStateVector;		Matrix3x3	mJGlobal;//Moment of inertia in global coordinates 						 //must update with orientation to be accurate, not really needed	Matrix3x3	mJInverseGlobal;	//This gets computed frequently, inverse mass moment of inertia	//I made these private because changing one always requires a change in the other	protected:	float		mMass;			//Mass in KG	float		mMassInv;			void	RK4Update(double);};StateVector GlobalRK4Update(StateVector &in,double DT);/*	y(t + dt) = y(x) + dt*y'(x) + [dt^2 / 2!]y''(x) + [dt^3/3!]y'''(x)	Linear:	y(x) is current position	y'(x) is velocity	y''(x) is acceleration , or Force/Mass or Force * iMass*//*	Need to be able to:	Convert Vector3D to a quaternion	(implemented this, haven't tested it)	Convert	quaternion to an orientation matrix	(implemented this, don't know if it works)	Transpose a 3x3 matrix (just swap rows and columns, should be easy)	(implemented this, tested it, seems to work fine)*//*  *A*:	To compute global inverse mass moment of inertia tensor:	-Integrate quaternion orientation	(see note B)	-Turn the quaternion orientation into a matrix orientation	(R, a 3x3 matrix)	-The inverse moment of inertia matrix in global space is:		R * mInvInertiaTensor * R^T	where R^T is the transpose of R	-The angular velocity is the inverse moment of inertia matrix multiplied by the angular		momentum vector    *B*:	The change in quaternion orientation with respect to time:	dq(t)	/	dt	=	(.5)(w(t))(q(t))		where w(t) is the quaternion which represents the angular velocity*/#endif

Share on other sites
I approach this problem differently than either of you. My Rigid Object class contains almost no physics logic code and is almost entirely data and accessors for that data. It doesn't even know that it is being integrated. I don't have near the level of customization that johnstanp's design has, but mine is designed to not be the ultimate physics object class (i've gone down that road before), and is focused on simplicity, no coupling to other objects, and a data-driven design.

A rigid object contains the following data members:
• scale - a scalar value which dynamically and uniformly scales all 3 dimensions of the object's shape (default is 1).
• center of mass - setable by the user as object-space or world space.
• position - the user's position on the object (usually for graphics), to which the shape is relatively positioned. This position and the center or mass are kept synchronized whenever the state of the object's transformation is changed.
• orientation - a 3x3 orthonormal matrix representing the rotation from world to object space.
• velocity - world space velocity vector.
• bias velocity - a temporary velocity zeroed out every frame and used in the constraint solvers to correct errors. (I need to get rid of this and find a better way that is less hack-ish)
• angular momentum - world space angular momentum vector. I should probably change this to angular velocity for consistency. It doesn't really matter which one, the constraint solvers require as many computations either way.
• bias angular momentum - same as bias velocity.
• force - world space mass-acceleration of the object.
• torque - world space torque vector of the object.
• mass, inverse mass - setable by the user, which can retrieve this data from the collision shape if desired.
• inverse body space inertia tensor - same semantics as mass.
• inverse world space inertia tensor - read-only externally and updated internally from the body space tensor whenever orientation is changed.
• a pointer to the object's deactivation state- the object is created internally and is NULL by default.
• a pointer to the object's collision shape instance - similar to the deactivation state, created from a collision shape template object and then given the object's transformation whenever the shape is retrieved from the object.
• a boolean for whether or not the object is static.

The reason why I opted for a simple, hard-coded design is that there is never really any need in physical simulations for exotic object types. I considered using a component-based system for the object, but there was no way that I could get it to work without seriously sacrificing performance. In order to allow objects to have extensions such as a deactivation, I have instead given them pointers to hard-coded object types which are NULL by default. This way, an object can have a shape or deactivation state, but doesn't have to if the resources are not needed or available. To give an object a deactivation state, one calls object->addDeactivationState().

All of the logic for collision detection/deactivation/integration is handled by subclasses of interfaces (CollisionWorld, Deactivator, Integrator), which allows the user to provide different implementations for integration or collision detection, etc if they want. The user adds the desired objects to the different physics systems and they will then operate on the objects if those objects have all of the data required for those systems. For instance, if an object doesn't have a collision shape yet is added to the collision world, no collision detection will be performed.

This is the 3rd or 4th revision of my physics engine which I've been working on for 1.5 years. I can discuss it in more detail if desired, especially about the collision detection and constraint subsystems.

Share on other sites
mJBody.mat[0][0]	=	Ixx;mJBody.mat[1][1]	=	Iyy;mJBody.mat[2][2]	=	Izz;if(Ixx != 0)	mJBodyInverse.mat[0][0]	=	1.0f/Ixx;else	mJBodyInverse.mat[0][0]	=	0.0f;if(Iyy != 0)	mJBodyInverse.mat[1][1]	=	1.0f/Iyy;else	mJBodyInverse.mat[1][1]	=	0.0f;if(Izz != 0)	mJBodyInverse.mat[2][2]	=	1.0f/Izz;else	mJBodyInverse.mat[2][2]	=	0.0f;

Funny code. You can leave out the else-statements completely or

mJBody.mat[0][0] = Ixx ? 1.0f/Ixx : 0.0f;mJBody.mat[1][1] = Iyy ? 1.0f/Iyy : 0.0f;mJBody.mat[2][2] = Izz ? 1.0f/Izz : 0.0f;

Share on other sites
This is way too complicated, over engineered and abstract in my opinion. Make your rigid body class as *simple* as possible with some helper functions, like e.g. transforming points and vectors into the local body space and vice versa, applying forces and impulses. There is really no need for ll these indirections.

Note that you also need some collision shape to identify contacts between shapes.

Share on other sites
Quote:
 Original post by DonDickieDThis is way too complicated, over engineered and abstract in my opinion. Make your rigid body class as *simple* as possible with some helper functions, like e.g. transforming points and vectors into the local body space and vice versa, applying forces and impulses. There is really no need for ll these indirections.Note that you also need some collision shape to identify contacts between shapes.

I agree.

Share on other sites
Quote:
 This is way too complicated, over engineered and abstract in my opinion. Make your rigid body class as *simple* as possible with some helper functions, like e.g. transforming points and vectors into the local body space and vice versa, applying forces and impulses. There is really no need for ll these indirections.

Interesting. I'd like to see your Rigid body definitions, because I actually think of mine as somewhat over-simplified! It also seems that the differences between Agessera's definitions and mine are ultimately pretty small (considering how vastly different they could have been). The only functions I have that might be considered 'physics logic' only has to do with differential equation solving, and they're ultimately just place holders at that!

Share on other sites
I was refering to the original design. I see no problems with Aressera implementation.

Note that your are not solving ODEs in *constrained* rigid body dynamics, but DAEs. As long as you don't solve for the constraint forces explicetely (like e.g. Baraff) you can really forget about integrators for ODEs.

Simplified the problem can be written as:

dx/dt = v
dv/dt = M^-1 * ( f_ext + JT * lambda )

C(x) = 0

Share on other sites
Again, I'd like to see your equivalent definitions. That'd probably help give us ideas.

• 48
• 12
• 10
• 10
• 9
• Forum Statistics

• Total Topics
631380
• Total Posts
2999673
×