A Stabilized Hybrid Discontinuous Galerkin Method for Nearly Incompressible Linear Elasticity Problem

ABSTRACT

In this work, a primal hybrid finite element method for nearly incompressible linear elasticity problem on triangular meshes is shown. This method consists of coupling local discontinuous Galerkin problems to the primal variable with a global problem for the Lagrange multiplier, which is identified as the trace of the displacement field. Also, a local post-processing technique is used to recover stress approximations with improved rates of convergence in H(div) norm. Numerical studies show that the method is locking free even using equal or different orders for displacement and stress fields and optimal convergence rates are obtained.

Keywords:
linear elasticity; discontinuous Galerkin method; stabilization hybrid method; locking free

RESUMO

Neste trabalho, um método de elementos finitos híbrido primal para o problema de elasticidade linear quase incompressível em malhas triangulares é mostrado. Este método consiste no acomplamento de problemas de Galerkin descontínuos locais para a variável primal com um problema global para o multiplicador de Lagrange, o qual é identificado como o traço do campo de deslocamento. Além disso, uma técnica de pós-processamento local é usada para recuperar as aproximações das tensões com melhores taxas de convergeˆncia na norma H(div). Estudos numéricos mostram que o método é “locking free” quando é usado ordens iguais ou diferentes para os campos de deslocamento e tensão e taxas ótimas de convergência são obtidas.

Palavras-chave:
elasticidade linear; método de Galerkin descontínuo; método híbrido de estabilização; locking free

1 INTRODUCTION

Historically, many finite elements have been constructed, in the last decades, for the incompressible and nearly incompressible linear elasticity materials, in order to avoid the “volumetric locking” phenomena. This is an undesirable unphysical stiff response to deformation exhibited when using the standard Galerkin method based upon the conventional displacement formulation with low order elements and the Poisson’s ratio tends to 1/2. Many successful strategies to avoid the locking effects based on mixed finite element methods^{1), (2), (3), (4}, higher-order method⁵, and discontinuous Galerkin finite element methods^{6), (7), (8), (9), (10} can be found in the literature.

Discontinuous Galerkin (DG) methods are naturally a suitable alternative for solving linear elasticity problems. Robustness, local conservation and flexibility for implementing h and p-adaptivity strategies are well known advantages of DG methods stemming from the use of finite element spaces consisting of discontinuous piecewise polynomials. A natural connection between DG formulations and hybrid methods have been exploited successfully in many problems^{11), (12), (13), (14), (15} and they are still being developed. These hybrid formulations have improved stability, robustness and flexibility of the DG methods with reduced complexity and computational cost. Hybrid finite element methods are characterized by the introduction of new unknowns variables, the Lagrange multipliers, defined at the edges of the elements to weakly impose continuity on the element interfaces. Their reduced complexity and computational cost is obtained by eliminating all degrees of freedom of the primal variables at the element level resulting a global system only with degrees of freedom of the multipliers.

In this paper, we extend the studies presented in^{16), (17), (18} considering a primal Stabilized Hybrid Discontinuous Galerkin (SHDG) method, for the nearly incompressible linear elasticity problem using triangular elements. The multiplier, identified with the trace of the displacement field, is interpolated discontinuously. The method is stable for any order of interpolation of the displacement field and the multiplier even for nearly incompressible materials. Stress approximations are recovered by a local post-processing for both displacement and stress field using the multiplier approximation obtained by the SHDG method. The residual form of the equilibrium equation at the element level is added in this local post-processing leading to optimal rates of convergence of the stress approximations in H(div) norm for compressible and nearly incompressible elasticity problems with uniform meshes.

In the next section the equations that define the mechanical problem including the condition of nearly compressibility is presented. In Section 3, we introduce the SHDG formulation for the nearly incompressible linear elasticity problem based on displacement fields and the multiplier. In Section 4, we show Linear elasticity, Discontinuous Galerkin method, Stabilization Hybrid method, Locking free aspects related to the implementation and introduce a local post-processing technique to recover the pressure field. In Section 5, some numerical examples showing the h-convergence rates for the displacement field and the result of the post-processing technique for the recovery of the stress field convergence rates. Finally we have some concluions in Section 6.

