When a multibody system reaches a singular position, one or more degrees of freedom appear instantaneously and the jacobian matrix of the constraint equations becomes rank-deficient. The classical kinematic formulation is based on the factorization of the jacobian and, therefore, fails in singular positions. In this paper we develop an efficient method, which uses a penalty and an augmented Lagrangian formulation, and successfully handles singular positions. This formulation automatically copes with redundant incompatible constraints and guarantees the stability of the constraints during numerical integration. Critical numerical examples are shown which corroborate these findings.