I'm hoping this thread will create a discussion regarding some of the various integration methods out there.
I'll get the ball rolling with Euler Integration.
I'm going to keep it pretty simple (yeah, right) by deliberately ignoring ROTATION and by only mentioning DYNAMICS equations (at least for now).
For those of you who managed to pick up some physics during your school days, you may recall the following formulae:
CALCULATE ACCELERATION
It should be clear enough that if we presume to KNOW the values for "f" and "m", then we can substitute the formula and calculate a=f/m , yes?
CALCULATE VELOCITY - A TASTE OF INTEGRATION
It follows that having calculated a=f/m, and presuming we KNOW deltaTIME, we can substitute again and calculate dV=a*dT.
We might choose at this moment to INTEGRATE the above two formulae, such that we can find dV=(f/m)*dT in a single step.
Also, we can integrate once more to find V
If I=initial velocity , to find new velocity V:
CALCULATE DISTANCE - INTEGRATE AGAIN?
Again following this logic, since we calculated dV=a*dT, and since we already know the value of dT, we can once more substitute to calculate dD=V*dT.
Let's take it to the next level.
By substitution, dD = (I + ((f/m)*dT)) * dT
We could if we wanted to integrate yet again..
If P = initial position, to find new position Q:
As long as we remember to grab the new Velocity out of the formula, theres no reason why we can't integrate all these calculations into one, and write a procedure which calculates new Velocity and Position at once, and call it UpdateObject or something like that.
We now have a basic set of formulae for calculating the new position and velocity of an object based on its old position, velocity, and summed input Force (includes Gravity), after some amount of Time has passed (known as the TimeStep).
What may not be so obvious is that we may use what we have learned to calculate the position and velocity of the object at any future (or past!) moment in Time.
This becomes critical when we begin attempting to implement a Collision Response. What we will be doing is calculating the exact moment when a collision occurs, then from that, calculating the position of the object at the moment of collision, then modifying the object velocity (reflecting the collision ray), and then repositioning the object to where it SHOULD have been, had it not passed through the wall it was meant to bounce off, ie move the object along the new r"reflected" ray by the same distance it had passed through the wall when it collided.
All we are really doing is throwing our current object values and various Time values at the formulae above.
One thing I note immediately with regard to gamecoding is that the value deltaTime is quite simply the elapsed time between frames.
Most of the sourcecode I've looked at uses "fixed timesteps" which create a flawed physical model simulation that is prone to "jerkyness", especially where the actual timestep fluctuates greatly.
Is this just another example of cut and paste coding, or is there some reason why a hardcoded timestep is a Good Thing? (I am led to believe that Runge-Kutta Integration performs poorly with fluctuating timesteps, but I really don't understand the benefits of RK either).
Who wants to shed some light on Runge-Kutta and Verlette integration, and who can describe something else, maybe something new?
I'll get the ball rolling with Euler Integration.
I'm going to keep it pretty simple (yeah, right) by deliberately ignoring ROTATION and by only mentioning DYNAMICS equations (at least for now).
For those of you who managed to pick up some physics during your school days, you may recall the following formulae:
f=ma (FORCE = MASS * ACCELERATION)a=dV/dT (ACCELERATION = deltaVELOCITY / deltaTIME)v=dD/dT (VELOCITY = deltaDISTANCE / deltaTIME)CALCULATE ACCELERATION
It should be clear enough that if we presume to KNOW the values for "f" and "m", then we can substitute the formula and calculate a=f/m , yes?
CALCULATE VELOCITY - A TASTE OF INTEGRATION
It follows that having calculated a=f/m, and presuming we KNOW deltaTIME, we can substitute again and calculate dV=a*dT.
We might choose at this moment to INTEGRATE the above two formulae, such that we can find dV=(f/m)*dT in a single step.
Also, we can integrate once more to find V
If I=initial velocity , to find new velocity V:
V = I + ((f/m)*dT)CALCULATE DISTANCE - INTEGRATE AGAIN?
Again following this logic, since we calculated dV=a*dT, and since we already know the value of dT, we can once more substitute to calculate dD=V*dT.
Let's take it to the next level.
By substitution, dD = (I + ((f/m)*dT)) * dT
We could if we wanted to integrate yet again..
If P = initial position, to find new position Q:
Q = P + ( I + ((f/m)*dT)) * dT )As long as we remember to grab the new Velocity out of the formula, theres no reason why we can't integrate all these calculations into one, and write a procedure which calculates new Velocity and Position at once, and call it UpdateObject or something like that.
We now have a basic set of formulae for calculating the new position and velocity of an object based on its old position, velocity, and summed input Force (includes Gravity), after some amount of Time has passed (known as the TimeStep).
What may not be so obvious is that we may use what we have learned to calculate the position and velocity of the object at any future (or past!) moment in Time.
This becomes critical when we begin attempting to implement a Collision Response. What we will be doing is calculating the exact moment when a collision occurs, then from that, calculating the position of the object at the moment of collision, then modifying the object velocity (reflecting the collision ray), and then repositioning the object to where it SHOULD have been, had it not passed through the wall it was meant to bounce off, ie move the object along the new r"reflected" ray by the same distance it had passed through the wall when it collided.
All we are really doing is throwing our current object values and various Time values at the formulae above.
One thing I note immediately with regard to gamecoding is that the value deltaTime is quite simply the elapsed time between frames.
Most of the sourcecode I've looked at uses "fixed timesteps" which create a flawed physical model simulation that is prone to "jerkyness", especially where the actual timestep fluctuates greatly.
Is this just another example of cut and paste coding, or is there some reason why a hardcoded timestep is a Good Thing? (I am led to believe that Runge-Kutta Integration performs poorly with fluctuating timesteps, but I really don't understand the benefits of RK either).
Who wants to shed some light on Runge-Kutta and Verlette integration, and who can describe something else, maybe something new?
Hello Homer,
Just some comments from me.
You mentioned Euler integration but I have not seen anything on it here. If I am not wrong Euler integration is a numerial method right?
I would like to point out that your definitions are wrong
a=dV/dT (ACCELERATION = deltaVELOCITY / deltaTIME)
should be
a=dV/dT (ACCELERATION = deeVELOCITY / deeTIME)
Your usage of delta is wrong. deltavelocity is completely different from deevelocity (deeVee). delta means a small increment or decrement while deeVee on the other hand means limits as delta approaches 0.
Anyway your integration does not look like integration at all. :)
PS: I just went through a series of lectures on mechanics (including particle under variable force) , so I guess now is the right time to dicuss :)
Regards,
Victor
Just some comments from me.
You mentioned Euler integration but I have not seen anything on it here. If I am not wrong Euler integration is a numerial method right?
I would like to point out that your definitions are wrong
a=dV/dT (ACCELERATION = deltaVELOCITY / deltaTIME)
should be
a=dV/dT (ACCELERATION = deeVELOCITY / deeTIME)
Your usage of delta is wrong. deltavelocity is completely different from deevelocity (deeVee). delta means a small increment or decrement while deeVee on the other hand means limits as delta approaches 0.
Anyway your integration does not look like integration at all. :)
PS: I just went through a series of lectures on mechanics (including particle under variable force) , so I guess now is the right time to dicuss :)
Regards,
Victor
a=dV/dT
I was taught many years ago that dV = change in velocity.
a=(d2-d1) / (t2-t1) is not correct?
You mention deeVelocity and also that it should be limited as it approaches zero. Limited in what way? Are you simply inferring that due to numerical inaccuracy that small values near zero should be replaced with zero, for example -.0001 = 0 = +.0001 ?
I wasn't aware this was being implied by the formulae, although I've seen such mechanisms in source, I was under the impression they were to prevent numerical drift.
As for Euler integration,
1. Set Time to its initial value.
2. Initialize the levels.
3. Compute the rates of change of the Levels at the current value of Time.
4. Use the rates of change to compute the Levels at Time + TIME STEP according to the formula:
LEVELTime+TIME STEP = LEVELTime + TIME STEP * RATETime
5. Add TIME STEP to Time.
6. Repeat steps 3-5 until Time is equal to FINAL TIME.
It's basically using the formulae I gave with incremental value of Time.
Correct?
AFAIK Euler Integration is a method of APPROXIMATING the "levels" of a closed system at any moment in time, but it assumes that the "rates of change" are CONSTANT for a given TimeStep. For example, that Acceleration at T1 is the same as Acceleration at T2, is the same as Acceleration at any moment between T1 and T2. This is simply incorrect and is the cause of the inaccuracy of Euler's method. The inaccuracy can be minimized by using very small TimeSteps, but this is simply not a viable solution in a time-critical application. For example, one source I saw divided each "elapsed frame" worth of Time into 150 TimeSteps, then applied 150 iterations of physics equations to each and every object, including 150 intersection tests between each object and everything else in the world... this is WAY too expensive.
I'm NOT a physics guru, and CERTAINLY NOT a grand master of simulation, which is why I chose to post this thread here in this forum.
I know a lot of you are currently studying this stuff, some of you are studying physics and CS at once, and more than likely have your own ideas to bring to the table.
By all means, do so :)
Thanks to Roticv for the feedback, much appreciated.
I was taught many years ago that dV = change in velocity.
a=(d2-d1) / (t2-t1) is not correct?
You mention deeVelocity and also that it should be limited as it approaches zero. Limited in what way? Are you simply inferring that due to numerical inaccuracy that small values near zero should be replaced with zero, for example -.0001 = 0 = +.0001 ?
I wasn't aware this was being implied by the formulae, although I've seen such mechanisms in source, I was under the impression they were to prevent numerical drift.
As for Euler integration,
1. Set Time to its initial value.
2. Initialize the levels.
3. Compute the rates of change of the Levels at the current value of Time.
4. Use the rates of change to compute the Levels at Time + TIME STEP according to the formula:
LEVELTime+TIME STEP = LEVELTime + TIME STEP * RATETime
5. Add TIME STEP to Time.
6. Repeat steps 3-5 until Time is equal to FINAL TIME.
It's basically using the formulae I gave with incremental value of Time.
Correct?
AFAIK Euler Integration is a method of APPROXIMATING the "levels" of a closed system at any moment in time, but it assumes that the "rates of change" are CONSTANT for a given TimeStep. For example, that Acceleration at T1 is the same as Acceleration at T2, is the same as Acceleration at any moment between T1 and T2. This is simply incorrect and is the cause of the inaccuracy of Euler's method. The inaccuracy can be minimized by using very small TimeSteps, but this is simply not a viable solution in a time-critical application. For example, one source I saw divided each "elapsed frame" worth of Time into 150 TimeSteps, then applied 150 iterations of physics equations to each and every object, including 150 intersection tests between each object and everything else in the world... this is WAY too expensive.
I'm NOT a physics guru, and CERTAINLY NOT a grand master of simulation, which is why I chose to post this thread here in this forum.
I know a lot of you are currently studying this stuff, some of you are studying physics and CS at once, and more than likely have your own ideas to bring to the table.
By all means, do so :)
Thanks to Roticv for the feedback, much appreciated.
For the sake of discussion of Particle Physics, I'd like to use the example of a Sphere object. Let's describe a Sphere using the following structure:
Let's assume that the Sphere begins life at some position (X,Y,Z) in space, with zero velocity, and Time=0. Let's further assume that the only external force acting on the Sphere is Gravity. Finally, let's assume that the Mass of the Sphere is 2.0
Now let's assume that some Time has passed, such that Time=0.500 (half a second has elapsed).
Let's start calculating :)
First, we wanna know the Acceleration in order to obtain Velocity.
The formula f=ma (or rather, a=f/m) is what we need to use.
Since we know the Mass in this case is 1, its quite simple:
a=f/2
Let's assume the Gravity force has a value of 0.86, represented as a Vector:
G=<0.0f, -0.86f, 0.0f>
We can now calculate the Acceleration a due to Gravity G:
a=<0.0f, -0.43f, 0.0f>
Having done this we can now calculate the change in velocity due to acceleration over Time:
a=dV/dT, or rather, dV=a*dT
Since we know that dT=0.500 seconds, it follows that
dV=<0.0f, -0.215f, 0.0f>
and thus we can calculate the New Velocity=OldVelocity+dV
(not shown).
Having done that, we can calculate the change in Position due to Velocity over Time, and then we have a new Position and Velocity.
Should we wish to impart further external forces on the object, we should simply sum their vectors into the "f" vector.
Now we have our Sphere object accelerating under Gravity.
Let's suppose that we have a flat Plane on the "floor" which we want our Sphere to collide with, and then deflect from (bounce off).
We will assume that the Surface has a Normal vector already calculated.
We will also assume we are capable of calculating the Distance from a given Plane to a given Point In Space.
We should measure the Distance from the Plane to the Sphere Center.
If the distance is greater than the Radius of the Sphere, then the Sphere is NOT intersecting the Plane, and everything is cool.
If the distance is <= the Sphere Radius, then the Sphere has intersected the Plane, and we need to "fix" the situation.
I need to make a diagram or something to illustrate this further, before I describe the formulae for resolving the collision. The beautiful part of the formulae I will be documenting is that the plane doesn't have to be on the "floor", it could be a ceiling, a wall, or part of a slope, the collision will still be resolved correctly no matter what is going on.
For the purposes of discussion I choose to describe a sphere falling onto a flat floor plane with its Normal facing UP, simply because its easy to visualize.
Now I'm off to make some crap art :)
Sphere struct
vPosition Vec3 <> ;Sphere Position in 3D space (CURRENT)
vOldPos Vec3 <> ;Sphere Position in 3D space (OLD)
vVelocity Vec3 <> ;Sphere Velocity
fRadius REAL8 ?
fMass REAL4 ?
Sphere ends
Let's assume that the Sphere begins life at some position (X,Y,Z) in space, with zero velocity, and Time=0. Let's further assume that the only external force acting on the Sphere is Gravity. Finally, let's assume that the Mass of the Sphere is 2.0
Now let's assume that some Time has passed, such that Time=0.500 (half a second has elapsed).
Let's start calculating :)
First, we wanna know the Acceleration in order to obtain Velocity.
The formula f=ma (or rather, a=f/m) is what we need to use.
Since we know the Mass in this case is 1, its quite simple:
a=f/2
Let's assume the Gravity force has a value of 0.86, represented as a Vector:
G=<0.0f, -0.86f, 0.0f>
We can now calculate the Acceleration a due to Gravity G:
a=<0.0f, -0.43f, 0.0f>
Having done this we can now calculate the change in velocity due to acceleration over Time:
a=dV/dT, or rather, dV=a*dT
Since we know that dT=0.500 seconds, it follows that
dV=<0.0f, -0.215f, 0.0f>
and thus we can calculate the New Velocity=OldVelocity+dV
(not shown).
Having done that, we can calculate the change in Position due to Velocity over Time, and then we have a new Position and Velocity.
Should we wish to impart further external forces on the object, we should simply sum their vectors into the "f" vector.
Now we have our Sphere object accelerating under Gravity.
Let's suppose that we have a flat Plane on the "floor" which we want our Sphere to collide with, and then deflect from (bounce off).
We will assume that the Surface has a Normal vector already calculated.
We will also assume we are capable of calculating the Distance from a given Plane to a given Point In Space.
We should measure the Distance from the Plane to the Sphere Center.
If the distance is greater than the Radius of the Sphere, then the Sphere is NOT intersecting the Plane, and everything is cool.
If the distance is <= the Sphere Radius, then the Sphere has intersected the Plane, and we need to "fix" the situation.
I need to make a diagram or something to illustrate this further, before I describe the formulae for resolving the collision. The beautiful part of the formulae I will be documenting is that the plane doesn't have to be on the "floor", it could be a ceiling, a wall, or part of a slope, the collision will still be resolved correctly no matter what is going on.
For the purposes of discussion I choose to describe a sphere falling onto a flat floor plane with its Normal facing UP, simply because its easy to visualize.
Now I'm off to make some crap art :)
Note that the Sphere has a field for "Old Position".
The order of operations might look something like this:
-Calculate new velocity, and thus new position
-Check for intersection(s)...
-If Intersection detected,
--Calculate T=Time Of Intersection
--Calculate Position and Velocity at T
--Set OLD POSITION to Collision Point
--Calculate Reflection as Force
--Apply Reflection as input Force
--Calculate Position and Velocity at end of current TimeStep
--Set CURRENT POSITION to the above values
-Else
--Set OLD POSITION = CURRENT POSITION
-EndIf
How does that sound?
The order of operations might look something like this:
-Calculate new velocity, and thus new position
-Check for intersection(s)...
-If Intersection detected,
--Calculate T=Time Of Intersection
--Calculate Position and Velocity at T
--Set OLD POSITION to Collision Point
--Calculate Reflection as Force
--Apply Reflection as input Force
--Calculate Position and Velocity at end of current TimeStep
--Set CURRENT POSITION to the above values
-Else
--Set OLD POSITION = CURRENT POSITION
-EndIf
How does that sound?
Hello Homer,
You are correct that a = dv/dt.
However dv/dt is defined as lim ?t -> 0 (v+?v)/(t+?t) where ?t = deltatime and ?v = deltavelocity (Just some elementary calculus).
a = (v1 - v2)/(t1 - t2)
Provided v is constant. If not you would have to solve using real integration.
I am not very sure about Euler's integration, I would read up more on it before commenting on that.
By the way collisions somehow reminds me of principle of conservation of momentum and resitutition. The kind of mechanics crap that I have to put up with for the rest of the year.
PS: I like to think of g as 9.81ms^-1. :-D
PSS: Off to do some mechanics homework :)
Sorry if I am not really helping you, I am just looking more at the mechanical perpective rather than how to code it etc etc.
You are correct that a = dv/dt.
However dv/dt is defined as lim ?t -> 0 (v+?v)/(t+?t) where ?t = deltatime and ?v = deltavelocity (Just some elementary calculus).
a = (v1 - v2)/(t1 - t2)
Provided v is constant. If not you would have to solve using real integration.
I am not very sure about Euler's integration, I would read up more on it before commenting on that.
By the way collisions somehow reminds me of principle of conservation of momentum and resitutition. The kind of mechanics crap that I have to put up with for the rest of the year.
PS: I like to think of g as 9.81ms^-1. :-D
PSS: Off to do some mechanics homework :)
Sorry if I am not really helping you, I am just looking more at the mechanical perpective rather than how to code it etc etc.
is euler integrating the same as just calculating the riemann sum? (which is part of the definition of integrating).
I have no idea what the Riemann Sum is about.
Euler integration assumes that rates of change are constant between the start and end times, and therefore that various values can be extrapolated at any moment BETWEEN the start and end times.
Since things like Acceleration are exceedingly unlikely to be linear (let alone constant), Euler's method only produces reasonably accurate results when the time intervals applied are very small, and is only perfectly accurate when the time intervals are infinitely small.
Therefore its really not much good for this particular application.
The resolving of a collision is indeed about conservation of energy (read momentum). In the example of a moving object colliding with a static planar surface, we see the linear and angular momentum of an object as both inputs to the collision and as outputs resulting from it. That is to say, depending on the angles involved in the collision, energy is transformed from one to the other - an object that is moving slowly but rotating quickly before it collides may be moving quickly but rotating slowly after the collision. The surface itself plays a small part in this because it has at the very least a "coefficient of restitution" which indicates how much energy is ABSORBED by the surface (is LOST to the object), and may well have other properties such as a friction coefficient.
Euler integration assumes that rates of change are constant between the start and end times, and therefore that various values can be extrapolated at any moment BETWEEN the start and end times.
Since things like Acceleration are exceedingly unlikely to be linear (let alone constant), Euler's method only produces reasonably accurate results when the time intervals applied are very small, and is only perfectly accurate when the time intervals are infinitely small.
Therefore its really not much good for this particular application.
The resolving of a collision is indeed about conservation of energy (read momentum). In the example of a moving object colliding with a static planar surface, we see the linear and angular momentum of an object as both inputs to the collision and as outputs resulting from it. That is to say, depending on the angles involved in the collision, energy is transformed from one to the other - an object that is moving slowly but rotating quickly before it collides may be moving quickly but rotating slowly after the collision. The surface itself plays a small part in this because it has at the very least a "coefficient of restitution" which indicates how much energy is ABSORBED by the surface (is LOST to the object), and may well have other properties such as a friction coefficient.
Ok heres that crap art I mentioned earlier..
Image #1 shows our Sphere at various Times.
Beginning at Time 0, our sphere is born into existance.
Time passes, and it begins to fall under the force of Gravity.
By Time 1 we can see it descending toward the Surface Plane.
However, at Time 2, the Sphere has penetrated the Surface considerably.
There must have been a collision at some moment between Time 1 and Time 2.
Image#2 shows our Sphere at Time 0 (same as before), at Time 1 (same as before) and also shows a new Time X (moment of impact) and Time 2 (where the bloody thing SHOULD be if everything is working).
This image also shows the sphere's Direction of Travel at each TimeStep using some arrows.
Having determined that an impact occured between T1 and T2 (as in Image#1) is obviously not enough information for accurate simulation.
What we really need is a way to determine the actual moment of impact (Time X), and the position, direction and speed of the Sphere at that moment.
I'm going to give you a procedure which is capable of A) determining if a collision has occured in the current TimeStep and B) returning the deltaTime (from T1) to the moment of impact, the position of the sphere at that moment, the direction, and the speed.
After that we'll take a look at how to resolve the collision.
Image #1 shows our Sphere at various Times.
Beginning at Time 0, our sphere is born into existance.
Time passes, and it begins to fall under the force of Gravity.
By Time 1 we can see it descending toward the Surface Plane.
However, at Time 2, the Sphere has penetrated the Surface considerably.
There must have been a collision at some moment between Time 1 and Time 2.
Image#2 shows our Sphere at Time 0 (same as before), at Time 1 (same as before) and also shows a new Time X (moment of impact) and Time 2 (where the bloody thing SHOULD be if everything is working).
This image also shows the sphere's Direction of Travel at each TimeStep using some arrows.
Having determined that an impact occured between T1 and T2 (as in Image#1) is obviously not enough information for accurate simulation.
What we really need is a way to determine the actual moment of impact (Time X), and the position, direction and speed of the Sphere at that moment.
I'm going to give you a procedure which is capable of A) determining if a collision has occured in the current TimeStep and B) returning the deltaTime (from T1) to the moment of impact, the position of the sphere at that moment, the direction, and the speed.
After that we'll take a look at how to resolve the collision.
As promised, here is a procedure which can not only determine whether a collision occurred during the current Frame, but also return the relative time when impact occurred, and the object position at the moment of impact. This allows us to "correct" the position of the object and then to "project" the position for (TimeStep-deltaTime) in order to put the bugger where it ought to be (at "T2", post-collision).
We can even determine if the object was ALREADY intersecting ("touching") the Surface, which will be important later when we want to apply STATIC FORCES to the object (ie make it slide/roll around on the ground).
Of course, an object should never be "already intersecting" for the example at hand. Our objects always begin life in mid-air.
This procedure is obviously FAR MORE USEFUL than one which simply returns TRUE or FALSE for a collision test, as it sets us well on track toward RESOLVING the collision.
We can even determine if the object was ALREADY intersecting ("touching") the Surface, which will be important later when we want to apply STATIC FORCES to the object (ie make it slide/roll around on the ground).
Of course, an object should never be "already intersecting" for the example at hand. Our objects always begin life in mid-air.
This procedure is obviously FAR MORE USEFUL than one which simply returns TRUE or FALSE for a collision test, as it sets us well on track toward RESOLVING the collision.
;===============================================
;This procedure tests for intersection of a given Plane and Sphere.
;If no collision is found, it returns FALSE.
;If a collision is found, it returns TRUE, and also returns
;the exact Time of the collision and Position of the Sphere at that moment.
;IN pNormal = The Plane given as Surface Normal, which we assume has been precalculated from the "first triangle" in the Quad.
;IN pSphere = ptr to Sphere struct
;OUT pFirstTouched = position of sphere('s origin) at moment of impact
;OUT pTouchTime = deltaTime to the moment of impact
;===============================================
SpherePlaneSweep proc pNormal:ptr Vec3, pSphere:ptr Sphere, pFirstTouched:ptr Vec3, pTouchTime:DWORD
local fd0:REAL8 ;<== distance from Plane to (previous position of sphere)
local fd1:REAL8 ;<== distance from Plane to (current position of sphere)
;Calculate distance from Plane to OLD POSITION of Sphere
mov ecx,pSphere
invoke PlaneDistance,pNormal, addr [ecx].Sphere.vOldPos
fstp fd0
;check if it was touching on previous "TimeStep"
;(ie is distance from Plane to OLD POSITION less than the Sphere's Radius?)
fld fd0
fabs
mov ecx,pSphere
fld [ecx].Sphere.fRadius
fcomp [ecx].Sphere.fRadius
__FJLE @F ;if( fabs(d0) <= r )
mov ebx,pFirstTouched ;Ci = C0
mov ecx,pSphere
fld [ecx].Sphere.vOldPos.X
fld [ecx].Sphere.vOldPos.Y
fld [ecx].Sphere.vOldPos.Z
fstp [ebx].Vec3.Z
fstp [ebx].Vec3.Y
fstp [ebx].Vec3.X
mov ebx,pTouchTime
fld r4_0_0 ;u = 0;
fstp REAL4 ptr[ebx]
return TRUE
@@:
mov ecx,pSphere
invoke PlaneDistance, pNormal, addr [ecx].Sphere.vPosition
fstp fd1
;//check if the sphere penetrated during this "TimeStep"
;if( d0>r && d1<r )
fld fd0 ;if( d0>r
fcomp [ecx].Sphere.fRadius
__FJG @F
fld fd1 ;&& d1<r)
fcomp [ecx].Sphere.fRadius
__FJL @F
;Message "Penetrated This TimeStep"
;u = (d0-r)/(d0-d1); //normalized time
mov ebx,pTouchTime
fld fd0
fsub [ecx].Sphere.fRadius
fld fd0
fsub fd1
fdiv
fstp dword ptr[ebx]
;Ci = (1-u)*C0 + u*C1; //point of first contact
mov ebx,pTouchTime
fld1
fsub REAL4 ptr[ebx]
fmul [ecx].Sphere.vOldPos.X
fld REAL4 ptr[ebx]
fmul [ecx].Sphere.vPosition.X
fadd
mov ebx,pFirstTouched
fstp [ebx].Vec3.X
mov ebx,pTouchTime
fld1
fsub REAL4 ptr[ebx]
fmul [ecx].Sphere.vOldPos.Y
fld REAL4 ptr[ebx]
fmul [ecx].Sphere.vPosition.Y
fadd
mov ebx,pFirstTouched
fstp [ebx].Vec3.Y
mov ebx,pTouchTime
fld1
fsub REAL4 ptr[ebx]
fmul [ecx].Sphere.vOldPos.Z
fld REAL4 ptr[ebx]
fmul [ecx].Sphere.vPosition.Z
fadd
mov ebx,pFirstTouched
fstp [ebx].Vec3.Z
return TRUE
@@:
;No penetration occured during this "TimeStep"
return FALSE
SpherePlaneSweep endp
There's just one problem with that procedure which I have noticed so far.
The returned "touchpoint" is the point of intersection ON THE PLANE.
It is NOT the "corrected position of the Sphere Origin".
We need to move the Sphere by an amount equal to its Radius and in the direction which is NEGATIVE to the pre-collision direction, ie, back the way it came from.
Since our Velocity contains both Speed and Direction, we must Normalize the Velocity in order to get the Direction, then negate it so we have "minus vDirection", then multiply it by the Sphere Radius to obtain a "positional offset", then add the Position , and store as Position.
I note a rather obvious optimisation here.
Rather than negating the xyz of the Direction vector, we can simply negate the Radius value instead and obtain the same result.
Not much of an optimisation, but every bit helps.
The returned "touchpoint" is the point of intersection ON THE PLANE.
It is NOT the "corrected position of the Sphere Origin".
We need to move the Sphere by an amount equal to its Radius and in the direction which is NEGATIVE to the pre-collision direction, ie, back the way it came from.
Since our Velocity contains both Speed and Direction, we must Normalize the Velocity in order to get the Direction, then negate it so we have "minus vDirection", then multiply it by the Sphere Radius to obtain a "positional offset", then add the Position , and store as Position.
I note a rather obvious optimisation here.
Rather than negating the xyz of the Direction vector, we can simply negate the Radius value instead and obtain the same result.
Not much of an optimisation, but every bit helps.
I know this code looks rather complex for a sphere that falls down onto a surface that faces up, so let me exclaim loudly, just one more time, "THIS STUFF WORKS FOR ALL CASES". The orientation of the surface and the direction of travel of the sphere are totally irrelevant. If we build a rectangular room with six surfaces, the sphere will bounce around inside the room as it should. If we tilt the room on 45 degrees, the sphere will STILL behave itself. Happy Day :)
It is more fun if you consider drag force (laminar) due to air acting on the sphere when undergoing free fall.
ma = W + kv
ma = mg + kv
a = g + kv/m
dv/dt = (mg + kv)/m
m/(mg+kv) dv/dt = 1
m/k ln = t if t=0,v=0
Then you would have a more interesting equation, one which can be solved by solving Differential equations.
Maybe I am spourting rubbish here. Please forgive me :-D
ma = W + kv
ma = mg + kv
a = g + kv/m
dv/dt = (mg + kv)/m
m/(mg+kv) dv/dt = 1
m/k ln = t if t=0,v=0
Then you would have a more interesting equation, one which can be solved by solving Differential equations.
Maybe I am spourting rubbish here. Please forgive me :-D
but wouldnt airfriction be unnesserary computional, below a certain speed?
for example a car, when hitting 70km/h, start calc airfriction?
for example a car, when hitting 70km/h, start calc airfriction?
Ultrano -
Resistance in the direction of travel due to Drag is easier to implement than that - all you do is multiply the current Velocity by a coefficient between 0 and 1, which is how I slow my Camera down.
However, should the density of the medium change (for example, an object collides with the surface of a lake) then things get more interesting.
We have to start treating objects as collections of particles and apply perparticle physics to obtain correct results in such conditions, although there's shortcuts for specific geometries, general cases require more math than specific cases.
Daydreamer -
Using the method I just described, the drag is dependant on the velocity - the faster you travel, the more drag you experience. Simple is good.
Resistance in the direction of travel due to Drag is easier to implement than that - all you do is multiply the current Velocity by a coefficient between 0 and 1, which is how I slow my Camera down.
However, should the density of the medium change (for example, an object collides with the surface of a lake) then things get more interesting.
We have to start treating objects as collections of particles and apply perparticle physics to obtain correct results in such conditions, although there's shortcuts for specific geometries, general cases require more math than specific cases.
Daydreamer -
Using the method I just described, the drag is dependant on the velocity - the faster you travel, the more drag you experience. Simple is good.
I'm getting some strange results from my existing code, and I'd appreciate some help.
Probably the single most important procedure I am calling is called PointAndPlane. It's a lot like "PlaneDistance" in that it compares a 3D point in space to a Plane and returns the distance between the point and the plane (positive = in front, negative = behind, 0 = on plane), but uses different formula to calculate it.
The formula is (X ? P) . N
where X is the point being tested, P is some point on the plane, and N is the surface normal of the plane.
This seems pretty straightforward, and I doubt there is any bug here.
Here is my code:
The Vec3Dot procedure leaves its result on the fpu, and thus the PointAndPlane procedure does as well.
I better post the code for Vec3Dot now.
Why am I posting this stuff?
Simple - I am getting obviously incorrect Distances from the above code.
This is screwing everything up. Look at this line taken from my debug logfile:
Frame: 294 Time = 4.8443883230969296 Seconds FPS = 60.68878
SphereCenter Distance to Plane = 1061.819 : Y=1.611078
I have a flat rectangular surface in XZ, at Y=0 height.
From the first of the two triangles I extract my SurfaceNormal.
When I test against this plane I give the first vertex of the first triangle as the "poiint on plane".
If my test point is roughly Y=1.6, how the hell can I be roughly 1062 units away from the plane? That's crazy.
I have another procedure for calculating planedistance as I mentioned, but its only valid in "world space", whereas this one is able to be applied in "object space" (not that I am doing this right now, it's just nice to have options, right?) so I'd really like to see this one working.
I can provide more information for example the calculation of SurfaceNormal and the vertices of the quad surface, but tell me, is there any obvious problem in what I have presented today?
Probably the single most important procedure I am calling is called PointAndPlane. It's a lot like "PlaneDistance" in that it compares a 3D point in space to a Plane and returns the distance between the point and the plane (positive = in front, negative = behind, 0 = on plane), but uses different formula to calculate it.
The formula is (X ? P) . N
where X is the point being tested, P is some point on the plane, and N is the surface normal of the plane.
This seems pretty straightforward, and I doubt there is any bug here.
Here is my code:
PointAndPlane proc pPointOnPlane, pNormalOfPlane, pPointToCheck
local lvec:Vec3
mov ebx, pPointOnPlane
mov ecx, pPointToCheck
fld [ecx].Vec3.X
fsub [ebx].Vec3.X
fstp lvec.X
fld [ecx].Vec3.Y
fsub [ebx].Vec3.Y
fstp lvec.Y
fld [ecx].Vec3.Z
fsub [ebx].Vec3.Z
fstp lvec.Z
invoke Vec3Dot,addr lvec,pNormalOfPlane
ret
PointAndPlane endp
The Vec3Dot procedure leaves its result on the fpu, and thus the PointAndPlane procedure does as well.
I better post the code for Vec3Dot now.
;This computes the dot product of 2 vectors
; ** WARNING ** The return value is on the FPU STACK
Vec3Dot proc uses esi edi pVector1, pVector2
; The dot product is this equation: V1.V2 = (V1.x * V2.x + V1.y * V2.y + V1.z * V2.z)
; In math terms, it looks like this: V1.V2 = ||V1|| ||V2|| cos(theta).
; It returns the distance of the projected vector, vVector2, from the start of vVector1.
mov esi,pVector1
mov edi,pVector2
fld [esi].Vec3.X
fmul [edi].Vec3.X
fld [esi].Vec3.Y
fmul [edi].Vec3.Y
fld [esi].Vec3.Z
fmul [edi].Vec3.Z
fadd
fadd
ret
Vec3Dot endp
Why am I posting this stuff?
Simple - I am getting obviously incorrect Distances from the above code.
This is screwing everything up. Look at this line taken from my debug logfile:
Frame: 294 Time = 4.8443883230969296 Seconds FPS = 60.68878
SphereCenter Distance to Plane = 1061.819 : Y=1.611078
I have a flat rectangular surface in XZ, at Y=0 height.
From the first of the two triangles I extract my SurfaceNormal.
When I test against this plane I give the first vertex of the first triangle as the "poiint on plane".
If my test point is roughly Y=1.6, how the hell can I be roughly 1062 units away from the plane? That's crazy.
I have another procedure for calculating planedistance as I mentioned, but its only valid in "world space", whereas this one is able to be applied in "object space" (not that I am doing this right now, it's just nice to have options, right?) so I'd really like to see this one working.
I can provide more information for example the calculation of SurfaceNormal and the vertices of the quad surface, but tell me, is there any obvious problem in what I have presented today?
Well, I can confirm that my code for calculating SurfaceNormals from Triangles is working fine. I suspected this already, because I have lighting enabled, and the surface is being lit appropriately.
Anyhow, here's another line from my LogFile:
SurfaceNormal: X-0.000000e+647 Y1 Z0
It really should be 0,1,0
but I think we can call -0.000000e+647 close enough to 0, yes?
So, the Floor surface is definitely pointing "upwards"...
For those who might be interested, here is the source for calculating SurfaceNormal from Triangle, along with all relevant math functions, all apparently working 100% :)
Anyhow, here's another line from my LogFile:
SurfaceNormal: X-0.000000e+647 Y1 Z0
It really should be 0,1,0
but I think we can call -0.000000e+647 close enough to 0, yes?
So, the Floor surface is definitely pointing "upwards"...
For those who might be interested, here is the source for calculating SurfaceNormal from Triangle, along with all relevant math functions, all apparently working 100% :)
SurfaceNormalFromTriangle proc pTriangle, pDest
local vVector1[3]:FLOAT
local vVector2[3]:FLOAT
local vNormal[3]:FLOAT
; Get 2 vectors from the polygon (2 sides), Remember the order!
; (Triangles are 3 floats, ie, 12 bytes each)
mov ebx,pTriangle
add ebx,(sizeof Vec3*2)
invoke VectorFrom2Points, ebx, pTriangle, addr vVector1http://www.asmcommunity.net/board/posting.php?mode=editpost&p=157683
mov ebx,pTriangle
add ebx,(sizeof Vec3)
invoke VectorFrom2Points, ebx,pTriangle, addr vVector2
; Take the cross product of our 2 vectors to get a perpendicular vector
invoke Vec3Cross, addr vVector2, addr vVector1, pDest
; Now we have a normal, but it's at a strange length, so let's make it length 1.
invoke Vec3Normalize, pDest
ret
SurfaceNormalFromTriangle endp
; This proc calculates the vector between 2 points
VectorFrom2Points proc uses ecx esi edi pPoint1, pPoint2, pDest
mov esi, pPoint1
mov edi,pPoint2
xor eax,eax
mov ecx,pDest
mov [ecx].Vec3.X,eax
mov [ecx].Vec3.Y,eax
mov [ecx].Vec3.Z,eax
; In order to get a vector from 2 points (a direction) we need to
; subtract the second point from the first point.
fld [esi].Vec3.X
fsub [edi].Vec3.X
fstp [ecx].Vec3.X
fld [esi].Vec3.Y
fsub [edi].Vec3.Y
fstp [ecx].Vec3.Y
fld [esi].Vec3.Z
fsub [edi].Vec3.Z
fstp [ecx].Vec3.Z
ret
VectorFrom2Points endp
;This proc calculates a perpendicular vector from 2 given vectors by taking the cross product.
;Returns in given buffer the cross product (Direction the polygon is facing - Normal)
Vec3Cross proc uses ecx esi edi pVector1, pVector2 , pDest
; If we are given 2 vectors (directions of 2 sides of a polygon)
; then we have a plane define. The cross product finds a vector that is perpendicular
; to that plane, which means it's point straight out of the plane at a 90 degree angle.
mov esi,pVector1
mov edi,pVector2
mov ecx,pDest
; The X value for the vector is: (V1.y * V2.z) - (V1.z * V2.y)
; vNormal.x = ((vVector1.y * vVector2.z) - (vVector1.z * vVector2.y));
fld [esi].Vec3.Y
fmul [edi].Vec3.Z
fld [esi].Vec3.Z
fmul [edi].Vec3.Y
fsub
fstp [ecx].Vec3.X
; The Y value for the vector is: (V1.z * V2.x) - (V1.x * V2.z)
; vNormal.y = ((vVector1.z * vVector2.x) - (vVector1.x * vVector2.z));
fld [esi].Vec3.Z
fmul [edi].Vec3.X
fld [esi].Vec3.X
fmul [edi].Vec3.Z
fsub
fstp [ecx].Vec3.Y
; The Z value for the vector is: (V1.x * V2.y) - (V1.y * V2.x)
; vNormal.z = ((vVector1.x * vVector2.y) - (vVector1.y * vVector2.x));
fld [esi].Vec3.X
fmul [edi].Vec3.Y
fld [esi].Vec3.Y
fmul [edi].Vec3.X
fsub
fstp [ecx].Vec3.Z
ret
Vec3Cross endp
; This returns a normalize vector (A vector exactly of length 1)
Vec3Normalize proc uses esi pNormal
local fmag
invoke Vec3Magnitude,pNormal ; Get the magnitude of our normal
fstp fmag
; Now that we have the magnitude, we can divide our normal by that magnitude.
; That will make our normal a total length of 1. This makes it easier to work with too.
mov esi,pNormal
fld [esi].Vec3.X
fdiv fmag
fstp [esi].Vec3.X
fld [esi].Vec3.Y
fdiv fmag
fstp [esi].Vec3.Y
fld [esi].Vec3.Z
fdiv fmag
fstp [esi].Vec3.Z
ret
Vec3Normalize endp
Daydreamer -
Using the method I just described, the drag is dependant on the velocity - the faster you travel, the more drag you experience. Simple is good.
Using the method I just described, the drag is dependant on the velocity - the faster you travel, the more drag you experience. Simple is good.
in cameraexample which is emulating drag in your eyemuscles and slow speeds
in example above with gravity speed=9.82^2 you can reach highspeeds vs as you describe it as linear drag increase approximation, while realworld air drag balances this falling with exponential increase
meaning parachutes wont work correctly, you get big craters from low altitudes of falling objects
and your ferrari is suddenly capable of 400mph instead of 200mph
planes get ridicoulous speeds
You probably missed the point - the drag coefficient has nothing to do with how velocity is attained - it just slows things down.
It's just a value somewhat less than 1 ( ie 0.8 ) which we use to reduce the speed of stuff (we multiply velocity by coefficient).
So if coefficient is 1.0, speed is unchanged.
If coefficient is 0, speed is zero.
If coefficient is something less than 1 and greater than 0, speed is reduced.
The faster something is travelling, the faster its speed will be reduced, simply because the numbers are bigger. This is a ratio thing.
I never desribed linearity.
Velocity is altered at each timestep by acceleration (summed forces).
Therefore velocity is not linear.
Actually, later, we'll apply two more such coefficients, called KR and KF.
KR is the "coefficient of restitution" which determines how "bouncy" a surface is. 1.0 means very bouncy (no loss), and 0.0 makes things get stuck on contact (vertical).
KF is the "coefficient of friction" and in all the examples I've seen, for some reason or other, it works in reverse, such that 1.0 makes things stick, and 0.0 lets things slide freely (horizontal).
It's just a value somewhat less than 1 ( ie 0.8 ) which we use to reduce the speed of stuff (we multiply velocity by coefficient).
So if coefficient is 1.0, speed is unchanged.
If coefficient is 0, speed is zero.
If coefficient is something less than 1 and greater than 0, speed is reduced.
The faster something is travelling, the faster its speed will be reduced, simply because the numbers are bigger. This is a ratio thing.
I never desribed linearity.
Velocity is altered at each timestep by acceleration (summed forces).
Therefore velocity is not linear.
Actually, later, we'll apply two more such coefficients, called KR and KF.
KR is the "coefficient of restitution" which determines how "bouncy" a surface is. 1.0 means very bouncy (no loss), and 0.0 makes things get stuck on contact (vertical).
KF is the "coefficient of friction" and in all the examples I've seen, for some reason or other, it works in reverse, such that 1.0 makes things stick, and 0.0 lets things slide freely (horizontal).
Rethinking the problem of finding the intersection of a sphere and a plane, I think it might end up faster all said and done to use a slightly different algorithm anyway.
I'm still miffed that I haven't got a reasonable answer as to why the code shown isn't right, but nonetheless:
There's an algorithm for finding the intersection of a plane and a ray.
A ray is two points in 3D space.
All the fancy algorithms I've looked at are ways of splitting time up imho.
All said and done, they test two points and a plane.
If theres an intersection, they go on to find the intersection point.
This algorithm finds the intersection point immediately:
From NeHe:
t= (Xn dot (PointOnPLANE - Raystart)) / (Xn dot Raydirection)
t represents the distance from the start until the intersection point along the direction of the ray. Therefore substituting t into the ray equation we can get the collision point. There are a few special cases though. If Xn dot Raydirection = 0 then these two vectors are perpendicular (ray runs parallel to plane) and there will be no collision. If t is negative the collision takes place behind the starting point of the ray along the opposite direction and again there is no intersection.
( don't sue please, we're learning stuff :) )
and the following C code:
int TestIntersionPlane(const Plane& plane,const TVector& position,const TVector& direction, double& lamda, TVector& pNormal)
{
double DotProduct=direction.dot(plane._Normal); // Dot Product Between Plane Normal And Ray Direction
double l2;
// Determine If Ray Parallel To Plane
if ((DotProduct<ZERO)&&(DotProduct>-ZERO))
return 0;
l2=(plane._Normal.dot(plane._Position-position))/DotProduct; // Find Distance To Collision Point
if (l2<-ZERO) // Test If Collision Behind Start
return 0;
pNormal=plane._Normal;
lamda=l2;
return 1;
}
We have to forgive the error due to our linear solution of a nonlinear problem, but here goes:
Our ray consists of the position of the sphere object's origin (extended by radius after time, see below), and the origin as it would be in the next timestep (ie, if influenced by summed forces, gravity, velocity, etc, in the next frame), and then extended by the length of the radius. This is so we are testing the point Tangent on the sphere at two moments, not its origin, which saves us some conditional logic work later, and prevents "cross frame error".
If the ray between those points intersects the plane, we have more work.
If not, we have nothing to worry about in the current frame.
More Work:
We have found an intersection point.
It indicates the point Tangent to the Sphere, ie, on its Boundary, due to we bothered to add the Radius into consideration already.
We now need to subtract the Radius back out.
I know this seems messy, but it beats all the conditional logic and cross-frame errors that can otherwise occur if we do it the traditional way.
So we have this tangent point, we know the vector that the sphere was moving in, all we need to do is reverse that vector, then scale it by enough (move the sphere back along the ray enough) to bring the sphere origin to "its radius" distance from the plane.
That all sounds a bit complex too, but it's really not.
I'll see if I can make some more crap art.
I'm still miffed that I haven't got a reasonable answer as to why the code shown isn't right, but nonetheless:
There's an algorithm for finding the intersection of a plane and a ray.
A ray is two points in 3D space.
All the fancy algorithms I've looked at are ways of splitting time up imho.
All said and done, they test two points and a plane.
If theres an intersection, they go on to find the intersection point.
This algorithm finds the intersection point immediately:
From NeHe:
t= (Xn dot (PointOnPLANE - Raystart)) / (Xn dot Raydirection)
t represents the distance from the start until the intersection point along the direction of the ray. Therefore substituting t into the ray equation we can get the collision point. There are a few special cases though. If Xn dot Raydirection = 0 then these two vectors are perpendicular (ray runs parallel to plane) and there will be no collision. If t is negative the collision takes place behind the starting point of the ray along the opposite direction and again there is no intersection.
( don't sue please, we're learning stuff :) )
and the following C code:
int TestIntersionPlane(const Plane& plane,const TVector& position,const TVector& direction, double& lamda, TVector& pNormal)
{
double DotProduct=direction.dot(plane._Normal); // Dot Product Between Plane Normal And Ray Direction
double l2;
// Determine If Ray Parallel To Plane
if ((DotProduct<ZERO)&&(DotProduct>-ZERO))
return 0;
l2=(plane._Normal.dot(plane._Position-position))/DotProduct; // Find Distance To Collision Point
if (l2<-ZERO) // Test If Collision Behind Start
return 0;
pNormal=plane._Normal;
lamda=l2;
return 1;
}
We have to forgive the error due to our linear solution of a nonlinear problem, but here goes:
Our ray consists of the position of the sphere object's origin (extended by radius after time, see below), and the origin as it would be in the next timestep (ie, if influenced by summed forces, gravity, velocity, etc, in the next frame), and then extended by the length of the radius. This is so we are testing the point Tangent on the sphere at two moments, not its origin, which saves us some conditional logic work later, and prevents "cross frame error".
If the ray between those points intersects the plane, we have more work.
If not, we have nothing to worry about in the current frame.
More Work:
We have found an intersection point.
It indicates the point Tangent to the Sphere, ie, on its Boundary, due to we bothered to add the Radius into consideration already.
We now need to subtract the Radius back out.
I know this seems messy, but it beats all the conditional logic and cross-frame errors that can otherwise occur if we do it the traditional way.
So we have this tangent point, we know the vector that the sphere was moving in, all we need to do is reverse that vector, then scale it by enough (move the sphere back along the ray enough) to bring the sphere origin to "its radius" distance from the plane.
That all sounds a bit complex too, but it's really not.
I'll see if I can make some more crap art.