Skip to content

Multi-point constraints (MPCs)

Multi-point constraints (MPCs) enforce linear relations among selected degrees of freedom (DOFs) in a finite element model. In numgeo, constraints introduced through *Equation are homogeneous and are imposed by the penalty method. The involved DOFs may belong to different nodes, as for the general and EqualDOF equation types, or to the same node, as for *Equation, linear.

Penalty formulation

Consider \(n_c\) homogeneous linear constraint equations. The \(q\)th equation is written as 1

\[ g_q(\boldsymbol{d})=\boldsymbol{a}_q^T\boldsymbol{d}=0, \]

where \(\boldsymbol{d}\) is the global solution vector and \(\boldsymbol{a}_q\) contains the coefficients of the DOFs involved in constraint \(q\). Collecting all constraints gives

\[ \boldsymbol{z}(\boldsymbol{d})=\boldsymbol{A}\boldsymbol{d}=\boldsymbol{0}, \]

where row \(q\) of \(\boldsymbol{A}\) is \(\boldsymbol{a}_q^T\). The current *Equation implementation supports homogeneous equations only; no constant right-hand-side term is introduced.

The unconstrained total potential \(\Pi_0\) is augmented by one quadratic penalty contribution per constraint:

\[ \Pi(\boldsymbol{d}) = \Pi_0(\boldsymbol{d}) + \frac{1}{2}\sum_{q=1}^{n_c}\alpha_q g_q^2(\boldsymbol{d}), \]

where \(\alpha_q>0\) is the effective penalty parameter of constraint \(q\). Differentiation with respect to \(\boldsymbol{d}\) gives the penalty contribution to the residual,

\[ \boldsymbol{r}^{p} = \sum_{q=1}^{n_c}\alpha_q\boldsymbol{a}_q g_q = \sum_{q=1}^{n_c}\alpha_q\boldsymbol{a}_q\boldsymbol{a}_q^T\boldsymbol{d}, \]

and the corresponding tangent contribution,

\[ \boldsymbol{K}^{p} = \sum_{q=1}^{n_c}\alpha_q\boldsymbol{a}_q\boldsymbol{a}_q^T. \]

Equivalently, with \(\boldsymbol{\Lambda}=\operatorname{diag}(\alpha_1,\ldots,\alpha_{n_c})\),

\[ \boldsymbol{r}^{p}=\boldsymbol{A}^T\boldsymbol{\Lambda}\boldsymbol{A}\boldsymbol{d}, \qquad \boldsymbol{K}^{p}=\boldsymbol{A}^T\boldsymbol{\Lambda}\boldsymbol{A}. \]

For a linear equilibrium problem \(\boldsymbol{K}\boldsymbol{d}=\boldsymbol{f}\), the augmented system is therefore

\[ \left(\boldsymbol{K}+\boldsymbol{K}^{p}\right)\boldsymbol{d}=\boldsymbol{f}. \]

In a nonlinear implicit analysis, numgeo subtracts \(\boldsymbol{r}^{p}\) from the correction right-hand side and adds \(\boldsymbol{K}^{p}\) to the current tangent matrix. For an individual constraint \(g=\boldsymbol{a}^T\boldsymbol{d}\), the assembled terms are consequently

\[ \Delta\boldsymbol{b}=-\alpha\boldsymbol{a}g, \qquad \Delta\boldsymbol{K}=\alpha\boldsymbol{a}\boldsymbol{a}^T, \]

where \(\Delta\boldsymbol{b}\) denotes the contribution to the correction right-hand side.

Penalty scaling in numgeo

The value prescribed through *Penalty is a dimensionless multiplier \(\alpha_0\), not the effective penalty stiffness itself. For each generated constraint \(q\), numgeo evaluates

\[ K_{\mathrm{local},q} = \max\left( K_{\min}, \max_{i\in\mathcal{I}_q}\left|K_{ii}\right| \right), \]

and sets

\[ \alpha_q=\alpha_0 K_{\mathrm{local},q}. \]

Here, \(\mathcal{I}_q\) is the set of global DOFs involved in constraint \(q\), \(K_{ii}\) are their current diagonal tangent entries and \(K_{\min}=10^{-2}\) prevents a vanishing penalty where all relevant diagonal entries are zero. The default multiplier is

\[ \alpha_0=10^3. \]

Because the local tangent is re-evaluated during the solution process, the effective penalty parameter may change during a nonlinear analysis.

For EqualDOF and the general two-set equation, the coefficient vector is \(\boldsymbol{a}=[1,-1]^T\). For *Equation, linear, the coefficients are the factors entered by the user. numgeo does not normalise these factors. A common scaling \(\boldsymbol{a}\mapsto c\boldsymbol{a}\) leaves the exact equation \(\boldsymbol{a}^T\boldsymbol{d}=0\) unchanged, but scales both the penalty residual and the penalty tangent by \(c^2\). Coefficients should therefore be scaled consistently and, where practical, chosen with comparable magnitudes.

Set-based linear equations

For *Equation, linear, one node set \(\mathcal{S}\) and \(n_t\geq 2\) DOF/factor pairs are specified. For every node \(m\in\mathcal{S}\), numgeo generates the independent same-node constraint

\[ g^{(m)} = \sum_{k=1}^{n_t}a_k d_k^{(m)} =0. \]

Thus, a set containing \(n_s\) nodes produces \(n_s\) separate equations. Nodes in the set are not coupled to one another by this equation type. For example, the two-dimensional relation

\[ n_1u_1^{(m)}+n_2u_2^{(m)}=0 \]

suppresses the displacement component in direction \(\boldsymbol{n}=[n_1,n_2]^T\) at every node of the set. Direction cosines provide a well-scaled coefficient vector for this application.

Numerical considerations

The penalty method introduces no additional unknowns, but the following aspects must be considered:

  • For any finite \(\alpha_q\), the constraint is satisfied only approximately. Increasing \(\alpha_0\) generally reduces the constraint violation.
  • Excessively large penalty multipliers can make the tangent matrix ill-conditioned and impair nonlinear convergence or the accuracy of the linear solve.
  • The rank-one matrix \(\alpha_q\boldsymbol{a}_q\boldsymbol{a}_q^T\) couples every pair of DOFs appearing in an equation. Constraints between distant DOFs may therefore increase matrix fill-in and solution cost.
  • Redundant or mutually inconsistent constraints should be avoided because they can worsen conditioning and may prevent a meaningful solution.
  • Scaling all coefficients of a linear equation changes its effective numerical penalty. This scaling must be considered together with the selected value of \(\alpha_0\).

Dynamic analyses

In an implicit dynamic analysis, the equation is assembled for the displacement solution. Velocities and accelerations are obtained from the displacements through the selected time-integration operator rather than being constrained independently. For constant linear MPCs under geometrically linear kinematics, the derived velocities and accelerations inherit the same relation up to the penalty and solver tolerances.

In an explicit dynamic analysis, numgeo evaluates the penalty force \(-\alpha\boldsymbol{a}g\) from the current solution, while no tangent contribution is assembled.

  1. Robert D. Cook and others. Concepts and applications of finite element analysis. John Wiley & Sons, 2007. ↩