2 THE MODEL PROBLEM

Let Ω in R² an open bounded domain with piecewise Lipschitz boundary Γ_D = ∂Ω of an elastic body subjected to external force f ∈ [L ²(Ω)]² where Γ_D is the Dirichlet boundary. The kinematical model of linear elasticity problem in two dimensions consists in finding a displacement vector field u satisfying

$\begin{array}{cc}-div\sigma \left(\mathit{u}\right)=\mathit{f}& in\Omega ,\\ \sigma \left(\mathit{u}\right)=\mathbb{D}\epsilon \left(\mathit{u}\right)& in\Omega ,\\ \mathit{u}=\mathit{g}& on{\Gamma }_D,\end{array}$ (2.1)

where σ (u) is the symmetric Cauchy stress tensor, ε(u) = $\frac{1}{2}$ (grad u+grad u_T ) is the linear strain tensor, D is the isotropic elasticity tensor and g is a given boundary displacement.

For linear, homogeneous and isotropic material σ (u) is given by σ (u) = Dε(u) = 2με(u) + λ(tr ε(u))I, where tr ε(u) = div u, I is the identity tensor and λ and μ are called the Lamé parameters which are given in terms of elasticity modulus E and Poisson ratio ν by

$\lambda =\frac{E\nu }{\left(1+\nu \right)\left(1-2\nu \right)}and\mu =\frac{E}{2\left(1+\nu \right)}.$

Nearly incompressible materials are modeled by Poisson ratio when ν → 1/2. It is well known that in this limit classical Galerkin FEM approximations lead to volumetric locking.

3 THE SHDG FORMULATION

In this section we present the Stabilized Hybrid Discontinuous Galerkin (SHDG) formulation for the linear elasticity problem based in¹⁶. Let T _h be a finite element mesh on Ω consisting of triangles ${\left\{K\right\}}_{K\in T_h}$, E _h = {e: e is an edge of K for all K ∈ T _h } denote the set of all edges of all elements K of the mesh T _h , and $E_h^0=\left\{e\in E_h:eisaninterioredge\right\}$ is the set of interior edges. In this primal hybrid formulation the multiplier λ is identified as the trace of the displacement field. That is, u: λ = u|e on each edge e ∈ E _h .

To construct the primal hybrid formulation we introduce the following broken function space for the displacement field

$V_h^k=\left\{v_h\in L^2\left(\Omega \right):v_h|K\in {\left[S_k\left(K\right)\right]}^2\forall K\in T_h\right\},$ (3.1)

where S _k (K ) = P _k (K) (the space of polynomial functions of degree at most k in both variables) and for the Lagrange multiplier the discontinuous space is

$M_h^l=\left\{\lambda \in L^2\left(E_h\right):\lambda {|}_e={\left[pl\left(e\right)\right]}^2,\forall e\in E_h^0\right\}$ (3.2)

with polynomial functions p _l (e), of degree at most l, on each edge e. In these finite dimension spaces, the SHDG method is formulated as:

Find the pair $\left[u_h,{\lambda }_h\right]\in V_h^k\times M_h^l$ such that, for all $\left[v_h,{\mu }_h\right]\in V_h^k\times M_h^l$

$$\begin{gathered} \sum _{K\in T_h}{\int }_K\mathbb{D}\epsilon \left(u_h\right):\epsilon \left(v_h\right)dx-\sum _{K\in T_h}{\int }_{\partial K}\mathbb{D}\epsilon \left(u_h\right)n_K·\left(v_h-{\mu }_h\right)ds \\ -\sum _{K\in T_h}{\int }_{\partial K}\mathbb{D}\epsilon \left(v_h\right)n_K·\left(u_h-{\lambda }_h\right)ds \\ +2\mu \sum _{K\in T_h}\beta {\int }_{\partial K}\left(u_h-{\lambda }_h\right)·\left(v_h-{\mu }_h\right)ds \\ =\sum _{K\in T_h}{\int }_{K}f·v_{h}dx, \end{gathered}$$ (3.3)

