Constraint equations arise in the dynamics of mechanical systems whenever there is the need to restrict kinematically possible motions of the system. In practical applications, constraint equations can be used to simulate complex, connected systems. If the simulation must be carried out numerically, it is useful to look for a formulation that leads straightforwardly to a numerical approximation. In this paper, we extend the methodology of our previous work to incorporate the dynamics of holonomically and nonholonomically constrained systems. The constraint equations are cast in a variatonal form, which may be included easily, in the time finite element framework. The development of the weak constraint equations and their associated "tangent" operators is presented. We also show that this approach to constraint equations may be employed to develop time finite elements using a quaternion parametrization of finite rotation. Familiarity with the notation and methodology of our previously presented work is assumed.