Coupling matrix for (frictional) contact problems

Dear all,

I would like to compute/access the coupling matrix that arises in the mixed formulation of contact problems, possibly with friction. By mixed I mean, of course, that the unknown is the couple (\mathbf u, \boldsymbol\lambda) where \boldsymbol \lambda is a vector Lagrange multiplier representing contact stresses (normal and tangential).

The semi-discrete version of the field equations in the case, say, of a linear problem is of the form

\mathbf M \ddot{\mathbf u} + \mathbf K \mathbf u + {\color{}\mathbf C}\boldsymbol\lambda = \mathbf f.

Is there a way I can compute the rectangular matrix \textcolor{}{\mathbf C} in GetFEM?

In my case, I use a \mathbb P_1 discretization for \mathbf u and a piecewise \mathbb P_0 discretization for \boldsymbol\lambda on the contact boundary.

Thanks for your help!

well, matrices and matrices…, I would like to ask you back what do we need the matrix C for?

GetFEM computes the entire tangent matrix for the global system including all variables that you have defined in a GetFEM model. This includes the matrix you are asking about. The important question is what do you need to do with it?

How it works:

md = gf.Model("real")
md.add_fem_variable("u", mfu)
md.add_fem_data("u_prev", mfu)
md.add_filtered_fem_variable("lambda", RG_SLAVE)
md.add_initialized_data("r", 10*E) # where E is the Young's module
md.add_initialized_data("mu", 0.1) # frict. coeff.
md.add_initialized_data("alpha", 1)
md.add_integral_contact_between_nonmatching_meshes_brick\
  (mim, "u", "u", "lambda", "r", "mu", RG_SLAVE, RG_MASTER, 1,
   "alpha", "u_prev", "u_prev")
  # not sure whether it is slave first and master second or the opposite
  # 1 is the choice of which small sliding contact algorithm to use

md.assembly("build_all")
K=md.tangent_matrix()

IIu = range(*md.interval_of_variable("u").cumsum())
IIλ = range(*md.interval_of_variable("lambda").cumsum())
C = gf.Spmat("copy", K, IIu, IIλ)

Hi Kostas,

Thanks for your prompt feedback. I would like to implement the active set method for persistent contact (to ensure energy conservation) problems in 3D. For this, I have to perform static condensation on the Lagrange multiplier, hence some manipulations on the coupling matrix are necessary, unless I am missing something.

It is rarely necessary to perform manipulations on the tangent matrices. This misconceptions comes from the fact that, historically, people have been working with linear problems only, and then they try to fit nonlinear problems into the limited linear framework that they got used to. It ends up going in circles to reach their goal instead of going straight.

The straight forward way is to setup your nonlinear residual equations in the correct way in the first place. The only thing that maters are the residual equations. That’s what you solve for.

Ok, I see.

So, let me put it this way: Is it possible to state nonlinear equations in terms of dofs-based unknowns (node-based, edge-based, or face-based primal or dual unknowns, i.e. displacements or stresses), instead of using GWFL, for which the unknowns are functions?

I see what your use case is. What you refer to is discretize-first and model afterwards, instead of model first and discretize afterwards.

There are not many cases where writing the (nonlinear) equation after the discretization makes much sense because the physical phenomena that matter for continuum mechanics, are non-discrete. That’s why it is called continuum mechanics.

Having said that, there are examples of discretize-first and model afterwards in GetFEM, but not as many as the opposite. Relevant methods:

model.add_nodal_contact_between_nonmatching_meshes_brick
model.add_nodal_contact_with_rigid_obstacle_brick
model.add_explicit_matrix
model.add_explicit_rhs

There are many techniques to mention, for example the GWFL based model.interpolation(...) function can be used as a tool for a discretize-first model-afterwarts purpose. Another trick to exploit GWFL is to use element-wise constant fem (with mf.set_classical_fem(0)) and then use an integration method with just 1 integration point (in the center of the element).

If you want to see a classical example of how a nonlinear equation is set-up on an already discretized system, I would recommend you to study in depth the file getfem_contact_and_friction_nodal.cc