Rigid body motion handled by Lagrange multipliers

I consider periodic 3D linear elasticity with only Neumann boundary conditions applied. Such a problem is not well defined as the body can freely move in three directions. I want to constrain these three DOFs. I implemented the penalty method, but it alters the results. I want to use Lagrange multipliers, but I am unsure what FE space I should choose because the multipliers are, in fact, only three global DOFs.

Shall I use VectorL2 space of order 0 for multipliers and periodic VectorH1 of order 2 for displacement? I tried this combination and it somehow gives a stiffer solution.


Hi Petr,

[with NGSolve I have so far only used the first approach below, so please take my suggestions rather as a possible guideline, which might or might not work because I could have overlooked some detail/aspect. I would be interested to hear about your success with any of the other method, if you try them…]

From my (engineering) experience, you only need a rather small penalty parameter to make your problem well-defined.
If the parameter is chosen small enough, the resulting perturbation of the solution can be negligible for practical purposes. I am not sure though whether this works well with iterative solvers. You might do a convergence study for decreasing the penalty parameter to get “a feeling” for that.

As a second options, you could try to let the FESpace parameter
“dirichlet_bbbnd” do the trick. It allows you to set Dirichlet BCs on nodes in 3d. How to name the corners of your geometry depends on the geometry backend (CSG, OCC…) that you use.

As a third option, you could manually create a “bitarray” based on fes.FreeDofs(), where you switch the bits corresponding to the displacements of three nodes to be constrained.
I hope this works also in the context of periodic dofs.

A fourth option is to constrain all 6 rigid body modes with one Lagrange multiplier per mode. The space for each multiplier would be a NumberFESpace. This is, however, without any linear algebra tricks, probably not the best approach performance-wise.


thanks for advising me.
Penalty method works suffitiently with small penalty parameter (relatively to Young’s modulus of elasticity). I tried this FE space combination:

V = ng.Periodic(ng.VectorH1(mesh, order=1) # for displacement
Q = ng.VectorValued(ng.NumberSpace(mesh, order=0)) # lagrange multiplier for translational rigid mode
fe_space = V * Q

It seemed to work. I will investigate it more closely w.r.t convergence. I found no explanation what ‘NumbefSpace’ does actually - is it a sort of a non-local FE space ? Or could you point me where I can find relevant information?


Hi Petr,

the NumberFESpace contains “only a number” for enforcing an equation. Therefore the system matrix will have one additional row/column, which is “dense” as the number couples with all other involved shape functions.

In the source code the implementation is in comp/numberfespace.hpp (-.cpp) if you are interested in the implementation details.