Tangent contact stiffness

Figure 134: Visualization of the tangential differential displacements

To find the tangent contact stiffness matrix, please look at Figure 134, part a). At the beginning of a concrete time increment, characterized by time $t_n$, the slave node at position $\boldsymbol{p}_n$ corresponds to the projection vector $\boldsymbol{q}_n$ on the master side. At the end of the time increment, characterized by time $t_{n+1}$ both have moved to positions $\boldsymbol{p}_{n+1}$ and $\boldsymbol{q}_{n+1}$, respectively. The differential displacement between slave and master surface changed during this increment by the vector $\boldsymbol{s}$ satisfying:

$$\boldsymbol{s}=(\boldsymbol{p}_{n+1} - \boldsymbol{q}_{n+1}) -(\boldsymbol{p}_n - \boldsymbol{q}_n).$$ (185)

Here, $\boldsymbol{q_{n+1}}$ satisfies

$$\boldsymbol{q}_{n+1}= \sum_j \varphi_j(\xi^m _n, \eta^m_n) \boldsymbol{q_j}_{n+1} .$$ (186)

Since (the dependency of $\varphi_j$ on $\xi$ and $\eta$ is dropped to simplify the notation)

$$\boldsymbol{p}_{n}= \boldsymbol{X} + \boldsymbol{u}_{n},$$ (187)

$$\boldsymbol{p}_{n+1}= \boldsymbol{X} + \boldsymbol{u}_{n+1},$$ (188)

$$\boldsymbol{q}_{n}= \sum_j \varphi_j \left[ \boldsymbol{X_j} + (\boldsymbol{u_j})_{n} \right ],$$ (189)

$$\boldsymbol{q}_{n+1}= \sum_j \varphi_j \left[ \boldsymbol{X_j} + (\boldsymbol{u_j})_{n+1} \right ],$$ (190)

this also amounts to

$$\boldsymbol{s}= \boldsymbol{u}_{n+1} - \sum_j \varphi_j (\boldsym... ...- \left [ \boldsymbol{u}_{n} - \sum_j \varphi_j (\boldsymbol{u_j})_{n} \right].$$ (191)

Notice that the local coordinates take the values of time $t_n$ (the superscript m denotes iteration m within the increment). The differential tangential displacement now amounts to:

$$\boldsymbol{t} \equiv \boldsymbol{t_{n+1}}=\boldsymbol{t_{n}} + \boldsymbol{s}- s \boldsymbol{n},$$ (192)

where

$$s=\boldsymbol{s} \cdot \boldsymbol{n}.$$ (193)

Derivation w.r.t. $\boldsymbol{u_i}$ satisfies (straightforward differentiation):

$$\frac{\partial \boldsymbol{t} }{\partial \boldsymbol{u_i} } = \fr... ...ldsymbol{u_i} } - s \frac{\partial \boldsymbol{n} }{\partial \boldsymbol{u_i} }$$ (194)

$$\frac{\partial s}{\partial \boldsymbol{u_i} } = \frac{\boldsymbol... ...bol{m}\Vert } \cdot \frac{\partial \boldsymbol{s} }{\partial \boldsymbol{u_i} }$$ (195)

and

$$\frac{\partial \boldsymbol{n} }{\partial \boldsymbol{u_i} } = \fr... ...bol{m } \otimes \frac{\partial \Vert \boldsymbol{m} \Vert }{\boldsymbol{u_i} }.$$ (196)

The derivative of $\boldsymbol{s}$ w.r.t. $\boldsymbol{u_i}$ can be obtained from the derivative of $\boldsymbol{r}$ w.r.t. $\boldsymbol{u_i}$ by keeping $\xi$ and $\eta$ fixed (notice that the derivative is taken at $t_{n+1}$, consequently, all derivatives of values at time $t_n$ disappear):

$$\frac{\partial \boldsymbol{s} }{\partial \boldsymbol{u_i} } = \delta _{ip} \boldsymbol{I} - \varphi_i (1 - \delta _{ip}) \boldsymbol{I}.$$ (197)

Physically, the tangential contact equations are as follows (written at the location of slave node p):

• an additive decomposition of the differential tangential displacement in stick $\boldsymbol{t^e}$ and slip $\boldsymbol{t^p}$ :

  $$\boldsymbol{t}= \boldsymbol{t^e} + \boldsymbol{t^p}.$$ (198)

  
  