where the stabilization parameter β depends on the mesh parameter h. Here, we have considered the following definition for this stabilization parameter β = $\frac{{\beta }_0}{h}$, ∀e ∈ E _h with β₀ > 0. The influence of β₀ on the stability and accuracy of this primal hybrid formulation is analyzed in¹⁶ for the compressible elasticity problem.

4 COMPUTATIONAL IMPLEMENTATION

Considering that v_h , belonging to the broken function space $V_h^k$, is defined independently on each element K ∈ T _h , we observe that equation (3.3) can be split into a set of local problems defined on each element K coupled to the global problem defined on E _h , as follow:

Local problems: Find u_h|K ∈ $V_h^k$ (K) = $V_h^k$|K, such that, for all v_h|K ∈ $V_h^k$ (K),

$$\begin{gathered} {\int }_K\mathbb{D}\epsilon \left(u_h\right):\epsilon \left(v_h\right)dx-{\int }_{\partial K}\mathbb{D}\epsilon \left(u_h\right)n_K·v_{h}ds \\ -{\int }_{\partial K}\mathbb{D}\epsilon \left(v_h\right)n_K·\left(u_h-{\lambda }_h\right)ds \\ +2{\int }_{\partial K}\beta \left(u_h-{\lambda }_h\right)·v_{h}ds={\int }_{K}f·v_{h}dx, \end{gathered}$$ (4.1)

Global Problem: Find λ_h ∈$M_h^l$, such that, for all μ_h ∈$M_h^l$,

$$\begin{gathered} \sum _{K\in T_h}{\int }_{\partial K}\mathbb{D}\epsilon \left(u_h\right)n_K·{\mu }_{h}ds \\ -\sum _{K\in T_h}2\mu {\int }_{\partial K}\beta \left(u_h-{\lambda }_h\right)·{\mu }_{h}ds=0,\forall {\mu }_h\in M_h^l. \end{gathered}$$ (4.2)

Given that the multiplier of the proposed hybrid formulation is identified with the trace of the primal variable u on the element edges, for appropriate choices of β, we can always eliminate the degrees-of-freedom of the primal variable u_h at the element level in favor of the degrees-of-freedom of the multiplier leading to a global system in the multiplier only. Note that we can adopt any (l) order of continuous or discontinuous interpolation functions to Lagrange multiplier λ_h independently of the (k) order adopted for the primal variable u_h . Here, we consider only discontinuous interpolations for the Lagrange multiplier.

4.1 Stress and displacement local post-processing

In most engineering applications, stresses are the variables of main interest. Classically, in displacement finite element formulation, stresses are computed indirectly using the displacement approximation and the constitutive equation only. With this classical approach, the stress approximation is given by

${\sigma }_h=\mathbb{D}\epsilon \left(u_h\right)$ (4.3)

which converges at best with the following rates in L² and H(div) norms:

${||\sigma -{\sigma }_h||}_{L^2}={||\mathbb{D}\epsilon \left(u\right)-\mathbb{D}\epsilon \left(u_h\right)||}_{L^2}=C{h}^k,$ (4.4)

${||\sigma -{\sigma }_h||}_{H\left(div\right)}=C{h}^{k-1}.$ (4.5)

As an improved alternative to compute stress approximations, we propose a local post-processing consisting in solving at each element K ∈ T _h the local problem:

$\begin{array}{l}-\text{div}\sigma \left(\text{u}\right)=\text{f in}K,\\ \text{𝔸}\sigma \left(\text{u}\right)=\epsilon \left(\text{u}\right)\text{in}K,\\ \text{u}={\lambda }_h\text{on}\partial K,\end{array}$ (4.6)

in stress and displacement fields, in which A = D^-1 and λ_h is given by the solution of the global problem. Stress and displacement approximations [σ_pp , u_pp ] for [σ, u], solution of (4.6), are obtained in the finite dimension spaces

