Energy loss problem with physics objects bouncing off each other.

Hi, can someone who knows about vector math tell me what I'm doing wrong here.

I've made up a physics demo, using some of Sswift's physics code, which has 100 balls bouncing around inside a cube with gravity applied.

It looks fairly realistic, for the most part, however there seems to be a problem with energy loss when the balls bounce off each other. They go from flying around so fast that they're barely visible, to lying on the floor of the cube barely moving, within a few minutes.

This is with a maximum elasticity value of 1.0. The problem was also present before I factored elasticity into the code. This only shows up if the balls are enabled to bounce off each other. If they're only bouncing off the walls then there's no apparent energy loss. My physics code is below.

Also, how would I go about, mathematically applying thrust to a physics object parellel to the plane of the surface it's standing on, so that it can be properly 'pushed' along an inclined surface (for a game character walking up and down ramps and hills, for example). {edit}I think multiplying the x, y, z, velocity vectors by the respective collision normals to align them to the space of the collision plane, applying the thrust to the x and y vextors, as required, and then multplying the three vectors by the normals calculated from the original velocity vectors to convert them back to global space, would work. I'll have to try it out.{/edit}

--- Physics code

; Note: Delta time has already been applied to the gravity# variable at the start of the main loop.

Function apply_basic_physics( the_entity_data.ENTITY_STRUCT )
the_entity = the_entity_data\entity_handle
vx# = the_entity_data\physics_data\vx#
vy# = the_entity_data\physics_data\vy#
vz# = the_entity_data\physics_data\vz#
;the_mass# = the_entity_data\physics_data\mass#
the_elasticity# = the_entity_data\physics_data\elasticity#
;the_friction# = the_entity_data\physics_data\friction#

For the_collision_index = 1 To CountCollisions ( the_entity )

the_entity_collided_with = CollisionEntity ( the_entity, the_collision_index )
the_target_entitiy_type = GetEntityType ( the_entity_collided_with )

; Get the data structure of the target collision object.
the_target_entity_data.ENTITY_STRUCT = Object.ENTITY_STRUCT( EntityName$ ( the_entity_collided_with ) )

; Get the collision normals.
coll_nx# = CollisionNX# ( the_entity, the_collision_index )
coll_ny# = CollisionNY# ( the_entity, the_collision_index )
coll_nz# = CollisionNZ# ( the_entity, the_collision_index )

; Compute the dot product of the avatars's motion vector and the normal of the surface collided with.
VdotN# = vx# * coll_nx# + vy# * coll_ny# + vz# * coll_nz#

; Check if the target collision object has physics properties.
If the_target_entity_data\physics_data <> Null

; Calculate the normal force to be applied to the source and target objects.
NFx# = -the_elasticity# * coll_nx# * VdotN#
NFy# = -the_elasticity# * coll_ny# * VdotN#
NFz# = -the_elasticity# * coll_nz# * VdotN#

; Apply the positive normal force to the motion vector of the source object.
vx# = vx# + NFx#
vy# = vy# + NFy#
vz# = vz# + NFz#

; Apply the negative normal force to the motion vector of the target object.
the_target_entity_data\physics_data\vx# = the_target_entity_data\physics_data\vx# - NFx#
the_target_entity_data\physics_data\vy# = the_target_entity_data\physics_data\vy# - NFy#
the_target_entity_data\physics_data\vz# = the_target_entity_data\physics_data\vz# - NFz#

Else

; Calculate and add the normal force to the motion vector of the source object.
vx# = vx# + ( -( 1.0 + the_elasticity# ) * coll_nx# * VdotN# )
vy# = vy# + ( -( 1.0 + the_elasticity# ) * coll_ny# * VdotN# )
vz# = vz# + ( -( 1.0 + the_elasticity# ) * coll_nz# * VdotN# )

EndIf

Next

vy# = vy# + gravity#

PositionEntity the_entity, EntityX# ( the_entity ) + ( vx# * Delta_Time# ), EntityY# ( the_entity ) + ( vy# * Delta_Time# ), EntityZ# ( the_entity ) + ( vz# * Delta_Time# )

the_entity_data\physics_data\vx# = vx#
the_entity_data\physics_data\vy# = vy#
the_entity_data\physics_data\vz# = vz#


DBG_total_energy# = DBG_total_energy# + Sqr# ( vx# * vx# + vy# * vy# + vz# * vz# ) ; DEBUG CODE :: Calculate the ball's energy.

End Function

I seem to have sorted this issue out. Most of the loss in energy was simply due to the opposing vectors cancelling each other out. The remainder was probably due to the collision not being recorded by a slower moving object moving away from a faster moving object that hits it.

I added some code to check for and factor in this type of collision, and now, after a while, the system settles down to an equilibrium state where no energy is being lost.