• a stick law ( $K_t\equiv \lambda a_p$, where $\lambda$ is the stick slope and $a_p$ the representative slave area for the slave node at stake) defining the tangential force exerted by the slave side on the master side at the location of slave node p:

  $$\boldsymbol{F_T} = K_t \boldsymbol{t^e}.$$ (199)

  
  

• a Coulomb slip limit:

  $$\Vert \boldsymbol{F_T} \Vert \le \mu \Vert \boldsymbol{F_N} \Vert$$ (200)

  
  

• a slip evolution equation:

  $$\boldsymbol{\boldsymbol{\dot{t}^p} }= \dot{\gamma} \frac{\boldsymbol{F_T} }{\Vert \boldsymbol{F_T} \Vert }$$ (201)

  
  

First, a difference form of the additive decomposition of the differential tangential displacement is derived. Starting from

$$\boldsymbol{t}= \boldsymbol{t^e} + \boldsymbol{t^p},$$ (202)

one obtains after taking the time derivative:

$$\boldsymbol{\dot{t}}= \boldsymbol{\dot{t}^e} + \boldsymbol{\dot{t}^p}.$$ (203)

Substituting the slip evolution equation leads to:

$$\dot{\gamma} \frac{\boldsymbol{F_T} }{\Vert \boldsymbol{F_T} \Vert }=\boldsymbol{\dot{t}} - \boldsymbol{\dot{t}^e},$$ (204)

and after multiplying with $K_t$:

$$K_t \dot{\gamma} \frac{\boldsymbol{F_T} }{\Vert \boldsymbol{F_T} \Vert }=K_t \boldsymbol{\dot{t}} - K_t \boldsymbol{\dot{t}^e}$$ (205)

Writing this equation at $t_{n+1}$, using finite differences (backwards Euler), one gets:

$$K_t \Delta \gamma_{n+1} \frac{\boldsymbol{F_T}_{n+1} }{\Vert\bold... ...lta \boldsymbol{t}_{n+1} - K_t \boldsymbol{t^e}_{n+1} + K_t \boldsymbol{t^e}_n,$$ (206)

where $\Delta \gamma_{n+1} \equiv \dot{\gamma} _{n+1} \Delta t$ and $\Delta \boldsymbol{t}_{n+1} \equiv \boldsymbol{t}_{n+1} - \boldsymbol{t}_n$. The parameter $K_t$ is assumed to be independent of time.

Now, the radial return algorithm will be described to solve the governing equations. Assume that the solution at time $t_n$ is known, i.e. $\boldsymbol{t^e}_n$ and $\boldsymbol{t^p}_n$ are known. Using the stick law the tangential forc $\boldsymbol{F_T}_n$ can be calculated. Now we would like to know these variables at time $t_{n+1}$, given the total differential tangential displacement $\boldsymbol{t_{n+1}}$. At first we construct a trial tangential force $\boldsymbol{F_T}_{n+1}^{trial}$ which is the force which arises at time $t_{n+1}$ assuming that no slip takes place from $t_n$ till $t_{n+1}$. This assumption is equivalent to $\boldsymbol{t^p}_{n+1}=\boldsymbol{t^p}_n$. Therefore, the trial tangential force satisfies (cf. the stick law):

$$\boldsymbol{F_T}_{n+1}^{trial} = K_t(\boldsymbol{t}_{n+1} - \boldsymbol{t^p}_n).$$ (207)

Now, this can also be written as:

$$\boldsymbol{F_T}_{n+1}^{trial} = K_t(\boldsymbol{t}_{n+1} - \boldsymbol{t}_n + \boldsymbol{t}_n - \boldsymbol{t^p}_n).$$ (208)

or

$$\boldsymbol{F_T}_{n+1}^{trial} = K_t \Delta \boldsymbol{t}_{n+1} + K_t \boldsymbol{t^e}_n.$$ (209)

Using Equation (206) this is equivalent to:

$$\boldsymbol{F_T}_{n+1}^{trial} = K_t \Delta \gamma_{n+1} \frac{\b... ...l{F_T}_{n+1} }{\Vert\boldsymbol{F_T}_{n+1} \Vert} + K_t \boldsymbol{t^e}_{n+1},$$ (210)

or

$$\boldsymbol{F_T}_{n+1}^{trial} = (K_t \Delta \gamma_{n+1} \frac{1 }{\Vert\boldsymbol{F_T}_{n+1} \Vert} + 1) \boldsymbol{F_T}_{n+1} .$$ (211)

From the last equation one obtains

$$\boldsymbol{F_T}_{n+1}^{trial} \parallel \boldsymbol{F_T}_{n+1}$$ (212)

and, since the terms in brackets in Equation (211) are both positive:

$$\Vert \boldsymbol{F_T}_{n+1}^{trial}\Vert = K_t \Delta \gamma_{n+1} + \Vert \boldsymbol{F_T}_{n+1} \Vert.$$ (213)

The only equation which is left to be satisfied is the Coulomb slip limit. Two possibilities arise:

1. $\Vert \boldsymbol{F_T}_{n+1}^{trial}\Vert \le \mu \Vert \boldsymbol{F_N}_{n+1} \Vert$.

   In that case the Coulomb slip limit is satisfied and we have found the solution:

   $$\boldsymbol{F_T}_{n+1} = \boldsymbol{F_T}_{n+1}^{trial} = K_t( \boldsymbol{t_{n+1}} - \boldsymbol{t^p}_n)$$ (214)

   
   

   and

   $$\frac{\partial \boldsymbol{F_T}_{n+1}}{\partial \boldsymbol{t}_{n+1} } = K_t \boldsymbol{I}.$$ (215)

   
   

   No extra slip occurred from $t_n$ to $t_{n+1}$.

2. $\Vert \boldsymbol{F_T}_{n+1}^{trial}\Vert > \mu \Vert \boldsymbol{F_N}_{n+1} \Vert$.

   In that case we project the solution back onto the slip surface and require $\Vert \boldsymbol{F_T}_{n+1}\Vert = \mu \Vert \boldsymbol{F_N}_{n+1} \Vert$. Using Equation (213) this leads to the following expression for the increase of the consistency parameter $\gamma$:

   $$\Delta \gamma _{n+1}= \frac{\Vert\boldsymbol{F_T}_{n+1}^{trial}\Vert- \mu \Vert \boldsymbol{F_N}_{n+1} \Vert }{K_t},$$ (216)

   
   

   which can be used to update $\boldsymbol{t^p}$ (by using the slip evolution equation):

   $$\Delta \boldsymbol{t^p}= \Delta \gamma_{n+1} \frac{\boldsymbol{F_... ...rac{\boldsymbol{F_T}_{n+1}^{trial} }{\Vert\boldsymbol{F_T}_{n+1}^{trial} \Vert}$$ (217)

   
   

   The tangential force can be written as:

   $$\boldsymbol{F_T}_{n+1}= \Vert \boldsymbol{F_T}_{n+1} \Vert \frac{... ...c{\boldsymbol{F_T}_{n+1}^{trial}}{\Vert \boldsymbol{F_T}_{n+1}^{trial} \Vert }.$$ (218)

   
   

   Now since

   $$\frac{\partial \Vert \boldsymbol{a} \Vert}{\partial \boldsymbol{b... ...ymbol{a}\Vert } \cdot \frac{\partial \boldsymbol{a} }{\partial \boldsymbol{b} }$$ (219)

   
   

   and

   $$\frac{\partial }{\partial \boldsymbol{b} } \left (\frac{\boldsymb... ...ght ) \right ] \cdot \frac{\partial \boldsymbol{a} }{\partial \boldsymbol{b} },$$ (220)

   
   

   where $\boldsymbol{a}$ and $\boldsymbol{b}$ are vectors, one obtains for the derivative of the tangential force:

   $$\frac{\partial \boldsymbol{F_T}_{n+1} }{\partial \boldsymbol{t}_{n+1} } = \mu \boldsymbol{\xi }_{n+1} \otimes \left [ \frac{\boldsymbol{F... ...t \frac{\partial \boldsymbol{F_N}_{n+1}}{\partial \boldsymbol{t}_{n+1}} \right]$$

   $$+ \mu \frac{\Vert\boldsymbol{F_N}_{n+1} \Vert}{\Vert\boldsymbol{F... ...frac{\partial \boldsymbol{F_T}_{n+1}^{trial} }{\partial \boldsymbol{t}_{n+1}, }$$ (221)

   
   

   where

   $$\boldsymbol{\xi }_{n+1} \equiv \frac{\boldsymbol{F_T}_{n+1}^{trial}}{\Vert \boldsymbol{F_T}_{n+1}^{trial} \Vert }.$$ (222)

   
   

   One finally arrives at (using Equation (215)

   $$\frac{\partial \boldsymbol{F_T}_{n+1} }{\partial \boldsymbol{u_i}_{n+1} } = \mu \boldsymbol{\xi }_{n+1} \otimes \left [ -\boldsymbol{n} \cdot \frac{\partial \boldsymbol{F_N}_{n+1}}{\partial \boldsymbol{u_i}_{n+1}} \right]$$

   $$+ \mu \frac{\Vert\boldsymbol{F_N}_{n+1} \Vert}{\Vert\boldsymbol{F... ...ot K_t \frac{\partial \boldsymbol{t}_{n+1} }{\partial \boldsymbol{u_i}_{n+1}. }$$ (223)

   
   

   All quantities on the right hand side are known now (cf. Equation (170) and Equation (194)).

   In CalculiX, for node-to-face contact, Equation (185) is reformulated and simplified. It is reformulated in the sense that $\boldsymbol{q}_{n+1}$ is assumed to be the projection of $\boldsymbol{p}_{n+1}$ and $\boldsymbol{q}_n$ is written as (cf. Figure 134, part b))

   $$\boldsymbol{q}_{n}= \sum_j \varphi_j(\xi^m _{n+1}, \eta^m_{n+1}) \boldsymbol{q_j}_{n} .$$ (224)

   
   

   Part a) and part b) of the figure are really equivalent, they just represent the same facts from a different point of view. In part a) the projection on the master surface is performed at time $t_n$, and the differential displacement is calculated at time $t_{n+1}$, in part b) the projection is done at time $t_{n+1}$ and the differential displacement is calculated at time $t_n$. Now, the actual position can be written as the sum of the undeformed position and the deformation, i.e. $\boldsymbol{p}=(\boldsymbol{X}+\boldsymbol{v})^s$ and $\boldsymbol{q}=(\boldsymbol{X}+\boldsymbol{v})^m$ leading to:

   $$\boldsymbol{s}=(\boldsymbol{X}+\boldsymbol{v})^s_{n+1} - (\boldsy... ...mbol{v})^s_n + (\boldsymbol{X}+\boldsymbol{v})^m_n(\xi^m _{n+1}, \eta^m_{n+1}).$$ (225)

   
   

   Since the undeformed position is no function of time it drops out:

   $$\boldsymbol{s}=\boldsymbol{v}^s_{n+1} - \boldsymbol{v}^m_{n+1}(\x... ...a^m_{n+1}) -\boldsymbol{v}^s_n + \boldsymbol{v}^m_n(\xi^m _{n+1}, \eta^m_{n+1})$$ (226)

   
   

   or:

   $$\boldsymbol{s} =\boldsymbol{v}^s_{n+1} - \boldsymbol{v}^m_{n+1}(\xi^m _{n+1}, \e... ...n+1}) -\boldsymbol{v}^s_n + \boldsymbol{v}^m_n(\xi^m _{n}, \eta^m_{n})\noindent$$ (227)

   $$+ \boldsymbol{v}^m_n(\xi^m _{n+1}, \eta^m_{n+1})- \boldsymbol{v}^m_n(\xi^m _{n}, \eta^m_{n})$$ (228)

   
   

   Now, the last two terms are dropped, i.e. it is assumed that the differential deformation at time $t_n$ between positions $(\xi^m _{n}, \eta^m_{n})$ and $(\xi^m _{n+1}, \eta^m_{n+1})$ is neglegible compared to the differential motion from $t_n$ to $t_{n+1}$. Then the expression for $\boldsymbol{s}$ simplifies to:

   $$\boldsymbol{s}=\boldsymbol{v}^s_{n+1} - \boldsymbol{v}^m_{n+1}(\x... ...\eta^m_{n+1}) -\boldsymbol{v}^s_n + \boldsymbol{v}^m_n(\xi^m _{n}, \eta^m_{n}),$$ (229)

   
   

   and the only quantity to be stored is the difference in deformation between $\boldsymbol{p}$ and $\boldsymbol{q}$ at the actual time and at the time of the beginning of the increment.