The Open Dynamics Engine (ODE) is a free, industrial quality library for simulating articulated rigid body dynamics. For example, it is good for simulating ground vehicles, legged creatures, and moving objects in VR environments. It is fast, flexible and robust, and it has built-in collision detection. ODE is being developed by Russell Smith.
If ``rigid body simulation'' does not make much sense to you, check out What is a Physics SDK?.
The is the user guide for the third prototype release of ODE (version 0.025). This release is to help continue testing the core algorithms and solidifying the API. There will probably only be a few (hopefully minor) changes to the API in the future.
ODE is best for simulating articulated rigid body structures. An articulated structure is created when rigid bodies of various shapes are connected together with joints of various kinds. Examples are ground vehicles (where the wheels are connected to the chassis) or legged creatures (where the legs are connected to the body).
ODE is designed to be used in interactive or real-time simulation. It is particularly good for simulating moving objects in changeable virtual reality environments. This is because it is fast, robust and stable, and the user has complete freedom to change the structure of the system even while the simulation is running.
ODE uses a highly stable integrator, so that the simulation errors should not grow out of control. The physical meaning of this is that the simulated system should not "explode" for no reason (believe me, this happens a lot with other simulators if you are not careful). ODE emphasizes speed and stability over physical accuracy.
ODE has hard contacts. This means that a special non-penetration constraint is used whenever two bodies collide. The alternative, used in many other simulators, is to use virtual springs to represent contacts. This is difficult to do right and extremely error-prone.
ODE has a built-in collision detection system. However you can ignore it and do your own collision detection if you want to. The current collision primitives are sphere, box and plane - more collision objects will come later, as well as fast identification of potentially intersecting objects.
Here are the features:
One current problem is that the mu parameter that controls Coulomb friction is not yet interpreted correctly - it specifies the maximum friction force that can be present at a contact, rather than specifying the maximum ratio of tangential to normal force. This will be fixed in the next release. [Note: I need lots more documentation about how the friction approximation works].
ODE is Copyright © 2001 Russell Smith.
This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation; either version 2.1 of the License, or (at your option) any later version.
This library is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.
Think you can help? Please write to me.
[Here is where I will write some background information about rigid body dynamics and simulation. But in the meantime, please refer to Baraff's excellent SIGGRAPH tutorial].
A rigid body has various properties from the point of view of the simulation. Some properties change over time:
The origin of this coordinate frame is the body's point of reference. Some values in ODE (vectors, matrices etc) are relative to the body coordinate frame, and others are relative to the global coordinate frame.
Note that the shape of a rigid body is not a dynamical properties (except insofar as it influences the various mass properties). It is only collision detection that cares about the detailed shape of the body.
The process of simulating the rigid body system through time is called integration. Each integration step advances the current time by a given step size, adjusting the state of all the rigid bodies for the new time value. There are two main issues to consider when working with any integrator:
Between each integrator step the user can call functions to apply forces to the rigid body. These forces are added to "force accumulators" in the rigid body object. When the next integrator step happens, the sum of all the applied forces will be used to push the body around. The forces accumulators are set to zero after each integrator step.
In real life a joint is something like a hinge, that is used to connect two objects. In ODE a joint is very similar: It is a relationship that is enforced between two bodies so that they can only have certain positions and orientations relative to each other. This relationship is called a constraint -- the words joint and constraint are often used interchangeably. The following picture shows three different constraint types:
The first is a ball and socket joint that constraints the ``ball'' of one body to be in the same location as the ``socket'' of another body. The second is a hinge joint that constraints the two parts of the hinge to be in the same location and to line up along the hinge axle. The third is a slider joint that constraints the ``piston'' and ``socket'' to line up, and additionally constraints the two bodies to have the same orientation.
Each time the integrator takes a step all the joints are allowed to apply constraint forces to the bodies they affect. These forces are calculated such that the bodies move in such a way to preserve all the joint relationships.
Each joint has a number of parameters controlling its geometry. An example is the position of the ball-and-socket point for a ball-and-socket joint. The functions to set joint parameters all take global coordinates, not body-relative coordinates. A consequence of this is that the rigid bodies that a joint connects must be positioned correctly before the joint is attached.
A joint group is a special container that holds joints in a world. Joints can be added to a group, and then when those joints are no longer needed the entire group of joints can be very quickly destroyed with one function call. However, individual joints in a group can not be destroyed before the entire group is emptied.
This is most useful with contact joints, which are added and remove from the world in groups every time step.
When a joint attaches two bodies, those bodies are required to have certain positions and orientations relative to each other. However, it is possible for the bodies to be in positions where the joint constraints are not met. This ``joint error'' can happen in two ways:
There is a mechanism to reduce joint error: during each simulation step each joint applies a special force to bring its bodies back into correct alignment. This force is controlled by the error reduction parameter (ERP), which has a value between 0 and 1.
The ERP specifies what proportion of the joint error will be fixed during the next simulation step. If ERP=0 then no correcting force is applied and the bodies will eventually drift apart as the simulation proceeds. If ERP=1 then the simulation will attempt to fix all joint error during the next time step. However, setting ERP=1 is not recommended, as the joint error will not be completely fixed due to various internal approximations. A value of ERP=0.1 to 0.8 is recommended (0.2 is the default).
A global ERP value can be set that affects most joints in the simulation. However some joints have local ERP values that control various aspects of the joint.
Traditionally the constraint equation for every joint has the form
where v is a velocity vector for the bodies involved, J is a ``Jacobian'' matrix with one row for every degree of freedom the joint removes from the system, and c is a right hand side vector. At the next time step, a vector lambda is calculated (of the same size as c) such that the forces applied to the bodies to preserve the joint constraint are
ODE adds a new twist. ODE's constraint equation has the form
where CFM is a square diagonal matrix. CFM mixes the resulting constraint force in with the constraint that produces it. By using CFM you can model tire slip, suspension, soft ground and other things.
Solving for lambda gives
Thus CFM simply adds to the diagonal of the original system matrix. This has the additional benefit of taking the system away from singularity and improving factorizer accuracy.
[write more...]
[There is a lot that needs to be written about collision handling.]
Collisions between bodies or between bodies and the static environment are handled as follows:
A typical simulation will proceed like this:
[Discuss the various methods and approximations that are used.]
The ODE library can be built to use either single or double precision floating point numbers. Single precision is faster and uses less memory, but the simulation will have more numerical error that can result in visible problems. You will get less accuracy and stability with single precision.
[must describe what factors influence accuracy and stability].
The floating point data type is dReal. Other commonly used types are dVector3, dVector4, dMatrix3, dMatrix4, dQuaternion.
There are various kinds of object that can be created:
All 3-vectors (x,y,z) supplied to ``set'' functions are given as individual x,y,z arguments.
All 3-vector result arguments to get() function are pointers to arrays of dReal.
Larger vectors are always supplied and returned as pointers to arrays of dReal.
All coordinates are in the global frame except where otherwise specified.
The ODE library is written in C++, but its public interface is made of simple C functions, not classes. Why is this?
The ODE library can be compiled in "debugging" or "release" mode. Debugging mode is slower, but function arguments are checked and many compile-time tests are done to ensure internal consistency. Release mode is faster, but no checking is done.
The world object is a container for rigid bodies and joints.
dWorldID dWorldCreate();
Create a new, empty world and return its ID number.
void dWorldDestroy (dWorldID);
Destroy a world and everything in it.
void dWorldStep (dWorldID, dReal stepsize);
Step the world.
Hmmm. dBodyGetRotation returns a 4x3 rotation matrix.
dJointGroupID dJointGroupCreate (int max_size);
Create a joint group.
max_size is in bytes.
Some joints, like hinge-2 need to be attached to two bodies to work.
A ball and socket joint looks like this:
void dJointSetBallAnchor (dJointID, dReal x, dReal y, dReal z);
Set the joint anchor point.
void dJointGetBallAnchor (dJointID, dVector3 result);
Get the joint anchor point.
A hinge joint looks like this:
When the hinge anchor or axis is set, the current position of the attached bodies is examined and that position will be the zero angle.
A slider joint looks like this:
void dJointSetSliderAxis (dJointID, dReal x, dReal y, dReal z);
Set the slider axis parameter.
void dJointGetSliderAxis (dJointID, dVector3 result);
Get the slider axis parameter.
When the axis is set, the current position of the attached bodies is examined and that position will be the zero position.
A hinge-2 joint looks like this:
The hinge-2 joint is the same as two hinges connected in series, with different hinge axes. An example, shown in the above picture is the steering wheel of a car, where one axis allows the wheel to be steered and the other axis allows the wheel to rotate.
The hinge-2 joint has an anchor point and two hinge axes. Axis 1 is specified relative to body 1 (this would be the steering axis if body 1 is the chassis). Axis 2 is specified relative to body 2 (this would be the wheel axis if body 2 is the wheel).
Axis 1 can have joint limits and a motor, axis 2 can only have a motor.
Axis 1 can function as a suspension axis, i.e. the constraint can be compressible along that axis.
When the anchor or axis is set, the current position of the attached bodies is examined and that position will be the zero angles.
The fixed joint maintains a fixed relative position and orientation between two bodies, or between a body and the static environment. Using this joint is almost never a good idea in practice, except when debugging. If you need two bodies to be glued together it is better to represent that as a single body.
Currently the fixed joint does not support a non-identity relative rotation between two bodies, it only supports a relative offset.
A contact joint looks like this:
The contact joint prevents body 1 and body 2 from inter-penetrating at the contact point. It does this by only allowing the bodies to have an ``outgoing'' velocity in the direction of the contact normal. Contact joints typically have a lifetime of one time step. They are created and deleted in response to collision detection.
Contact joints can simulate friction at the contact by applying special forces in the two friction directions that are perpendicular to the normal.
When a contact joint is created, a dContact structure must be supplied. This has the following definition:
struct dContact {
dSurfaceParameters surface;
dContactGeom geom;
dVector3 fdir1;
};
|
fdir1 is a "first friction direction" vector that defines a direction along which frictional force is applied. It must be of unit length and perpendicular to the contact normal (so it is typically tangential to the contact surface). It should only be defined if the dContactFDir1 flag is set in surface.mode. The "second friction direction" is a vector computed to be perpendicular to both the contact normal and fdir1.
surface is a substructure that is set by the user. Its members define the properties of the colliding surfaces. It has the following members:
| dContactMu2 | If not set, use mu for both friction directions. If set, use mu for friction direction 1, use mu2 for friction direction 2. |
| dContactFDir1 | If set, take fdir1 as friction direction 1, otherwise compute friction direction 1 to be perpendicular to the contact normal (in which case its resulting orientation is not defined). |
| dContactBounce | If set, the contact surface is bouncy, in other words the bodies will bounce off each other. The exact amount of bouncyness is controlled by the bounce parameter. |
| dContactSoftERP | If set, the error reduction parameter of the contact normal can be set with the soft_erp parameter. This is useful to make surfaces soft. |
| dContactSoftCFM | If set, the constraint force mixing parameter of the contact normal can be set with the soft_cfm parameter. This is useful to make surfaces soft. |
| dContactMotion1 | If set, the contact surface is assumed to be moving independently of the motion of the bodies. This is kind of like a conveyor belt running over the surface. When this flag is set, motion1 defines the surface velocity in friction direction 1. |
| dContactMotion2 | The same thing as above, but for friction direction 2. |
| dContactSlip1 | Slip, friction direction 1. [this needs more description]. |
| dContactSlip2 | Slip, friction direction 2. [this needs more description]. |
The joint geometry parameter setting functions should only be called after the joint has been attached to bodies, and those bodies have been correctly positioned, otherwise the joint may not be initialized correctly. If the joint is not already attached, these functions will do nothing.
For the parameter getting functions, if the system is out of alignment (i.e. there is some joint error) then the anchor/axis values will be correct with respect to body 1 only (or body 2 if you specified body 1 as zero in the dJointAttach() function).
The default anchor for all joints is (0,0,0). The default axis for all joints is (1,0,0).
When an axis is set it will be normalized to unit length. The adjusted axis is what the axis getting functions will return.
When measuring a joint angle or position, a value of zero corresponds to the initial position of the bodies relative to each other.
Note that there are no functions to set joint angles or positions (or their rates) directly, instead you must set the corresponding body positions and velocities.
When a joint is first created there is nothing to prevent it from moving through its entire range of motion. For example a hinge will be able to move through its entire angle, and a slider will slide to any length.
This range of motion can be limited by setting stops on the joint. The joint angle (or position) will be prevented from going below the low stop value, or from going above the high stop value. Note that a joint angle (or position) of zero corresponds to the initial body positions.
As well as stops, many joint types can have motors. A motor applies a torque (or force) to a joint's degree(s) of freedom to get it to pivot (or slide) at a desired speed. Motors have force limits, which means they can apply no more than a given maximum force/torque to the joint.
Motors have two parameters: a desired speed, and the maximum force that is available to reach that speed. This is a very simple model of real life motors, engines or servos. However, is it quite useful when modeling a motor (or engine or servo) that is geared down with a gearbox before being connected to the joint. Such devices are often controlled by setting a desired speed, and can only generate a maximum amount of power to achieve that speed (which corresponds to a certain amount of force available at the joint).
Motors can also be used to accurately model dry (or Coulomb) friction in joints. Simply set the desired velocity to zero and set the maximum force to some constant value - then all joint motion will be impeded by that force.
The alternative to using joint stops and motors is to simply apply forces to the affected bodies yourself. Applying motor forces is easy, and joint stops can be emulated with restraining spring forces. However applying forces directly is often not a good approach and can lead to severe stability problems if it is not done carefully.
Consider the case of applying a force to a body to achieve a desired velocity. To calculate this force you use information about the current velocity, something like this:
This has several problems. First, the parameter k must be tuned by hand. If it is too low the body will take a long time to come up to speed. If it is too high the simulation will become unstable. Second, even if k is chosen well the body will still take a few time steps to come up to speed. Third, if any other ``external'' forces are being applied to the body, the desired velocity may never even be reached (a more complicated force equation would be needed, which would have extra parameters and its own problems).
Joint motors solve all these problems: they bring the body up to speed in one time step, provided that does not take more force than is allowed. Joint motors need no extra parameters because they are actually implemented as constraints. They can effectively see one time step into the future to work out the correct force. This makes joint motors more computationally expensive than computing the forces yourself, but they are much more robust and stable, and far less time consuming to design with. This is especially true with larger rigid body systems.
Similar arguments apply to joint stops.
Here are the functions that set stop and motor parameters (as well as other kinds of parameters) on a joint:
| dParamLoStop | Low stop angle or position. Setting this to -dInfinity (the default value) turns off the low stop. For rotational joints, this stop must be greater than -Pi to be effective. |
| dParamHiStop | High stop angle or position. Setting this to dInfinity (the default value) turns off the high stop. For rotational joints, this stop must be less than Pi to be effective. If the high stop is less than the low stop then both stops will be ineffective. |
| dParamVel | Desired motor velocity (this will be an angular or linear velocity). |
| dParamFMax | The maximum force or torque that the motor will use to achieve the desired velocity. This must always be greater than or equal to zero. Setting this to zero (the default value) turns off the motor. |
| dParamFudgeFactor | The current joint stop/motor implementation has a small problem: when the joint is at one stop and the motor is set to move it away from the stop, too much force may be applied for one time step, causing a ``jumping'' motion. This fudge factor is used to scale this excess force. It should have a value between zero and one (the default value). If the jumping motion is too visible in a joint, the value can be reduced. Making this value too small can prevent the motor from being able to move the joint away from a stop. |
| dParamBounce | The bouncyness of the stops. This is a restitution parameter in the range 0..1. 0 means the stops are not bouncy at all, 1 means maximum bouncyness. |
| dParamStopERP | The error reduction parameter (ERP) used by the stops. |
| dParamStopCFM | The constraint force mixing (CFM) value used by the stops. Together with the ERP value this can be used to get spongy or soft stops. Note that this is intended for unpowered joints, it does not really work as expected when a powered joint reaches its limit. |
| dParamSuspensionERP | Suspension error reduction parameter (ERP). Currently this is only implemented on the hinge-2 joint. |
| dParamSuspensionCFM | Suspension constraint force mixing (CFM) value. Currently this is only implemented on the hinge-2 joint. |
If a particular parameter is not implemented by a given joint, setting it will have no effect.
These parameter names can be optionally followed by a digit (currently only 2) to indicate the second set of parameters, e.g. for the second axis in a hinge-2 joint.
Rigid body orientations are represented with quaternions. A quaternion is four numbers [cos(theta/2) sin(theta/2)*u] where theta is a rotation angle and `u' is a unit length rotation axis.
Every rigid body also has a 3x3 rotation matrix that is derived from the quaternion. The rotation matrix and the quaternion always match.
Some information about quaternions:
The following are utility functions for dealing with rotation matrices and quaternions.
void dRSetIdentity (dMatrix3 R);
Set R to the identity matrix (i.e. no rotation).
void dQSetIdentity (dQuaternion q);
Set q to the identity rotation (i.e. no rotation).
void dQtoR (const dQuaternion q, dMatrix3 R);
Convert quaternion q to rotation matrix R.
void dRtoQ (const dMatrix3 R, dQuaternion q);
Convert rotation matrix R to quaternion q.
The mass parameters of a rigid body are described by a dMass structure:
typedef struct dMass {
dReal mass; // total mass of the rigid body
dVector4 c; // center of gravity position in body frame (x,y,z)
dMatrix3 I; // 3x3 inertia tensor in body frame, about POR
} dMass;
|
The following functions operate on this structure:
void dMassSetZero (dMass *);
Set all the mass parameters to zero.
[There are quite a lot of these, but they're not standardized enough to document yet].
[Document these later].
This chapter describes the built-in collision detection system of ODE. Using ODE's collision detection is optional - an alternative collision detection system can be used as long as it can supply the right kinds of information to ODE.
If two bodies touch, or if a body touches a static feature in its environment, the contact is represented by one or more "contact points". Each contact point has a corresponding dContactGeom structure:
struct dContactGeom {
dVector3 pos; // contact position
dVector3 normal; // normal vector
dReal depth; // penetration depth
};
|
depth is the depth to which the two bodies inter-penetrate each other. If the depth is zero then the two bodies have a grazing contact, i.e. they "only just" touch. However, this is rare - the simulation is not perfectly accurate and will often step the bodies too far so that the depth is nonzero.
normal is a unit length vector that is, generally speaking, perpendicular to the contact surface.
The convention is that if body 1 is moved along the normal vector by a distance depth (or equivalently if body 2 is moved the same distance in the opposite direction) then the contact depth will be reduced to zero. This means that the normal vector points "in" to body 1.
In real life, contact between two bodies is a sophisticated thing. Representing contacts by contact points is only an approximation. Contact "patches" or "surfaces" might be more physically accurate, but representing these things in high speed simulation software is a challenge.
Each extra contact point added to the simulation will slow it down some more, so sometimes we are forced to ignore contact points in the interests of speed. For example, when two boxes collide many contact points may be needed to properly represent the geometry of the situation, but we may choose to keep only the best three. Thus we are piling approximation on top of approximation.
To use ODE's collision detection, geometry objects must be associated with the rigid bodies. A geometry object represents a rigid shape in space. Geometry objects are distinct from rigid bodies in that a geometry object has geometrical properties (size, shape, position and orientation) but no dynamical properties (such as velocity or mass).
Every geometry object is an instance of a class, such as sphere, plane, or box. You can define your own classes as well.
Every geometry object has a position vector and a 3*3 rotation matrix. If a geometry object is associated with a rigid body then its position and rotation is actually the position and rotation of that body. The point of reference for the standard classes usually corresponds to their centers of mass. This makes them particularly easy to connect to dynamics objects. If other points of reference are required, composite objects should be used to encapsulate the primitives.
dReal dGeomSphereGetRadius (dGeomID sphere);
Return the radius of the given sphere.
void dGeomSetPosition (dGeomID, dReal x, dReal y, dReal z);
void dGeomSetRotation (dGeomID, const dMatrix3 R);
Set the position vector and rotation matrix of the geometry object.
These functions are analogous to dBodySetPosition() and
dBodySetRotation().
If the geometry object is attached to a body, the body's position / rotation
will also be changed.
const dReal * dGeomGetPosition (dGeomID);
const dReal * dGeomGetRotation (dGeomID);
Return pointers to the geometry object's position vector and rotation matrix.
The returned values are pointers to internal data structures, so the vectors
are valid until any changes are made to the geometry object.
If the geometry object is attached to a body, the body's position / rotation
pointers will be returned, i.e. the result will be identical to calling
dBodyGetPosition() or dBodyGetRotation().
void dGeomDestroy (dGeomID);
Destroy a geometry object, removing it from the space first.
flags specifies how contacts should be generated if the objects touch. Currently the lower 8 bits of flags specifies the maximum number of contact points to generate. If this number is zero, this function just pretends that it is one - in other words you can not ask for zero contacts. All other bits in flags must be zero. In the future the other bits may be used to select other contact generation strategies.
contact must point to an array of contact geometry information that can hold at least the maximum number of contacts. Note: the elements of the contact array do not necessarily have to be contiguous. skip is the number of bytes between each dContactGeom structure in the contact array. If skip is sizeof(dContactGeom) then contact points to a "normal" (C-style) contact array. If skip is larger than this, then the dContactGeom structures are embedded in some other larger structures. It is an error for skip to be smaller than sizeof(dContactGeom).
If the objects touch, this returns the number of contact points generated (and updates the contact array), otherwise it returns 0 (and the contact array is not touched).
You can define your own geometry classes using the functions in this section. The standard geometry classes do not have any special access to the internals of ODE, they use the public functions exactly as you would.
Every geometry class has a unique number (0,1,2, etc...). A new geometry class (call it `X') must provide the following to ODE:
typedef int dColliderFn (dGeomID o1, dGeomID o2, int flags,
dContactGeom *contact, int skip);
|
typedef dColliderFn * dGetColliderFnFn (int num); |
This function is called infrequently - the return values are cached and reused.
typedef void dGetAABBFn (dGeomID g, dReal aabb[6]); |
Here are the functions used to manage custom classes:
Here is the definition of the dGeomClass structure:
struct dGeomClass {
int bytes; // bytes of custom data needed
dGetColliderFnFn *collider; // collider function
dGetAABBFn *aabb; // bounding box function
};
|
dGeomID dCreateGeom (int classnum);
Create a geometry object of the given class number.
The custom data block will initially be set to 0.
This object can be added to a space using dSpaceAdd().
int dGeomGetClass (dGeomID);
Given a geometry object, this returns its class number.
When you implement a new class you will usually write a function that does the following:
A space is a container for geometry objects. It is similar to the rigid body "world", except for collision instead of dynamics.
The space does high level collision culling, which means that it can identify which pairs of geometry objects are potentially touching. You can safely call dCollide() for only those pairs, instead of having to call dCollide() for every object-object pair. This can save a huge amount of time. Various collision culling algorithms will be available, however only the first is currently implemented:
Here are the functions used for spaces:
dSpaceID dSpaceCreate();
Create a space.
void dSpaceCollide (dSpaceID space, void *data,
dNearCallback *callback);
Call a callback function one or more times, for all potentially intersecting
objects in the space.
The callback function is of type dNearCallback, which is defined as:
The data variable is passed from dSpaceCollide() to the
callback function. Its meaning is user defined.
The o1 and o2 arguments are the geometry objects that may be near
each other.
The callback function can call dCollide() on o1 and o2,
perhaps first determining whether to collide them at all based on other
information.
typedef void dNearCallback (void *data, dGeomID o1, dGeomID o2);
Future interface enhancements:
[just notes for now]
What factors does execution speed depend on? Each joint removes a number of degrees of freedom (DOFs) from the system. For example the ball and socket removes three, and the hinge removes five. For each separate group of bodies connected by joints, where:
ODE currently relies on factorization of a ``system'' matrix that has one row/column for each DOF removed (this is where the O(m23) comes from). In a 10 body chain that uses ball and socket joints, roughly 30-40% of the time is spent filling in this matrix, and 30-40% of the time is spent factorizing it.
Thus, to speed up your simulation you might consider:
Use dJointAttach with arguments (body,0) or (0,body).
No. ODE is a computational engine, and is completely independent of any graphics library. However the examples that come with ODE use OpenGL, and most interesting uses of ODE will need some graphics library to make the simulation visible to the user. But that's your problem.
Sometimes when rigid bodies collide without restitution, they appear to inter-penetrate slightly and then get pushed apart so that they only just touch. The problem gets worse as the time step gets larger. What is going on?
The contact joint constraint is only applied after the collision is detected. If a fixed time step is being used, it is likely that the bodies have already penetrated when this happens. The error reduction mechanism will push the bodies apart, but this can take a few time steps (depending on the value of the ERP parameter).
This penetration and pushing apart sometimes makes the bodies look like they are bouncing, although it is completely independent of whether restitution is on or not.
Some other simulators have individual rigid bodies take variable sized timesteps to make sure bodies never penetrate much. However ODE takes fixed size steps, as automatically choosing a non-penetrating step size is problematic for an articulated rigid body simulator (the entire ARB structure must be stepped to account for the first penetration, which may result in very small steps).
There are three fixes for this problem:
In other words, how can you create a body that doesn't move, but that interacts with other bodies? The answer is to create a geometry object only, without the corresponding rigid body object. The geometry object is associated with a rigid body ID of zero. Then in the contact callback when you detect a collision between two geometry objects with a nonzero body ID and a zero body ID, you can simply pass those two IDs to the dJointAttach() function as normal. This will create a contact between the rigid body and the static environment.
Don't try to get the same effect by setting a very high mass/inertia on the ``motionless'' body and then resetting it's position/orientation on each time step. This can cause unexpected simulation errors.
[only notes for now]
Why don't I implement a proper friction pyramid or friction cone (e.g. Baraff's version) ?} Because I have to factor non-symmetric (and possibly indefinite) matrices, for either static or dynamic friction. Speed was considered more important - the current friction approximation only needs a symmetric factorization, which is twice as fast.
Matrix operations like factorization are expensive, so we must store the data in a way that is most useful to the matrix code. I want to do 4-way SIMD optimizations later, so the format is this: store the matrix by rows, and each row is rounded up to a multiple of 4 elements. The extra "padding" elements at the end of each row/column must be set to 0. This is called the "standard format". Hopefully this decision will remain good in the future, as more and more processors have 4-way SIMD (especially for fast 3D graphics).
The exception: matrices that have only one column or row (vectors), are always stored as consecutive elements in standard row format, i.e. there is no interior padding, only padding at the end.
Thus: all 3x1 floating point vectors are stored as 4x1 vectors: (x,x,x,0).
The dx prefix is used for internal structures that should never be visible externally. The d prefix is used for structures that are part of the public interface.