Here's my new code.

Function apply_basic_physics( the_entity_data.ENTITY_STRUCT ) ; the_collision_data.COLLISION_LIST_STRUCT )
; Note:
; The source collision object here is always a physics object. The target collision object may or may not be.
;the_entity_data.ENTITY_STRUCT = the_collision_data\source_entity_data
the_physics_data.PHYSICS_STRUCT = the_entity_data\physics_data
the_entity = the_entity_data\entity_handle
vx# = the_physics_data\vx#
vy# = the_physics_data\vy#
vz# = the_physics_data\vz#
;the_mass# = the_physics_data\mass#
the_elasticity# = the_physics_data\elasticity#
;the_friction# = the_physics_data\friction#

;the_target_entity_data.ENTITY_STRUCT = the_collision_data\target_entity_data

For the_collision_index = 1 To CountCollisions ( the_entity )

the_entity_collided_with = CollisionEntity ( the_entity, the_collision_index )
;the_target_entitiy_type = GetEntityType ( the_entity_collided_with )

; Get the data structure of the target collision object.
the_target_entity_data.ENTITY_STRUCT = Object.ENTITY_STRUCT( EntityName$ ( the_entity_collided_with ) )
the_target_physics_data.PHYSICS_STRUCT = the_target_entity_data\physics_data

; Get the collision normals.
coll_nx# = CollisionNX# ( the_entity, the_collision_index )
coll_ny# = CollisionNY# ( the_entity, the_collision_index )
coll_nz# = CollisionNZ# ( the_entity, the_collision_index )

; Compute the dot product of the avatars's motion vector and the normal of the surface collided with.
VdotN# = vx# * coll_nx# + vy# * coll_ny# + vz# * coll_nz#

; Check if the target collision object has physics properties.
If the_target_physics_data <> Null

; Calculate the normal force to be applied to the source and target objects.
NFx# = -the_elasticity# * coll_nx# * VdotN#
NFy# = -the_elasticity# * coll_ny# * VdotN#
NFz# = -the_elasticity# * coll_nz# * VdotN#

; Apply the positive normal force to the motion vector of the source object.
vx# = vx# + NFx#
vy# = vy# + NFy#
vz# = vz# + NFz#

; Apply the negative normal force to the motion vector of the target object.
the_target_physics_data\vx# = the_target_physics_data\vx# - NFx#
the_target_physics_data\vy# = the_target_physics_data\vy# - NFy#
the_target_physics_data\vz# = the_target_physics_data\vz# - NFz#

If Sgn ( vz# * coll_nz# ) = Sgn ( the_target_physics_data\vz# * coll_nz# )

VdotN# = the_target_physics_data\vx# * coll_nx# + the_target_physics_data\vy# * coll_ny# + the_target_physics_data\vz# * coll_nz#

; Calculate the normal force to be applied to the source and target objects.
NFx# = -the_target_physics_data\elasticity# * coll_nx# * VdotN#
NFy# = -the_target_physics_data\elasticity# * coll_ny# * VdotN#
NFz# = -the_target_physics_data\elasticity# * coll_nz# * VdotN#

; Apply the positive normal force to the motion vector of the target object.
the_target_entity_data\physics_data\vx# = the_target_entity_data\physics_data\vx# + NFx#
the_target_entity_data\physics_data\vy# = the_target_entity_data\physics_data\vy# + NFy#
the_target_entity_data\physics_data\vz# = the_target_entity_data\physics_data\vz# + NFz#

; Apply the negative normal force to the motion vector of the source object.
vx# = vx# - NFx#
vy# = vy# - NFy#
vz# = vz# - NFz#

EndIf

Else

; Calculate and add the normal force to the motion vector of the source object.
vx# = vx# + ( -( 1.0 + the_elasticity# ) * coll_nx# * VdotN# )
vy# = vy# + ( -( 1.0 + the_elasticity# ) * coll_ny# * VdotN# )
vz# = vz# + ( -( 1.0 + the_elasticity# ) * coll_nz# * VdotN# )

EndIf

Next

vy# = vy# + gravity#

PositionEntity the_entity, EntityX# ( the_entity ) + ( vx# * Delta_Time# ), EntityY# ( the_entity ) + ( vy# * Delta_Time# ), EntityZ# ( the_entity ) + ( vz# * Delta_Time# )

DBG_total_energy# = DBG_total_energy# + Sqr# ( vx# * vx# + vy# * vy# + vz# * vz# ) ; DEBUG CODE :: Calculate the ball's energy.

the_physics_data\vx# = vx#
the_physics_data\vy# = vy#
the_physics_data\vz# = vz#

End Function

Not really your problem, but something that might be of interest ...

Long time ago when I were a lad, some of the guys in the Physics lab at the university where I used to work did a lot of very long physical simulations. I seem to recall that they kept track of the total energy in their system. Due to small floating point errors, it was common for the energy to gradually drift, so they applied a fudge whereby any excess energy was redistributed evenly (or any deficit was rectified). This only needed to be done infrequently.

Yes, I was thinking of something similar at one point. Only really worth doing, though, if you've got a lot of physics collisions going on over a long period of time.

The code above has an elasticity constant, but 1.0 is not the correct value for perfectly inelastic collsions. (Those where no energy is lost)

NFx# = -the_target_physics_data\elasticity# * coll_nx# * VdotN#

See that line? Well in a perfectly inelastic collision, that elasticity value there would be 2.0. So if you set it to 1.0, then you SHOULD be losing almost all your energy instantly, because that would mean a perfectly elastic collision.

This line SHOULD say:

NFx# = -(1.0 + the_target_physics_data\elasticity#) * coll_nx# * VdotN#

Then a value of 1.0 for elasticity will represent a bouncy ball and 0.0 will represent a pool ball.

Sort of.

Those physics don't take into account mass, or how much energy is transferred to the other object. So 1.0 means (I think it's called) an elastic collision where the ball KEEPS all of it's energy when it collides with the ground, and 0.0 represents an inelastic collision where the ball TRANSFERS all of it's energy to the other object when they collide.