${\mathbb{W}}_h^k\left(K\right)=\left\{{\tau }_{ij}\in S_k\left(K\right),{\tau }_{ij}={\tau }_{ji},i,j=1,2\right\}$ (4.7)

and

$V_h^k\left(K\right)=\left\{v_i\in S_k\left(K\right),i=1,2\right\},$ (4.8)

respectively, considering the following residual form on each element K ∈ T _h .

Given λ_h , find $\left[{\sigma }_{pp|K},u_{pp|K}\right]\in {\mathbb{W}}_h^k\left(K\right)\times V_h^k\left(K\right)$, such that

$a_{pp}\left(\left[{\sigma }_{pp},u_{pp}\right],\left[{\tau }_h,v_h\right]\right)=f_{pp}\left(\left[{\tau }_h,v_h\right]\right)\forall \left[{\tau }_{h|K},v_{h|K}\right]\in {\mathbb{W}}_h^k\left(K\right),$ (4.9)

with

$\begin{array}{c}{a}_{pp}\left(\left[{\sigma }_{pp},u_{pp}\right],\left[{\tau }_h,v_h\right]\right)={\int }_K\mathbb{A}{\sigma }_{pp}:{\tau }_{h}dx+{\int }_K·div{\tau }_{h}dx\\ +{\int }_{K}div{\sigma }_{pp}·v_{h}dx+\frac{\delta }{2\mu }{\int }_{K}div{\sigma }_{pp}·div{\tau }_{h}dx\\ +2\mu {\int }_{\partial K}\beta u_{pp}·v_{h}ds\end{array}$ (4.10)

$\begin{array}{c}{f}_{pp}\left(\left[{\tau }_h,v_h\right]\right)={\int }_{\partial K}{\lambda }_h·{\tau }_h{n}_{K}ds-\frac{\delta }{2\mu }{\int }_{K}f·div{\tau }_{vh}dx\\ -{\int }_{K}f·v_{h}dx+2\mu {\int }_{\partial K}\beta {\lambda }_h·v_{h}ds\end{array}$ (4.11)

For appropriate choices of the stabilization parameters δ, we have observed the following convergence rate for the post-processed stress:

${||\sigma -{\sigma }_{pp}||}_{H\left(div\right)}=C{h}^k$ (4.12)

which is one order higher than that observed for |σ - σ _h |_H(div).

5 NUMERICAL RESULTS

In this section the behavior of the proposed formulation is tested to plane-strain problem, defined on square domain Ω = (0, 1) × (0, 1) with homogeneous boundary conditions, considering the elasticity modulus E = 1 and forcing term:

$f_1\left(x,y\right)=\left[2\nu \left(2\mu +\lambda \right)-\left(\mu +\lambda \right)\right]\sin \left(\pi x\right)\cos \left(\pi y\right),$ (5.1)

$f2\left(x,y\right)=\left[2\nu \left(2\mu +\lambda \right)-\left(3\mu +\lambda \right)\right]\sin \left(\pi y\right)\cos \left(\pi x\right),$ (5.2)

The approximate solutions have been obtained using triangular elements P _k - p _l where k and l denote the degree of Lagrangian polynomial space for the displacement field and the Lagrange multiplier, respectively. The post-processed stress approximations are recovered with Lagrangian polynomials of degree k + 1. Next, we present the incompressible elasticity problem, using uniform meshes, comparing the stress approximations σ _pp obtained from the local post-processing (4.9), with σ_h , obtained by the constitutive equation (2.1).

5.1 h-convergence of displacement field

Figures (1) and (2) present h-convergence studies for u_h and λ_h of SHDG approximations with equal order (l = k) in L ²-norm. In these experiments, we consider the parameters β = 20 and ν = 0.3, 0.49, 0.499, 0.4999, 0.49999 (nearly incompressible) with triangular elements, P ^{1 -} p ¹ in Figure (1) and P ^{2 -} p ² in Figure (2). Optimal rates of convergence are observed in L ²(Ω)-norm, with O(h ²) to P ^{1 -} p ¹, and O(h ³) to P ^{2 -} p ² with little loss of precision for nearly incompressible case.

