An outline of our semi-implicit simulation algorithm is as follows:
Initial calculations at t=0:
, 
and 
.
from
column heights and the index array MAX .
For each time step 
:
and the effective load
matrix 
as prescribed in Section 
.
, where 
is diagonal
and 
is lower triangular.
:
for 
using forward-substitution.
for using
back-substitution. Note that since is diagonal, the
computation of its inverse 
is trivial.
.
We have verified experimentally that the semi-implicit algorithm is
sufficiently stable for our simulation purposes. We use a time step
of 
in our simulator (the stability limit is reached
with a time step around 
). Our algorithm is much
more stable than an explicit Euler time integration
implementation. Given the same numerical settings, experiments show
that the largest time step allowed (before the system becomes
unstable) when using our algorithm is typically 100 to 150 times
larger than is possible when using the explicit method. On the other
hand, our algorithm is simpler than implicit methods described by
Bathe and Wilson Bathe76, such as the Houbolt method or
the Newmark method. In particular, the approximations we use for 
(Eq. ) and 
(Eq. ) are simple and require only knowledge of the
position matrix 
and the velocity matrix 
at
time t while other methods (for example the Houbolt method) requires
positions and velocities from several previous time steps.
Furthermore, fully implicit methods require that the external load
matrix 
in Eq. () be calculated
implicitly as well, which would be difficult in our case.
| Xiaoyuan Tu | January 1996 |