Here's the kicker though. Because I'm not transferring energy, if you have anything other than an elasticity of 1.0, then each time energy should be transferred from one ball to another (or from a ball to the floor) the energy is simply "lost". And that violates the laws of physics. Of course, if the balls collided with the floor presumably it would have so much mass anyhow that if the balls were not completely inelastic (0.0) then when they hit the floor it would get the energy but only barely budge. But of course you don't even budge the floor you treat it like an object with infinite mass. Either way, an even slightly inelastic collision with the floor will cause the balls to eventually slow down because the floor will just store the energy if anything at all, and not give it back.

Imagine for example pool balls colliding with a wall in a space station. Most likely they'd hit the wall, and then stick there. The wall gets all the energy and doesn't give it back. The ship moves a tiny bit but that is imperceptible because the difference in mass is so great. Either way, if the collision is inelastic, then the balls come to a rest.

Now if two inelastic balls of the same SIZE hit eachother in space, well that's a different story. They'll "bounce" off one another. One will gain the velocity and direction which the other previously had. Like when you have those steel balls on string and you let two go at the same time, they both bounce back, over and over again. Or if you release three, then three pop out on the other side.

I'm not really sure though why three pop out on the other side, rather than just one moving at three times the speed. Hm... I suppose it may have something to do with harmonics or something of the metal. Or not. I'll have to ask a physicist. :-)

Oh and what Steve said also applies. Floating point numbers aren't perfectly accurate. Energy can be lost. You'd have to be really careful with your calculations to make sure that the velocity A plus the veloicity of B before the collision equals the velocity of A plus the velocity of B AFTER the collision. And even then you might not be safe.

Here's some physics for you of two objects with variable mass and different coefficients of friction colliding with one another. Still not 100% proper physics though. I didn't know what to do when the coefficients of friction (elasticity) are different.

Be aware that you must make sure that you do this only once for each pair of colliding objects. If you just detect every object that collides and then determine what it collided with and do this code, then you'll UNDO the energy transfer that takes place when you do it the second time.

