A nonlinear finite element connector for the simulation of bolted assemblies

Abstract

Fine scale computations of bolted assemblies are generally too costly and hardly tractable within an optimization process. Thus, finite elements (FE) connectors or user-elements are commonly used in FE commercial codes by engineers as substitutes for bolts. In this paper, a non-linear FE connector with its identification methodology is proposed to model the behaviour of a single-bolt joint. The connector model is based on design parameters (bolt prestress, friction coefficient, bolt characteristics...). The axial behaviour of the connector reflects the preload effect and the axial bolt stiffness. The tangential connector behaviour accounts for frictional phenomena that occur in the bolt’s vicinity due to preload thanks to an elasto-plastic analogy for friction. Tangential and normal behaviours identification is performed on a generic elementary single bolt joint. The connector has been implemented in ABAQUS through a user-element subroutine. Comparisons of the quasistatic responses between full fine scale 3D computations and 3D simulations with connectors on various bolted assemblies are proposed. Results are in good agreement and a significant gain in terms of CPU time is obtained.

Appendices

Appendix A: Radial return mapping algorithm

• 1. Increment of tangential gap within the time step $\varDelta t_{n+1}$

  As specified in (5), the expression of the tangential displacement jump $\varvec{g}_{T}$ is given by

  $$\begin{aligned} \varDelta g_{T,(n+1)} = (\varvec{I}_3- \varvec{n}_{(0)}\otimes \varvec{n}_{(0)})\left( \varDelta u_{2,(n+1)}-\varDelta u_{1,(n+1)}\right) \end{aligned}$$

  (23)

  in the case of small pertubations.

• 2. Computation of the elastic trial state from (7) and evaluation of the slip criterion (8) at time $t_{n+1}$. The trial state consists in assuming that the considered increment is purely elastic, and then correcting it according to the value of the obtained threshold function. Since the total slip $\varvec{g}_{T,(n+1)} = \varvec{g}_{T,(n)} + \varDelta \varvec{g}_{T,(n+1)}$ is decomposed into an elastic and a plastic part, one has:

  $$\begin{aligned} \begin{array}{rcl} \varvec{t}_{T,(n+1)}^{tr}&{} = &{} c_{T}\left( \varvec{g}_{T,(n+1)}-\varvec{g}^s_{T,(n)}\right) \\ &{}=&{} \varvec{t}_{T,(n)} + c_T\,\varDelta \varvec{g}_{T,(n+1)}\\ \end{array} \end{aligned}$$

  (24)

  Here the vector $\varvec{t}_{T,(n)} = c_T\left( \varvec{g}_{T,(n)}-\varvec{g}^s_{T,(n)}\right)$ is the tangential force at time $t_n$. A value of the slip criterion which fulfills $f_{s,(n+1)}^{tr}\leqslant 0$ indicates stick. For $f_{s,(n+1)}^{tr}>0$ sliding occurs in the tangential direction, and a return mapping of the trial tractions to the slip surface has to be performed. For Coulomb’s model, one simply has:

  $$\begin{aligned} f_{s,(n+1)}^{tr} = \Vert \varvec{t}^{tr}_{T,(n+1)} \Vert - \mu \,p_{N,(n+1)} \end{aligned}$$

  (25)

  where $p_{N,(n+1)} = c_N \vert g_{N,(n+1)}\vert$ (since $g_{N,(n+1)}<0$ when contact occurs) and where $c_N$ is the normal connector stiffness.

