After initialization, the center-of-mass variables are integrated using the three-step Verlet method:
Friction terms are embedded in the force 
.
The integration of rotational variables uses a different central difference
expansion. The algorithm begins with the advancement of the Euler-Cayley
parameters by 
, using the following four steps:
Here 
represents the angular momentum of a rigid body transformed
to the body frame, and 
represents the transformed
angular velocity; q is the four-vector (
,
,
,
).
Finally, the torques are determined by
Again, all frictional terms are embedded in the force 
.
These half-step calculations are necessary for the final, full-step update:
The velocities of individual atoms, which are necessary for the friction evaluation, are determined using the angular velocities of the rigid bodies. The velocity due to rotation alone is determined in the body frame, then transformed into the lab frame and added to the translational velocity of the center of mass:
Finally, the atom positions are updated for the next energy evaluation using the new center-of-mass positions and the new Euler-Cayley parameters:
The rigid-body dynamics routine can be run in one of three modes: FREE, LANGevin, or TCOUpling. The FREE mode simply allows frictionless dynamics, the LANGevin mode adds friction- and bath-dependent random forces, and TCOUpling uses the method of Berendsen et al. (1984). The coordinates (and velocities) are affected by using a nonzero friction coefficient in the Langevin dynamics algorithm with zero random forces. The friction coefficient is computed by the program from the equation
where 
is given by the FBETa atom property.
Note that FBETa is by default zero; that is, FBETa has to be
specified in order to use TCOUpling.
The target temperature 
is specified by the TBATh parameter.
If the simulation involves light rigid groups such as water molecules, it is recommended that a time step on the order of .25 fsec be used. For heavier groups, time steps as large as 20 or 30 fsec are adequate.