I couldn't think of a way to allow you to do the collisoon code just for one of the objects, transfer the energy, and then do the other. Like if two balls collide with eachother head on, if you just did one first and then the other after, with the first one, you'd just cancel out it's momentum and end up with 0 movement, and then have nothing to transfer when you get to the second. They have to transfer their energy simulatanously. I suppose you could maintain a before and after state for all objects though and operate on the before state and set the after. Hm...

     ; Convert the player's velocities and directions into vectors.
	    	        MagDirToVector(ThisBike)
					MagDirToVector(HitPlayer)

					;DebugLog("Before:")
					;DebugLog(HitPlayer\VECTOR_X#)
					;DebugLog(HitPlayer\VECTOR_Y#)
					;DebugLog(HitPlayer\VECTOR_Z#)
					;DebugLog("")

					; Transfer momentum from one entity to the other.

						; First, get caclulate a vector from the center of mass of the first player,
						; to the center of mass of the second.
						
							CNx# = EntityX#(ThisBike\Pivot, True) - EntityX#(HitPlayer\Pivot, True)
							CNy# = EntityY#(ThisBike\Pivot, True) - EntityY#(HitPlayer\Pivot, True)
							CNz# = EntityZ#(ThisBike\Pivot, True) - EntityZ#(HitPlayer\Pivot, True)
						
						; Normalize the vector.
		
							CNl# = Sqr(CNx#*CNx# + CNy#*CNy# + CNz#*CNz#)
							CNx# = CNx# / CNl#
							CNy# = CNy# / CNl#
							CNz# = CNz# / CNl#

						; Find the length of each of the movement vectors along the normal.

							;V1dotN# = ThisBike\VECTOR_X#*CNx#  + ThisBike\VECTOR_Y#*CNy#  + ThisBike\VECTOR_Z#*CNz#
							;V2dotN# = HitPlayer\VECTOR_X#*CNx# + HitPlayer\VECTOR_Y#*CNy# + HitPlayer\VECTOR_Z#*CNz#

						; Get the mass of the objects involved in the collision.

							Mass1# = 1.0
							Mass2# = 1.0

						; Get the coefficient of restitution of each object.
						
							; The coefficient of restitution is a value between 0 and 1 which defines how elastic
							; the objects that collide are.  A metal ball is inelastic, and would have a coefficient
							; of 0... all energy would be transferred from the first object to the second during
							; a collision.  A basketball on the other hand has a coefficient of 0.6...  So it bounces
							; when it hits the floor rather than transferring all it's energy to the floor.  
							; A baseball, btw, has a coefficient of 0.55.
													
							E1# = 0.25
							E2# = 0.25
							
						; Calculate the coefficient of restitution to use for the collision.
							
							; Note:
							; This is probably not the correct physical way to model a collision when the two
							; objects have different coefficients of friction, but it should work well enough.
																										
							; If the first object is softer than the second object...
							If E1# > E2#
								E# = E1#		; Use the first object's coefficient.
							Else
								E# = E2#		; Use the second object's coefficient.
							EndIf

						; Caclulate the new velocities of the two objects.
							
							C# = CNx# * (ThisBike\VECTOR_X# - HitPlayer\VECTOR_X#) + CNy# * (ThisBike\VECTOR_Y# - HitPlayer\VECTOR_Y#) + CNz# * (ThisBike\VECTOR_Z# - HitPlayer\VECTOR_Z#)
							
							P1# = ((Mass2# * C#) / (Mass1# + Mass2#)) * (1.0 + E#)
							
							ThisBike\VECTOR_X# = ThisBike\VECTOR_X# - P1# * CNx#
							ThisBike\VECTOR_Y# = ThisBike\VECTOR_Y# - P1# * CNy#
							ThisBike\VECTOR_Z# = ThisBike\VECTOR_Z# - P1# * CNz#

							P2# = ((Mass1# * C#) / (Mass1# + Mass2#)) * (1.0 + E#)
							
							HitPlayer\VECTOR_X# = HitPlayer\VECTOR_X# + P2# * CNx#
							HitPlayer\VECTOR_Y# = HitPlayer\VECTOR_Y# + P2# * CNy#
							HitPlayer\VECTOR_Z# = HitPlayer\VECTOR_Z# + P2# * CNz#
	
						; Convert the movement vector back into a magnitude and direction.
							VectorToMagDir(ThisBike,  ThisBike\VECTOR_X#,  ThisBike\VECTOR_Y#,  ThisBike\VECTOR_Z#)
							
						; Store the new movement vector for the player collided with.
						
						    HitPlayer\Vx# = HitPlayer\VECTOR_X#
						    HitPlayer\Vy# = HitPlayer\VECTOR_Y#
							HitPlayer\Vz# = HitPlayer\VECTOR_Z#

Sswift, the line:

NFx# = -the_target_physics_data\elasticity# * coll_nx# * VdotN#

in my code is for collisions where physics objects are bouncing off each other, so the resulting normal force is added to the source object and subtracted from the target object which means its accumulative value should be equivalent to:

NFx# = -(1.0 + the_target_physics_data\elasticity#) * coll_nx# * VdotN#


Regarding my second set of code above, it seems to be a mistake to factor in the unregistered collision of a slower ball being hit from behind by a faster ball as this seems to result in the normal force being cancelled out (unless there's an error in my code somewhere). A more efficient check for missed collisions using:

----

t = 1
the_number_of_target_collisions = CountCollisions ( the_entity_collided_with )

For the_target_collision_index = 1 To the_number_of_target_collisions
If CollisionEntity ( the_entity_collided_with, the_target_collision_index ) = the_entity Then t = 0 : Exit
Next

If t

----

instead of:

----

If Sgn ( vz# * coll_nz# ) = Sgn ( the_target_physics_data\vz# * coll_nz# )

----

made the problem more obvious. The faster ball would simply move around the slower ball like it's climbing over it. There were a few other bugs in this code as well.


Here's a link to the compiled program, anyway. Everything seems to be behaving as it should apart from the energy loss problem.

http://www.brightrealm.com/downloads/demos/simple_demos/Bouncing_balls.zip