Best way to make a physics object hover

Started by
5 comments, last by oliii 17 years, 7 months ago
Hi, If I'm using a physics engine whats the best way to make an object hover? I've looked all over google and its not much help - just abstracted heavily mathematical theory on hovercrafts which isn't useful and such. The approach I've tried is to get the height and at every step, apply an upthrust obeying the inverse square law, assuming that the height will stabilise out. The problem is, when it hits the ground it will suddenly jump high up (K/0.01*0.01 is a large number), so I used clipping. Then another problem is that there won't be enough force above a certain height, so it will float gently down then jump up, when the floors upthrust tips the critical force to outweigh the gravity. Then I tried switching off gravity - no luck. I also tried more force when its falling, which seems a mediocre fix- the object will jut around at a fixed height as it re-adjusts itself. I wouldn't mind if it slowly cycled its height around a certain height, but its the slow falling and sudden jump upward as opposed to moving like a sine wave around its height. I've tried all kinds of laws like exponential, linear .etc And the working solution is very hackish- can anyone impart any advice? If it helps I'm using ODE. Thanks for your time.
Advertisement
In these kind fo systems (let it be A.I. sensors, or how much a A.I. car should apply force to correct it's steering), I usually look 'forward' in time, instead of the instantaneous inputs.

Say you use a sort of spring, like you do.

F = k / l(t)^2

instead of looking at l(t), I try to find out what the length will be after say, one second, or 1/2 second, or whatever looks good.

So if the craft is falling looking ahead in time, then you will have to apply a strong force, and has the equilibrium point is reached, the force will counteract the movement gently and naturally.

It's kind of a damping system if you will. There could be a slight sway, but as you said, that's not a problem, and that's tweakable depending on how much you look ahead in time. It's kind of a negative feedback system.

Everything is better with Metal.

Hmm, well, I've implemented a hovercraft (air cushion vehicle) simulation, but you say you aren't interested in the heavy math aspect of it.

Unfortunately, for anything that is realtime, interactive, and true 3D, you are going to have an extremely difficult time implementing this without the difficult theory!

I honestly don't believe that looking ahead in time, or any other time of coherence implementation, will work, except for in extremely idealized situations (i.e. sometime that is limited in what you can do with the hovercraft, or something that isn't particularly interactive).

I implement the air cushion by defining a certain number of sample points arranged in a circular fashion underneath the chassis of the hovercraft (the skirt and plenum chamber area). I run an RK4 solver to predict the new position and orientation, based on the inputs from the air cushion. Each sample point behaves like a spring, but a more detailed fluid flow simulation could be used instead of basic distance finding (which would require more sample points; finite element analysis).

Here's the force and torque function (the 'derivative' function when reading about numerical methods of problem solving). For each sample point, it measures the distance to the ground and applies a force at that point based on the relative velocity between the sample point and the ground, and a dampening factor is also added in there. The force per sample point is the total weight of the hovercraft, divided by the number of thrusters, at a given height (like a spring).