Figure 1:
h-convergence study for u_h and λ_h of primal SHDG approximations with ν = 0.3, 0.49, 0.499, 0.4999, 0.49999.

Figure 2:
h-convergence study for u _h and λ_h of primal SHDG approximations with ν = 0.3, 0.49, 0.499, 0.4999, 0.49999.

5.2 Convergence rates for the stress tensor

Figures (3)-(6) present the convergence rates in L ² and H (div) norm for the stress tensor approximations when calculated using the constitutive equation, σ_h , or by the post-processing technique, σ _pp , as presented in Section 4.1, respectively. In this study, we consider the parameters β = 20, ν = 0.3, 0.49, 0.499, 0.4999, 0.49999 (nearly incompressible) and the local stabilization parameter δ = 1 for all elements K ∈ T _h . We use triangular elements P _k - p _k with (k = 1, 2), adopting same order of interpolations k for displacement field and the Lagrange multiplier. For the stress field, we use k +1. As expected for k = 1, σ_h , calculated from the constitutive equation does not converge in H(div) norm and for k = 2 converges with order O(h ¹). Improved rates of convergence O(h ¹), for k = 1, O(h ²), for k = 2 are obtained for σ _pp , calculated by the proposed post-processing technique. Convergence rates are observed for all approximations σ _pp , independent of ν, using triangular elements. Same accuracies are obtained for ν = 0.49, 0.499, 0.4999, 0.49999.

Figure 3:
h-convergence for σ _pp and σ_h aproximations, in L ² norm using k + 1 = 2 and ν = 0.3, 0.49, 0.499, 0.4999, 0.49999.

Figure 4:
h-convergence for σ _pp and σ_h aproximations, in H (div) norm using k + 1 = 2 and ν = 0.3, 0.49, 0.499, 0.4999, 0.49999.

Figure 5:
h-Convergence for σ _pp and σ_h approximations in L ² norm using k + 1 = 3 and ν = 0.3, 0.49, 0.499, 0.4999, 0.49999.

Figure 6:
h-Convergence for σ _pp and σ_h approximations in H (div) norm using k + 1 = 3 and ν = 0.3, 0.49, 0.499, 0.4999, 0.49999.

6 CONCLUDING REMARKS

We have presented a Stabilized Hybrid Discontinuous Galerkin finite element formulation for incompressible linear elasticity problems, based in¹⁶. The method was developed from a primal hybrid formulation, with the Lagrange multiplier identified with the trace of the displacement field on the edges of the elements, leading to a set of local problems defined at the element level and a global problem in the multiplier only. The formulation preserves the main properties of the corresponding DG method, but with reduced computational complexity. Numerical experiments confirm the locking free ability of the method for nearly incompressible problems. Accurate approximations and optimal rates of convergence in L ² norm are obtained for triangular elements on uniform meshes.

A local post-processing based on residual forms of the constitutive and equilibrium equations at the element level, and using the multiplier as a kinematic boundary condition, is presented to recover stress approximations with improved stability, accuracy and robustness. Numerical experiments show that, compared to the standard post-processing based on the constitutive law, it retrieves optimal convergence rates in H (div) norm and improves accuracy using triangular elements.

ACKNOWLEDGEMENTS

The authors thank the Brazilian Research Council (CNPq) and the Rio de Janeiro State Foundation (FAPERJ) for their financial support to the development of this work.

REFERENCES

• ¹
  L.P. Franca, T.J.R. Hughes, A.F.D. Loula & I. Miranda. A new family of stable elements for nearly incompressible elasticity based on a mixed Petrov-Galerkin finite element formulation. Numer. Math., 53 (1988), 123-141.
• ²
  L.P. Franca & R. Stenberg. Error analysis of some Galerkin least squares methods for the elasticity equations. SIAM J. Numer. Anal., 28(6) (1991), 1680-1697.
