In this section, we shortly recapitulate interior penalty discontinuous Galerkin (DG) methods for non-linear solid mechanics. In contrast to the previously considered weakly conforming method, with DG methods the solution space is not constrained, but the strong force equilibrium is directly applied to bro-ken function spaces, see [18]. Applying integration by parts, we obtain a non-symmetric and singular variational setting, which can be stabilized by an additional non-symmetric term or by penalty, see [11, 12]. The symmetric version is applied to linear models in solid mechanics in [65, 97, 100, 146, 226].
Further formulations are applied to plasticity in [2, 147], for finite elasticity in [124, 211] and for failure in [159]. The use of discontinuous Petrov-Galerkin (DPG) methods for linear elasticity is studied in [37, 86, 126]. Hybrid methods for linear models are considered in [199], and divergence-free hybrid methods
-0.6 -0.5 -0.4 -0.3 -0.2 -0.1 0
100 1000 10000 100000
displacement, u2(A) [mm]
number of degrees of freedom Q1
Q2 WCM Q3-P1
Figure 5.12: Convergence plot of tip-displacement for Cook’s membrane.
are analyzed in [143]. An extensive overview on existing DG methods is given in [64, 184].
At first, purely primal methods are shown and then a hybrid DG method as introduced in [235] is presented. We directly introduce the methods for the case of non-linear elasticity, but, for simplicity of the presentation, we only consider the two-dimensional case on quadrilateral elements with homogeneous Dirichlet data ub = 0. For an H(div)-conforming stress and a discontinuous test function vh, the 1st Piola–Kirchhoff stress tensor P = P(Id + Gradu) satisfies
X
τh∈Th
Z
τh
P: Gradvh dx=Z
Γ
{P}N>+[vh] dγ+Z
ΓN
bt>vh dγ, (5.7) where Γ = SF∈FhF¯ denotes the skeleton of the triangulation and the averaged stress
{P}=
(P++P−)/2 for F 6⊂∂Ω P, for F ⊂∂Ω,
is defined analogously to the average of the normal stress (5.2). We remark that [vh] =0on ΓNby definition of the jump (5.1). Based on this identity, we present different variants of the DG method.
5.3.1 The discontinuous Galerkin method
As with the weakly conforming method, for DG methods we choose local poly-nomial ansatz spacesVh(τh), and define the broken (discontinuous) ansatz and
test space
VDGh =uh ∈L2(Ω)2 : uh|τh ∈Vh(τh) . (5.8) Then, based on the weak force equilibrium (5.7), we define the consistent non-linear DG form
aDGh (u,v) = X
τh∈Th
Z
τh
P(Id + Gradu) : Gradvdx
−
Z
Γ
{P(Id + Gradu)}N>+[v] dγ,
which is linear in the second argument, and for the exact solutionu we have aDGh (u,vh) = f(vh), vh ∈VDGh .
Since the linearization of the DG form is not positive definite for discontinuous ansatz functions, a stabilization is needed. Now, we discuss several possibilities (noted by an index ∗) to achieve stability by choosing an additional possibly non-linear form s∗h(·,·), linear in the second argument, so that in the case of small deformations
aDGh (uh,uh) +s∗h(uh,uh)>0, uh ∈VDGh ,uh 6=0.
Then, for sufficiently smallbt, there exists a unique solution uh ∈ VDGh , satis-fying
aDGh (uh,vh) +s∗h(uh,vh) = f(vh), vh ∈VDGh . (5.9) In order to preserve asymptotic consistency, s∗h(u,vh) → 0 for h → 0 is re-quired for the solutionuof the continuous problem and for discontinuous test functions vh. In general for DG methods even full consistency s∗h(u,vh) = 0 is achieved.
The convergence properties and the robustness of DG methods depend on the choice of the ansatz space and the stabilization. The most popular sta-bilization is the Nitsche-type penalty, which yields the incomplete interior penalty method (IIPDG)
sIIPDGh (uh,vh) =Z
Γ
θ
h[uh]>[vh] dγ, (5.10) introduced in [167] with the mesh-independent parameter θ >0. Note that in general this parameter depends on the polynomial degree of the local ansatz spaces Vh(τh). While this choice preserves the asymptotic order of the
ap-proximation as well as the condition number, good numerical results have also been observed with alternative h-scalings of the penalty parameter. Hence for some examples we allow θ to beh-dependent.
The symmetric interior penalty method (SIPDG) is obtained by minimizing the energy (see also Section 1.2.4), extended with consistency and stabilization terms, i.e., In the variational form this additionally yields the adjoint consistency term
sSIPDGh (uh,vh) =sIIPDGh (uh,vh)−sadjh (uh,vh), (5.11) where
sadjh (uh,vh) = Z
Γ[uh]>{∂FP(Id + Graduh) Gradvh}N+ dγ.
In the case of linear elasticity, the system matrix is symmetric and in this case one obtains optimal L2(Ω) convergence, see [11].
The case ofsadjh without any Nitsche-type stabilization is considered in [18], where we observe robustness and obtain a non-symmetric linearization in (5.9).
To analyze robustness for linear materials in the nearly incompressible case, penalty terms adapted to the large value of λLam´e are considered in [97, 100].
In case of bilinear ansatz spacesVh(τh), a reduced integration on the element faces together with a simple Nitsche-type penalty was introduced in [20], which shows improved convergence in case of locking. We consider the bilinear form
aRIDGh (uh,vh) = X where RRI denotes the mid-point rule approximating the integral on every element face.
While robustness and even a full convergence analysis is provided for many different DG variants in the linear setting, it remains important to reduce the computational cost. For all these methods the system matrix has the same structure with entries connecting all degrees of freedom of neighboring
ele-ments. Hybrid methods aim to reduce the numerical expense while preserving the robustness properties as we have already observed for the weakly conform-ing method.
Remark 5.3.1. Optimal L2 norm error estimates can be shown for suffi-ciently regular linear problems under the assumption of adjoint consistency, which follows from the symmetry of the methods, see, e.g., [64, Chapter 4.2.4].
While this works for SIPDG, for the further methods examples exist with a sub-optimal L2 norm, see [184, Chapter 2.8].
Also L∞ norm error estimates are, to the best of our knowledge, currently only available for the symmetric interior penalty method, see [48].
5.3.2 A low-order hybrid discontinuous Galerkin method with conforming traces
Here we shortly recapitulate the hybrid DG method as proposed in [235]. For the hybrid DG method we use a second discrete space
VΓ,h ⊂nuΓ,h∈C0(Γ)2 : uΓ,h=0 on ΓD
o,
for the approximation of the displacement vector on the skeleton and use traces of conforming bilinear ansatz functions, i.e., face-wise linear skeleton functions in VΓ,h. In the interior of the elements, we use the discontinuous space VDGh with linear ansatz functions Vh(τh), see (5.8). For a given linear skeleton function uΓ,h ∈ VΓ,h, we consider the linear volume approximation uh = Πlinh uΓ,h ∈ VDGh , defined in every quadrilateral element by the linear projection
Z h
∂τh
Πlinh u>Γ,hvh =Z h
∂τh
uΓ,h> δuh, vh ∈Vh(τh),
where the boundary integral on ∂τh is approximated by a trapezoidal quadra-ture rule
Z h
∂τh
f =X4
j=1
wjf(Xτh,j)
using the element corners Xτh,1,Xτh,2,Xτh,3,Xτh,4 and weights wj > 0. Since the strain approximation Fh = Id + Graduh is constant in τh, we obtain DivP(Fh) =0 and
Z
τh
P(Fh) : Grad(vh) dx=Z
∂τh
P(Fh)N>vh dγ≈
Z h
∂τh
P(Fh)N>vh.
This defines the non-symmetric hybrid form aHDGh (uΓ,h,vΓ,h) = X
τh∈Th
Z h
∂τhPId + Grad Πlinh uΓ,h
N>vΓ,h.
Adapting the standard Nitsche-type penalty term (5.10) to the hybridization yields
sHDGh (uΓ,h,vΓ,h) = X
τh∈Th
θ h
Z h
∂τh
Πlinh uΓ,h−ulinΓ,h>vΓ,h. Then, the hybrid solution uΓ,h ∈VΓ,h is computed by
aHDGh (uΓ,h,vΓ,h) +sHDGh (uΓ,h,vΓ,h) =Z
ΓN
bt>vΓ,h dγ, vΓ,h∈VΓ,h. The system matrix is non-symmetric with the sparsity pattern of conforming bilinear elements and all computations can be performed on element level which makes it convenient to include it in standard finite element codes.
We note that the hybrid DG method is based on the incomplete inte-rior penalty method, hence straightforward L2 estimates based on an Aubin–
Nitsche argument cannot be expected, see Remark 5.3.1.