Friction
When two bodies are in contact with each other, they can produce a frictional force.
Friction force is calculated by the normal force between the two bodies multiplied by a coefficient of friction. The coefficient may be static or kinetic depending on the relative speed between the bodies.
If the bodies are stationary to each other, then the friction is "static"; if the bodies are sliding across each other, then the friction is "kinetic". The static coefficient is usually higher than the kinetic coefficient (it's harder to get something moving than to keep something moving).
The frictional force is opposite of the direction of travel. and the friction doesn't exceed the amount of force needed to keep the bodies locked together.
https://en.wikipedia.org/wiki/Friction
Normal and Surface Vector
The actual friction equation is very simple, but you need vectors normal and parallel to the surface to complete the calculations.
Since everything in my simulation is still a sphere, those values are pretty straightforward to determine.
Taking the unit vector of the relative locations between the two bodies gives me a unit normal vector. (Only for two spheres in contact.)
I came up with a couple ways to get the velocity component along the surface. (I'm sure there are many more!) I'm fairly happy with my second approach (please feel free to suggest alternatives!).
After getting the normal vector, I use the dot product of the relative velocity and the normal vector to get the speed of the bodies into each other. The resulting speed is a scalar (magnitude only, no direction).
I can multiple that speed scalar with the normal vector to get a velocity vector normal to contact.
Finally I can subtract the total relative velocity by the normal relative velocity to get the relative velocity along the surface.
Take the unit vector of that, and we get the motion of travel along the surface!
With the unit normal and surface vectors, we can finish the friction calculations.
Take the dot product of the relative force and the unit normal vector to get the normal force.
The magnitude of the friction force is the normal force * the coefficient of friction. And then multiply that value by the unit vector along the surface to get the friction force vector.
Finally add that frictional force vector to the sum of all the other forces acting upon that body.
Orbital Aero Code - update()
To recap from a previous post, there is an update() function in the OrbitalAeroModel class. I've made some cleanup and optimization changes to it compared to previous posts. The new update function combines Indirect and Direct forces into a single function call.
// Primary function to process the physics and kinematics all bodies.
bool OrbitalAeroModel::update(double timeStep_seconds)
{
// It's important all forces are calculated before all positions are updated before all collisions
// are checked, etc. That way the order in which bodies are looped through does not affect how they
// affect other bodies in space.
// Calculate Motor Thrust (From control system and motor abilities.)
// Sum all the forces acting on every body, including forces from thrust, gravity, atmospheric drag,
// and friction.
processForces();
// Can now translate and rotate each body.
processKinematics(timeStep_seconds);
// After moving the body, check for collisions and adjust positions/velocities as appropriate.
processCollisions(timeStep_seconds);
return true;
}
Orbital Aero Code - processForces()
The previous processIndirectForces function has been renamed processForces and it now includes the frictional forces.
// Function to loop through all bodies and sum the forces acting on them, including forces from
// thrust, gravity, atmospheric drag, and friction.
bool OrbitalAeroModel::processForces()
{
// Loop through all bodies and calculate their forces upon each other.
// The outer loop is for the primary body that we're calculating forces for.
map <uint64, BodyModel>::iterator iteratorPrimaryBody = tableBodies.begin();
while (iteratorPrimaryBody != tableBodies.end())
{
BodyModel &primaryBody = iteratorPrimaryBody->second;
// Only process the primary body if it's active.
if (primaryBody.active)
{
// Create a vector to contain the sum of all forces.
Vector sumForces_newtons;
// Start with the thrust force of any onboard engines.
sumForces_newtons = primaryBody.forceThrust_newtons;
// The inner loop is for all other bodies that influence the primary.
map <uint64, BodyModel>::iterator iteratorSecondaryBody = tableBodies.begin();
while (iteratorSecondaryBody != tableBodies.end())
{
BodyModel &secondaryBody = iteratorSecondaryBody->second;
// Only contribute the secondary body's influence if it's active and not the primary body.
if (secondaryBody.active && (&primaryBody != &secondaryBody))
{
// Get the relative distance between the two bodies.
Vector relativeDistance_meters = secondaryBody.location_meters - primaryBody.location_meters;
double relativeDistanceMagnitude_meters = relativeDistance_meters.magnitude();
if (relativeDistanceMagnitude_meters > 0.0)
{
.
// Are we within the atmosphere of the secondary body?
.
.
All the gravity and atmospheric drag code has remained unchanged. The friction code comes into play next...
// Are the bodies in contact?
double combinedRadii_meters = primaryBody.radius_meters + secondaryBody.radius_meters;
if (relativeDistanceMagnitude_meters < (combinedRadii_meters + Constants::relativeDistanceForContact_meters))
{
// Apply friction forces.
// Friction force is determined by the friction coefficient * the normal force between the bodies.
// Getting the normal vector for spheres is simple, we can just take the unit vector
// of our relative distance.
Vector vectorNormalToContact = relativeDistance_meters.unit();
// We also need a vector along the surface slope in our direction of travel to know the
// direction the friction force is acting in. We can use the relative velocity between
// the bodies and the normal vector to get the surface vector.
// Get the relative velocity between the two bodies.
Vector relativeVelocity_metersPerSecond = secondaryBody.velocity_metersPerSecond - primaryBody.velocity_metersPerSecond;
// Calculate the magnitude of velocity going normal to the surface.
double relativeSpeedNormalToContact_metersPerSecond = Vector::dotProduct(relativeVelocity_metersPerSecond, vectorNormalToContact);
// Calculate the velocity vector that is normal to the surface.
Vector relativeVelocityNormalToContact_metersPerSecond = vectorNormalToContact * relativeSpeedNormalToContact_metersPerSecond;
// Now that we know the full relative velocity and the velocity component normal
// to the surface, we can subtract those two to get the velocity component along
// the surface. This vector is opposite of the direction of travel.
Vector relativeVelocityAlongSurface_metersPerSecond = relativeVelocity_metersPerSecond - relativeVelocityNormalToContact_metersPerSecond;
// And the unit vector of our surface velocity gives us the surface direction.
Vector vectorAlongSurface = relativeVelocityAlongSurface_metersPerSecond.unit();
// Determine coefficient of friction based on the speed along surface.
double relativeSurfaceSpeed_meters = relativeVelocityAlongSurface_metersPerSecond.magnitude();
double coefficientOfFriction = 0.0;
if (abs(relativeSurfaceSpeed_meters) > Constants::relativeSpeedForKineticFriction_metersPerSecond)
{
// Body is moving, kinetic friction.
coefficientOfFriction = Constants::defaultCoefficientKineticFriction;
}
else
{
// Body is not moving, static friction.
coefficientOfFriction = Constants::defaultCoefficientStaticFriction;
}
// Friction force is based on the force pushing the two bodies together.
// Gravitational force pulling them together plus the thrust of both bodies.
// It is not necessary to take into account the gravitational forces of all the other
// nearby bodies, as they will have an equal influence on both bodies.
Vector relativeForce_newtons = secondaryGraviationalForce_newtons + primaryBody.forceThrust_newtons - secondaryBody.forceThrust_newtons;
// Get the force magnitude by taking the dot product of the normal vector with the relative force.
double normalForceMagnitude_newtons = Vector::dotProduct(relativeForce_newtons, vectorNormalToContact);
// Friction Force = coefficientOfFriction * normalForce
double frictionForceMagnitude_newtons = normalForceMagnitude_newtons * coefficientOfFriction;
// Friction force is applied in the opposite direction of movement along the surface.
Vector frictionForce_newtons = vectorAlongSurface * frictionForceMagnitude_newtons;
// Apply friction force to the sum.
sumForces_newtons += frictionForce_newtons;
// Frictional force should be limited to the amount of force necessary to bring the body to a stop.
// As implemented, this will cycle between continuously between the direction of travel to apply
// friction. (It never allows the object to come to a complete stop.)
// Friction should also cause the bodies to spin.
} // (relativeDistanceMagnitude_meters < (combinedRadii_meters + Constants::relativeDistanceForContact_meters))
} // relativeDistanceMagnitude_meters > 0.0
} // ((secondaryBody.active) && (primaryBody.index != secondaryBody.index))
iteratorSecondaryBody++;
} // (iteratorSecondaryBody != tableBodies.end())
primaryBody.forceIndirect_newtons = sumForces_newtons;
} // (primaryBody.active)
iteratorPrimaryBody++;
} // (iteratorPrimaryBody != tableBodies.end())
return true;
}
// Function to loop through all bodies and sum the forces acting on them, including forces from
// thrust, gravity, atmospheric drag, and friction.
bool OrbitalAeroModel::processForces()
{
// Loop through all bodies and calculate their forces upon each other.
// The outer loop is for the primary body that we're calculating forces for.
map <uint64, BodyModel>::iterator iteratorPrimaryBody = tableBodies.begin();
while (iteratorPrimaryBody != tableBodies.end())
{
BodyModel &primaryBody = iteratorPrimaryBody->second;
// Only process the primary body if it's active.
if (primaryBody.active)
{
// Create a vector to contain the sum of all forces.
Vector sumForces_newtons;
// Start with the thrust force of any onboard engines.
sumForces_newtons = primaryBody.forceThrust_newtons;
// The inner loop is for all other bodies that influence the primary.
map <uint64, BodyModel>::iterator iteratorSecondaryBody = tableBodies.begin();
while (iteratorSecondaryBody != tableBodies.end())
{
BodyModel &secondaryBody = iteratorSecondaryBody->second;
// Only contribute the secondary body's influence if it's active and not the primary body.
if (secondaryBody.active && (&primaryBody != &secondaryBody))
{
// Get the relative distance between the two bodies.
Vector relativeDistance_meters = secondaryBody.location_meters - primaryBody.location_meters;
double relativeDistanceMagnitude_meters = relativeDistance_meters.magnitude();
if (relativeDistanceMagnitude_meters > 0.0)
{
// Calculate the force of gravity.
..
// Are we within the atmosphere of the secondary body?
.
.
// Are the bodies in contact?
double combinedRadii_meters = primaryBody.radius_meters + secondaryBody.radius_meters;
if (relativeDistanceMagnitude_meters < (combinedRadii_meters + Constants::relativeDistanceForContact_meters))
{
// Apply friction forces.
// Friction force is determined by the friction coefficient * the normal force between the bodies.
// Getting the normal vector for spheres is simple, we can just take the unit vector
// of our relative distance.
Vector vectorNormalToContact = relativeDistance_meters.unit();
// We also need a vector along the surface slope in our direction of travel to know the
// direction the friction force is acting in. We can use the relative velocity between
// the bodies and the normal vector to get the surface vector.
// Get the relative velocity between the two bodies.
Vector relativeVelocity_metersPerSecond = secondaryBody.velocity_metersPerSecond - primaryBody.velocity_metersPerSecond;
// Calculate the magnitude of velocity going normal to the surface.
double relativeSpeedNormalToContact_metersPerSecond = Vector::dotProduct(relativeVelocity_metersPerSecond, vectorNormalToContact);
// Calculate the velocity vector that is normal to the surface.
Vector relativeVelocityNormalToContact_metersPerSecond = vectorNormalToContact * relativeSpeedNormalToContact_metersPerSecond;
// Now that we know the full relative velocity and the velocity component normal
// to the surface, we can subtract those two to get the velocity component along
// the surface. This vector is opposite of the direction of travel.
Vector relativeVelocityAlongSurface_metersPerSecond = relativeVelocity_metersPerSecond - relativeVelocityNormalToContact_metersPerSecond;
// And the unit vector of our surface velocity gives us the surface direction.
Vector vectorAlongSurface = relativeVelocityAlongSurface_metersPerSecond.unit();
// Determine coefficient of friction based on the speed along surface.
double relativeSurfaceSpeed_meters = relativeVelocityAlongSurface_metersPerSecond.magnitude();
double coefficientOfFriction = 0.0;
if (abs(relativeSurfaceSpeed_meters) > Constants::relativeSpeedForKineticFriction_metersPerSecond)
{
// Body is moving, kinetic friction.
coefficientOfFriction = Constants::defaultCoefficientKineticFriction;
}
else
{
// Body is not moving, static friction.
coefficientOfFriction = Constants::defaultCoefficientStaticFriction;
}
// Friction force is based on the force pushing the two bodies together.
// Gravitational force pulling them together plus the thrust of both bodies.
// It is not necessary to take into account the gravitational forces of all the other
// nearby bodies, as they will have an equal influence on both bodies.
Vector relativeForce_newtons = secondaryGraviationalForce_newtons + primaryBody.forceThrust_newtons - secondaryBody.forceThrust_newtons;
// Get the force magnitude by taking the dot product of the normal vector with the relative force.
double normalForceMagnitude_newtons = Vector::dotProduct(relativeForce_newtons, vectorNormalToContact);
// Friction Force = coefficientOfFriction * normalForce
double frictionForceMagnitude_newtons = normalForceMagnitude_newtons * coefficientOfFriction;
// Friction force is applied in the opposite direction of movement along the surface.
Vector frictionForce_newtons = vectorAlongSurface * frictionForceMagnitude_newtons;
// Apply friction force to the sum.
sumForces_newtons += frictionForce_newtons;
// Frictional force should be limited to the amount of force necessary to bring the body to a stop.
// As implemented, this will cycle between continuously between the direction of travel to apply
// friction. (It never allows the object to come to a complete stop.)
// Friction should also cause the bodies to spin.
} // (relativeDistanceMagnitude_meters < (combinedRadii_meters + Constants::relativeDistanceForContact_meters))
} // relativeDistanceMagnitude_meters > 0.0
} // ((secondaryBody.active) && (primaryBody.index != secondaryBody.index))
iteratorSecondaryBody++;
} // (iteratorSecondaryBody != tableBodies.end())
primaryBody.forceIndirect_newtons = sumForces_newtons;
} // (primaryBody.active)
iteratorPrimaryBody++;
} // (iteratorPrimaryBody != tableBodies.end())
return true;
}
Results
In this example, a body is just above the Earth falling towards it. It does a short couple bounces (I'll cover that code later), and then friction brings it to a stop.
I also tried adding a thrust force down when the body was in contact with the ground; stops sooner.
When it comes to the friction code, it could be greatly improved upon in the future. Bodies should have different properties (smooth vs rough), friction should cause bodies to spin, and I'm applying the full static force even when stopped (their friction forces constantly twitch back and forth when stationary). For my purposes now it gives nice results for when objects scrape along the ground.
Copyright (c) 2016 Clinton Kam
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
Copyright (c) 2016 Clinton Kam
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