• ³
  R. Kouhia & R. Stenberg. A linear nonconforming finite element method for nearly incompressible elasticity and Stokes flow. Comput. Methods Appl. Mech. Engrg., 124 (1995), 195-212.
• ⁴
  M. Chiumenti, M. Cervera & R. Codina. A mixed three-field FE formulation for stress accurate analysis including the incompressible limit. Computer Methods in Applied Mechanics and Engineering, 283 (2015), 1095-1116.
• ⁵
  M. Vogelius. An analysis of the p-version of the finite element method for nearly incompressible materials. Numer. Math., 41 (1983), 39-53.
• ⁶
  J.W.B. Cockburn & D. Schötzau. Discontinuous Galerkin methods for incompressible elastic materials. Comput. Methods Appl. Mech. Engrg., 195 (2006), 3184-3204.
• ⁷
  P. Hansbo & M.G. Larson. Discontinuous Galerkin methods for incompressible and nearly incompressible elasticity by Nitsche’s method. Comput. Methods Appl. Mech. Engrg., 191 (2002), 1895-1908.
• ⁸
  T. Wihler. Locking-free adaptive discontinuous Galerkin fem for linear elasticity problems. Mathematics of Computation, 75(255) (2006), 1087-1102.
• ⁹
  Y. Chen, J. Huang, X. Huang & Y. Xu. On the Local Discontinuous Galerkin method for linear elasticity. Mathematical Problems in Engineering, (2010), 20 pp.
• ¹⁰
  D.A.D. Pietro & S. Nicaise. A locking-free discontinuous Galerkin method for linear elasticity in locally nearly incompressible heterogeneous media. Applied Numerical Mathematics, 63 (2013), 105-116.
• ¹¹
  B. Cockburn, J. Gopalakrishnan & R. Lazarov. Unified hybridization of discontinuous Galerkin, mixed, and continuous Galerkin methods for second order elliptic problems. SIAM J. Numer. Anal., 47(2) (2009), 1319-1365.
• ¹²
  S.C. Soon, B. Cockburn & H.K. Stolarski. A hybridizable discontinuous galerkin method for linear elasticity. International Journal of Numerical Methods in Engineering, 80 (2010), 1058-1092.
• ¹³
  N.C.B. Arruda, A.F.D. Loula & R.C. Almeida. Locally discontinuous but globally continuous Galerkin methods for elliptic problems. Comput. Methods Appl. Mech. Engrg., 255 (2013), 104-120.
• ¹⁴
  E.G. do Carmo, M.T.C.A. Fernandes & W.J. Mansur. Continuous/discontinuous Galerkin methods stabilized through transfer functions applied to the incompressible elasticity and to the Stokes problem. Computer Methods in Applied Mechanics and Engineering, 283 (2015), 806-840.
• ¹⁵
  G.J.B. dos Santos, A.F.D. Loula & A.J.B. dos Santos. A stabilized hybrid-mixed finite element formulation for the elasticity problems. In: Proceedings Series of the Brazilian Society of Applied and Computational Mathematics, 5 (2017).
• ¹⁶
  C.O. Faria, A.F.D. Loula & A.J.B. Santos. Primal stabilized hybrid and DG finite element methods for the linear elasticity problem. Computers & Mathematics with Applications, 68(4) (2014), 486-507.
• ¹⁷
  C.O. Faria, A.F.D. Loula & A.J.B. Santos. Stabilized hybridized finite element method for incompressible and nearly incompressible elasticity. In: Proceedings of the XXXIV Iberian Latin American Congress on Computational Methods in Engineering (CILAMCE), (Pirenópolis, Goiás, Brasil), (2013).
• ¹⁸
  C.O. Faria, A.J.B. Santos & A.F.D. Loula. Locking-free hybridized fem for incompressible linear elasticity problem. In: Proceedings of the XXXVI Iberian Latin American Congress on Computational Methods in Engineering (CILAMCE), (Rio de Janeiro, RJ, Brasil), (2015).

Publication Dates

• Publication in this collection
  Sep-Dec 2017

History

• Received
  25 Nov 2016
• Accepted
  06 Sept 2017