• 3. Return mapping procedure is derived from the time integration algorithm. In case the implicit Euler scheme is applied to approximate the evolution equation (9), one obtains:

  $$\begin{aligned} \varvec{g}^s_{T,(n+1)} = \varvec{g}^s_{T,(n)} + \lambda \,\varvec{n}_{T,(n+1)} \quad \end{aligned}$$

  (26)

  with

  $$\begin{aligned} \varvec{n}_{T,(n+1)} =\dfrac{\varvec{t}_{T,(n+1)} }{\Vert \varvec{t}_{T,(n+1)} \Vert } \end{aligned}$$

  (27)

  With the standard arguments regarding the projection schemes [36], one obtains:

  $$\begin{aligned} \left\{ \begin{array}{rcl} \varvec{t}_{T,(n+1)} &{} = &{} \varvec{t}^{tr}_{T,(n+1)}-\lambda \,c_T\,\varvec{n}_{T,(n+1)} \\ \varvec{n}_{T,(n+1)} &{} = &{} \varvec{n}_{T,(n+1)}^{tr} \quad \end{array}\right. \end{aligned}$$

  (28)

  with

  $$\begin{aligned} \varvec{n}_{T,(n+1)}^{tr} =\dfrac{\varvec{t}_{T,(n+1)}^{tr} }{\Vert \varvec{t}_{T,(n+1)}^{tr} \Vert } \end{aligned}$$

  (29)

  The multiplication of (28) by $\varvec{n}_{T,(n+1)}$ yields the condition from which $\lambda$ may be computed:

  $$\begin{aligned} \kappa (\lambda ) = \Vert \varvec{t}^{tr}_{T,(n+1)} \Vert - \mu p_{N,(n+1)} - c_T\,\lambda = 0 \end{aligned}$$

  (30)

  since $\Vert \varvec{t}_{T,(n+1)} \Vert - \mu p_{N,(n+1)} = 0$ on the yield surface. Thus, one gets:

  $$\begin{aligned} \lambda = \dfrac{1}{c_T}\left( \Vert \varvec{t}^{tr}_{T,(n+1)} \Vert - \mu \,p_{N,(n+1)} \right) \end{aligned}$$

  (31)

  Once $\lambda$ is known, one obtains the explicit results for Coulomb’s model:

  $$\begin{aligned} \begin{array}{rcl} \varvec{t}_{T,(n+1)} &{} = &{} \mu \,p_{N,(n+1)}\, \varvec{n}_{T,(n+1)}^{tr}\\ \varvec{g}^s_{T,(n+1)} &{} =&{} \varvec{g}^s_{T,(n)} + \lambda \ \varvec{n}_{T,(n+1)}^{tr} \end{array} \end{aligned}$$

  (32)

  This update completes the local integration algorithm for the frictional interface law.

Appendix B: Elasto-plastic tangent modulus

Let us consider the yield surface described in (25):

$$\begin{aligned} f_{s}(\varvec{t}_{T}, p_N) = \Vert \varvec{t}_{T} \Vert - \mu p_{N} \end{aligned}$$

(33)

Thus:

$$\begin{aligned} \begin{array}{rcl} \dot{f}_{s}(\varvec{t}_{T}, p_N) &{}=&{} \dfrac{\partial f_s}{\partial \varvec{t}_T}\dot{\varvec{t}}_T + \dfrac{\partial f_s}{\partial p_N}\dot{p}_N\\ &{}=&{} \dfrac{\partial f_s}{\partial \varvec{t}_T} \, c_T\left( \dot{\varvec{g}}_{T}-\dot{\varvec{g}}^s_{T}\right) + \dfrac{\partial f_s}{\partial p_N}\dot{p}_N\\ &{}=&{} \dfrac{\partial f_s}{\partial \varvec{t}_T} \, c_T\left( \dot{\varvec{g}}_{T}-\lambda \dfrac{\partial f_s(\varvec{t}_T)}{\partial \varvec{t}_T}\right) + \dfrac{\partial f_s}{\partial p_N}\dot{p}_N \end{array} \end{aligned}$$

(34)

The consistency condition ($\lambda \dot{f}_s(\varvec{t}_{T}, p_N) = 0$ if $f_s(\varvec{t}_T) = 0$) enables one to determine $\lambda$:

$$\begin{aligned} \lambda = \dfrac{c_T\dfrac{\partial f_s}{\partial \varvec{t}_T}\dot{\varvec{g}}_{T}+\dfrac{\partial f_s}{\partial p_N}\dot{p}_N }{c_T\dfrac{\partial f_s}{\partial \varvec{t}_T}\dfrac{\partial f_s}{\partial \varvec{t}_T}} \end{aligned}$$

(35)

Since $\dot{\varvec{g}}^s_T=\dot{\gamma }\,\dfrac{\partial f_s(\varvec{t}_T)}{\partial \varvec{t}_T} = \dot{\gamma } \,\varvec{n}_T=\lambda \,\varvec{n}_T$, one deduces that:

$$\begin{aligned} \lambda = \dfrac{c_T\,\varvec{n}_T^T\,\dot{\varvec{g}}_{T}+\dfrac{\partial f_s}{\partial p_N}\dot{p}_N }{c_T\,\varvec{n}_T^T\varvec{n}_T} = \dfrac{c_T\,\varvec{n}_T^T\,\dot{\varvec{g}}_{T}+\dfrac{\partial f_s}{\partial p_N}\dot{p}_N }{c_T} \end{aligned}$$

(36)

For the Coulomb’s friction model, one obtains:

$$\begin{aligned} \lambda = \dfrac{c_T\,\varvec{n}_T^T\,\dot{\varvec{g}}_{T}-\mu \dot{p}_N }{c_T} \end{aligned}$$

(37)

This provides the expression for $\dot{\varvec{t}}_{T}$ as a function of $\dot{\varvec{g}}_T$ and $\dot{g}_N$:

$$\begin{aligned} \begin{array}{rcl} \dot{\varvec{t}}_{T} &{}=&{} c_T\left( \dot{\varvec{g}}_{T}-\dot{\varvec{g}}^s_{T}\right) \\ &{} = &{} c_T\, \dot{\varvec{g}}_{T}-c_T\,\lambda \,\varvec{n}_T \\ &{}=&{} c_T\, \dot{\varvec{g}}_{T}-\dfrac{c_T\,\varvec{n}_T\left( c_T\,\varvec{n}_T^T\,\dot{\varvec{g}}_{T}\right) -\mu \,c_T \,\dot{p}_N\,\varvec{n}_T}{c_T}\\ &{}=&{} c_T\left( \mathbf {I_3} - \varvec{n}_T\otimes \varvec{n}_T\right) \dot{\varvec{g}}_{T} + \mu \, c_N\,\mathrm {sign}\left( g_{N}\right) \dot{g}_N\,\varvec{n}_T \end{array} \end{aligned}$$

(38)

where $\varvec{n}_T =\dfrac{\varvec{t}_T}{\Vert \varvec{t}_T\Vert }$. Note that $\text {sign}(g_N)=-1$.

Note that the frictional term (second term in right-hand side) makes of the matrix which links $\dot{\varvec{t}}_{T}$ to the normal gap increment $\dot{g}_N$, non-symmetric. This is because the Coulomb’s law of friction can be viewed as a non-associative constitutive model.

Tangential slip increment and gap increment may be written in matrix form. One can express the variation $\delta g_{N}\,\varvec{n}_T$ and $\delta \varvec{g}_{T}$ in terms of connector nodes displacement vector variation $\delta \varvec{g}$ as expressed in (5) as:

$$\begin{aligned} \begin{array}{rcl} \delta g_{N}\varvec{n}_T&{}=&{}\left( \left( \delta \varvec{x}_2-\delta \varvec{x}_1\right) \cdot \varvec{n}\right) \varvec{n}_T\\ &{}=&{} \begin{pmatrix} -\varvec{n}_T\otimes \varvec{n}\ \ &{} \varvec{n}_T\otimes \varvec{n}\end{pmatrix} \begin{pmatrix} \delta \varvec{x}_1\\ \delta \varvec{x}_2\end{pmatrix}\\ &{}=&{} \begin{pmatrix} -\varvec{n}_T\otimes \varvec{n}\ \ &{} \varvec{0}_{3\times 3}\ \ &{} \varvec{n}_T\otimes \varvec{n}\ \ &{} \varvec{0}_{3\times 3} \end{pmatrix} \delta \varvec{g}\\ &{}=&{}\varvec{N}_T \,\delta \varvec{g}\end{array} \end{aligned}$$

(39)

where $\varvec{0}_{n\times m}$ represents a n by m null matrix. With the definition of tangential gap given in (5), one gets:

$$\begin{aligned} \begin{array}{rcl} \delta \varvec{g}_{T}&{}=&{} (\mathbf {I_3} - \varvec{n}\otimes \varvec{n})\left( \delta \varvec{x}_2-\delta \varvec{x}_1\right) \\ &{}=&{} \begin{pmatrix} -\left( \mathbf {I_3} - \varvec{n}\otimes \varvec{n}\right) \ \ &{} \mathbf {I_3} - \varvec{n}\otimes \varvec{n}\end{pmatrix} \begin{pmatrix} \delta \varvec{x}_1\\ \delta \varvec{x}_2\end{pmatrix}\\ &{}=&{} \begin{pmatrix} -\left( \mathbf {I_3} - \varvec{n}\otimes \varvec{n}\right) \ \ &{} \varvec{0}_{3\times 3}\ \ &{} \mathbf {I_3} - \varvec{n}\otimes \varvec{n}\ \ &{} \varvec{0}_{3\times 3} \end{pmatrix} \delta \varvec{g}\\ &{}=&{}\varvec{M} \,\delta \varvec{g}\end{array} \end{aligned}$$

(40)

One deduces from (38) that:

$$\begin{aligned} \begin{array}{rcl} \varDelta \varvec{t}_T &{}=&{} c_T\left( \mathbf {I_3} -\varvec{n}_T\otimes \varvec{n}_T\right) \varDelta \varvec{g}_{T} +\mu \,c_N\,\mathrm {sign}\left( g_{N}\right) \varDelta {g}_N\,\varvec{n}_T\\ &{}=&{} \left[ c_T\left( \mathbf {I_3} -\varvec{n}_T\otimes \varvec{n}_T\right) \varvec{M}+\mu \,c_N\,\mathrm {sign}\left( g_{N}\right) \varvec{N_T} \right] \varDelta \varvec{d}\end{array} \end{aligned}$$

(41)

and, finally:

$$\begin{aligned}&\varDelta \varvec{t}_T\cdot \delta \varvec{g}_{T}= \delta \varvec{g}_{T}^T\varDelta \varvec{t}_T = \nonumber \\&\delta \varvec{d}^T \left( \varvec{M}^{T} \left[ c_T\left( \mathbf {I_3} -\varvec{n}_T\otimes \varvec{n}_T\right) \varvec{M}+ \mu \,c_N\,\mathrm {sign}\left( g_{N}\right) \varvec{N}_T \right] \right) \varDelta \varvec{d}\nonumber \\&=\delta \varvec{d}^T \varvec{K}_T \,\varDelta \varvec{d}\end{aligned}$$

(42)

Note that matrix $\varvec{K}_T$ is nonsymmetric.

Appendix C: Summary of the connector integration algorithm

As a reminder, a quantity a denoted by $a_{(n)}$ corresponds to the value of a evaluated at time increment n. However, to alleviate the notations, the normal $\mathbf {n}_{(n+1)}$ is noted only by $\mathbf {n}$ and stays equal to $\mathbf {n}_{(0)}$ under the small perturbations assumption.

The different steps of the connector model integration are summarized in Algorithm 1. $\mathbf {P}_T$ represents the matrix for passing from the global frame to the local frame of the connector.

For a two-dimensional problem, the stiffness matrix of a Timoshenko beam $\mathbf {K}_{bolt}$ is expressed as a function of the normal stiffness $c_N$, the bending stiffness $c_{bolt}$, its length L and the shear coefficient $\varPhi$ depending on the section geometry as:

$$\begin{aligned} \mathbf {K}_{bolt} {=} \begin{pmatrix} c_N &{} 0 &{} 0 &{} -c_N &{} 0 &{} 0 \\ 0 &{} c_{bolt} &{} c_{bolt}\frac{L}{2} &{} 0 &{} -c_{bolt} &{} c_{bolt}\frac{L}{2} \\ 0 &{} c_{bolt}\frac{L}{2} &{} c_{bolt}\frac{(4+\varPhi )L^2}{12} &{} 0 &{} -c_{bolt}\frac{L}{2} &{} c_{bolt}\frac{(2-\varPhi )L^2}{12} \\ -c_N &{} 0 &{} 0 &{} c_N &{} 0 &{} 0 \\ 0 &{} -c_{bolt} &{} -c_{bolt}\frac{L}{2} &{} 0 &{} c_{bolt} &{} -c_{bolt}\frac{L}{2} \\ 0 &{} c_{bolt}\frac{L}{2} &{} c_{bolt}\frac{(2-\varPhi )L^2}{12} &{} 0 &{} -c_{bolt}\frac{L}{2} &{} c_{bolt}\frac{(4+\varPhi )L^2}{12} \end{pmatrix} \end{aligned}$$