Here's this function:
void	Hovertank::Force_and_Torque(Vector3 *output_force,Vector3 *output_torque,Vector3 & position,Quaternion & quat_orientation,									Vector3 & linear_momentum,Vector3 & angular_momentum,Matrix3x3 & mat_orientation,									Vector3 & linear_velocity,Vector3 & angular_velocity){	*output_force  += Vector3(0.0f,-gpPhysicsManager->mGravAccel*mMass,0.0f);	//applies weight		Air_Cushion	transformed_cushion	=	mAirCushion;	Matrix4x4	cushion_sample_transformation;	cushion_sample_transformation.SetTranslation(&position);	cushion_sample_transformation = mat_orientation;//CurrentState.mMat_Orientation;		for(int i = 0; i < transformed_cushion.num_sample_points; i++)	{		transformed_cushion = cushion_sample_transformation * transformed_cushion;	}	transformed_cushion.orientation	=	mat_orientation * transformed_cushion.orientation;		Vector3	start = position;	Vector3	end = position + (transformed_cushion.orientation * transformed_cushion.rest_height);#if!	PRE_FIND_TOUCHED_LEAFS	std::vector<int>	leafs;	MoveTouchNode(gpBSP,0,&start,&end,mRadius,leafs);	//collect leafs#endif	Vector3	total_cushion_force(0.0f,0.0f,0.0f);			Vector3	total_cushion_moment(0.0f,0.0f,0.0f);		for(i = 0; i < transformed_cushion.num_sample_points; i++)	{		end	=	transformed_cushion + (transformed_cushion.orientation * 			(transformed_cushion.rest_height));	//FIXME: transformed_cushion.rest_height - radius ???		Vector3	pos = transformed_cushion;			HitInfo	CGTrace	=	#if	PRE_FIND_TOUCHED_LEAFS			TraceSphereToBSPFaces(gpBSP,transformed_cushion,end,transformed_cushion.sample_radius,gpPhysicsManager->physics_epsilon,mvBSPLeafs);#else			TraceSphereToBSPFaces(gpBSP,transformed_cushion,end,transformed_cushion.sample_radius,gpPhysicsManager->physics_epsilon,leafs);#endif		if(CGTrace.Time < 1.0f)		{			Vector3	local_pos	=	transformed_cushion - position;	//reflects rotation only#if	USE_REAL_DIST			float	dist_to_point	=	(CGTrace.HitPoint - transformed_cushion).GetLength();#else			float	dist_to_point	=	(DotProduct(&CGTrace.HitPoint,&transformed_cushion.orientation) - 										DotProduct(&transformed_cushion,&transformed_cushion.orientation));#endif			float	compression		=	fabs(dist_to_point - transformed_cushion.rest_height);					float	basic_force_magnitude		=	compression * transformed_cushion.weight_per_thruster;				if(basic_force_magnitude > transformed_cushion.max_sample_force)				basic_force_magnitude = transformed_cushion.max_sample_force;			Vector3	sample_force = transformed_cushion.orientation * basic_force_magnitude * -1.0f;			float	total_rel_vel(0.0f),rel_vel_fromlinear(0.0f),rel_vel_fromangular(0.0f);			rel_vel_fromlinear	=	DotProduct(&transformed_cushion.orientation,&linear_velocity);			rel_vel_fromangular	=	DotProduct(&transformed_cushion.orientation,				&CrossProduct(&angular_velocity,&local_pos));			total_rel_vel	=	rel_vel_fromlinear + rel_vel_fromangular;						sample_force	+=	transformed_cushion.orientation * total_rel_vel * -transformed_cushion.dampen_coefficient;			total_cushion_force	+=	sample_force;			total_cushion_moment += CrossProduct(&local_pos,&sample_force);		}	}		*output_force	+=	total_cushion_force;	*output_torque	+=	total_cushion_moment;}


I dont' know how helpful this is, but seeing as how I've already implemented this it may help to give some ideas.


...and do not wildly extrapolate. Just because Saddam Hussein gassed Kurds in 1990 doesn't mean he eats babies' brains.
Quote:Original post by oliii
In these kind fo systems (let it be A.I. sensors, or how much a A.I. car should apply force to correct it's steering), I usually look 'forward' in time, instead of the instantaneous inputs.


What you're describing is known as a "PID controller". Genjix, Googling for that term should turn up plenty of info on what is a very common problem.
tthibault:

Didn't say I was interested in heavy math. Just real-life math I cannot simulate ;)

Yours is a very interesting code which I'll have to check out.

Sneftel and oliii:

Looking from the tutorials this looks great. I found this on the net:
typedef struct{  double dState;      	// Last position input  double iState;      	// Integrator state  double iMax, iMin;  	  // Maximum and minimum allowable integrator state  double	iGain,    	// integral gain        	pGain,    	// proportional gain         	dGain;     	// derivative gain} SPid;double UpdatePID(SPid * pid, double error, double position){  double pTerm, dTerm, iTerm;  pTerm = pid->pGain * error;     // calculate the proportional term// calculate the integral state with appropriate limiting  pid->iState += error;  if (pid->iState > pid->iMax) pid->iState = pid->iMax;  else if (pid->iState < pid->iMin) pid->iState = pid->iMin;  iTerm = pid->iGain * iState;  // calculate the integral term  dTerm = pid->dGain * (position - pid->dState);  pid->dState = position;  return pTerm + iTerm - dTerm;}


I assume in my context error would mean the difference in height between my position and the desired position?

Would this roughly translate to:
float celPcHover::UpthrustForce (){  float pval, dval, ival;  float error = fabs (WorldToLocal (position - end_pos). Length ())  // calculate the proportional term  pval = error * prop_constant;  // calculate the integral term  istate += error;  ival = istate * int_constant;  // calculate the differential term  dval = pcmechobj->GetLinearVelocity () * vel_constant;  return pval + ival + dval;}


Am I correct in my reasoning? Would the integrator part be the total amount of excess above the desired value over the past?

Thanks.
Quote:
Why do you do this?

cushion_sample_transformation.SetTranslation(&position);

cushion_sample_transformation = mat_orientation;//CurrentState.mMat_Orientation;


Then you loop through all the sample points and transform them by that matrix. And then you again transform the cushion itself (why?)- wouldn't that be a double transform on each of the points?


I've defined the air cushion to have the following properties:
-a bunch of sample points laid out circularly *in object space*.
-An orientation (which points downward).

I realize now there are a few things I've left out which might confuse you, namely that the index operator of the air cushion structure returns the position of the indexed sample point, e.g.:

air_cushion returns the position of sample point at index i.

In order to use the air cushion sample points they must be rotated out to reflect the translation and orientaiton of the hovercraft ( if the hovercraft starts to bank at an angle of 30 degrees, the position of the sample points rotates along with the chassis of the hovercraft). The direction that the sample point fires in is simply rotated.

>>Is weight_per_thruster = weight / n_thruster ?

Yes

>>Then it looks like you clip the force. Whats sample_force ?

Sample force is the force contribution of the current sample force.

Quote:
Then you subtract sample_force from our earlier calculated force. Why is this?

transformed_cushion.orientation * total_rel_vel * - transformed_cushion.dampen_coefficient


Hmm, I'm not entirely sure what you mean by the first part of this question, but the second line is how you model a dampened spring. According to hook's law, a spring exerts a force somewhat proportional to a constant times the compression length (in this case, how far the sample point is being compressed from the rest height). In order to dampen the spring so that the system loses energy and can come to a halt, an added force must be added that acts against the relative velocity (hence the negative). dampen_coefficient is basically just an ad hoc number that works. The larger the dampen coefficient is, the more quickly the system losts energy...the lower the value, the more the hovercraft will 'bounce' (which looks cool).

>>Sorry for all the questions, but this is the first promising thing I've found for weeks. Are you using ODE per chan

Don't apologize for asking questions!!! It took me about a month to get this to actually work. I don't use ODE, I try to write everything by myself (but what I've implemented is far less than any physics library).

Note that I implemented the force_and_torque function, but was using it with an euler solver, and therefore the logic was right but the hovercraft would eventually topple over due to numerical instability. That's what actually took me the longest to get working properly. Here's the RK4 solver that calls force_and_torque. This is essentially taken from Game Physics by David Eberly. Note that force_and_torque is a virtual function defined by the rigid body, and so you can write a new force_and_torque function to simulate cars, boats, airplanes, etc without having to change the structure of the RK4Update function.

If you approach the theory in the correct direction an RK4 solver is actually quite intuitive, but the end result certainly doesn't look it! You should start off with reading about an RK2 solver (which is essentially an improved euler method) and ignore matrices of partial derivatives (they only come into play when you take into account the initial time at the frame, I simply take into account the change in time).

//Only need to normalize quaternions at the last stepvoid	RigidBody::RK4Update(double	DT){	//part1	double	halfdt = 0.5f * DT, sixthdt = DT / 6.0f;		Vector3	XN, PN, LN, VN, WN;	Quaternion QN;	Matrix3x3 RN;	/*	//Step1	*/	Vector3	A1DXDT = CurrentState.mLinearVelocity;	Quaternion A1DQDT = (Quaternion(CurrentState.mAngularVelocity) * CurrentState.mQuat_Orientation) * 0.5f;		Vector3	A1DPDT(0.0f,0.0f,0.0f);	Vector3	A1DLDT(0.0f,0.0f,0.0f);	Force_and_Torque(&A1DPDT,&A1DLDT,CurrentState.mPosition,CurrentState.mQuat_Orientation,CurrentState.mLinearMomentum,		CurrentState.mAngularMomentum,CurrentState.mMat_Orientation,CurrentState.mLinearVelocity,CurrentState.mAngularVelocity);		XN = CurrentState.mPosition + (A1DXDT * halfdt);	QN = CurrentState.mQuat_Orientation + (A1DQDT * halfdt);		PN = CurrentState.mLinearMomentum + (A1DPDT * halfdt);	LN = CurrentState.mAngularMomentum + (A1DLDT * halfdt);	Convert(QN,PN,LN,RN,VN,WN);	/*	//step 2	*/	Vector3 A2DXDT = VN;	Quaternion A2DQDT = (Quaternion(WN) * QN) * .5f;		Vector3 A2DPDT(0.0f,0.0f,0.0f);	Vector3 A2DLDT(0.0f,0.0f,0.0f);	Force_and_Torque(&A2DPDT,&A2DLDT,XN,QN,PN,LN,RN,VN,WN);	XN = CurrentState.mPosition + (A2DXDT * halfdt);	QN = CurrentState.mQuat_Orientation + (A2DQDT * halfdt);		PN = CurrentState.mLinearMomentum + (A2DPDT * halfdt);	LN = CurrentState.mAngularMomentum + (A2DLDT * halfdt);	Convert(QN,PN,LN,RN,VN,WN);	/*	//Step 3	*/	Vector3	A3DXDT = VN;	Quaternion A3DQDT = (Quaternion(WN) * QN) * .5f;	Vector3 A3DPDT(0.0f,0.0f,0.0f);	Vector3 A3DLDT(0.0f,0.0f,0.0f);	Force_and_Torque(&A3DPDT,&A3DLDT,XN,QN,PN,LN,RN,VN,WN);	XN = CurrentState.mPosition + (A3DXDT * DT);	QN = CurrentState.mQuat_Orientation + (A3DQDT * DT);		PN = CurrentState.mLinearMomentum + (A3DPDT * DT);	LN = CurrentState.mAngularMomentum + (A3DLDT * DT);	Convert(QN,PN,LN,RN,VN,WN);	/*	//Step 4	*/	Vector3	A4DXDT = VN;	Quaternion A4DQDT = (Quaternion(WN) * QN) * 0.5f;	Vector3 A4DPDT,A4DLDT;	PredictedStateVector.mPosition = CurrentState.mPosition + (A1DXDT + ((A2DXDT + A3DXDT) * 2.0f) + A4DXDT) * sixthdt;	PredictedStateVector.mQuat_Orientation = CurrentState.mQuat_Orientation + (A1DQDT + ((A2DQDT + A3DQDT) * 2.0f) + A4DQDT) * sixthdt;	PredictedStateVector.mQuat_Orientation.Normalize();	CurrentState.mLinearMomentum   = CurrentState.mLinearMomentum + (A1DPDT + ((A2DPDT + A3DPDT) * 2.0f) + A4DPDT) * sixthdt;	CurrentState.mAngularMomentum  = CurrentState.mAngularMomentum + (A1DLDT + ((A2DLDT + A3DLDT) * 2.0f) + A4DLDT) * sixthdt;	Convert(PredictedStateVector.mQuat_Orientation,CurrentState.mLinearMomentum,CurrentState.mAngularMomentum,		PredictedStateVector.mMat_Orientation,CurrentState.mLinearVelocity,CurrentState.mAngularVelocity);	PredictedStateVector.mMat_Orientation = PredictedStateVector.mQuat_Orientation.ToMatrix();	PredictedStateVector.mMat_Orientation.ToRadians(&PredictedStateVector.mEuler_Orientation.x,		&PredictedStateVector.mEuler_Orientation.y,&PredictedStateVector.mEuler_Orientation.z);		mJInverseGlobal	=	CurrentState.mMat_Orientation * mJBodyInverse * CurrentState.mMat_Orientation.Transpose();		//FIXME: this only works when the hulls are always spheres!	CurrentState.mEuler_Orientation  =  PredictedStateVector.mEuler_Orientation;	CurrentState.mMat_Orientation	  = PredictedStateVector.mMat_Orientation;	CurrentState.mQuat_Orientation	  = PredictedStateVector.mQuat_Orientation;}


Note that this is insanely expensive. The code is geared towards accuracy, not speed. An lower order solver would probably be a better balance. I also apologize for not being particularly eloquent, so feel free to ask me even more questions.
...and do not wildly extrapolate. Just because Saddam Hussein gassed Kurds in 1990 doesn't mean he eats babies' brains.
Quote:Original post by Sneftel
Quote:Original post by oliii
In these kind fo systems (let it be A.I. sensors, or how much a A.I. car should apply force to correct it's steering), I usually look 'forward' in time, instead of the instantaneous inputs.


What you're describing is known as a "PID controller". Genjix, Googling for that term should turn up plenty of info on what is a very common problem.


Ah cool, I'll try to remember the name for the next time [grin]

Everything is better with Metal.

This topic is closed to new replies.

Advertisement