The principle of virtual work is used to derive the Euler-Lagrange equations of motion in order to describe the dynamics of multibody android systems. The constrained variational equations are in fact differential-algebraic equations of high index and are cast as ordinary differential equations through differentiation of the constraint equations. The integration routine LSODAR and the fourth-order Runge-Kutta method are used to compute the generalized coordinates, their time derivatives and the body forces of two android models. The graphs of the constraint forces reveal the whiplash effect on the neck and that the stiffness of both multibody systems is due to large magnitude impulsive forces experienced by many bodies simultaneously.