Contact analysis is performed within each iteration until the equilibrium condition is satisfied. The residual vector that results from the global equilibrium is computed before the contact analysis starts. This vector is used to compute the contact forces resulting from the contact inequality constraints applied on the system. The constraint matrix is constructed and assembled to the left-hand-side coefficient matrix in each iteration. This would require the entire coefficient matrix to be reformulated and decomposed after each iteration, which is computationally expensive. Alternatively, using Crout factorization technique and isolating the Lagrange multiplier degrees of freedom at the end of the left-hand-side coefficient matrix, makes this action unnecessary. Only the part of the coefficient matrix corresponding to the Lagrange multiplier degrees of freedom is assembled and decomposed in each iteration without having to reconstruct and decompose the entire coefficient matrix. This technique saves a tremendous amount of computation.
The overall contact analysis procedure is performed, for each node on the contactor surface, as follows: