Local adaptive refinement method applied to solid mechanics

Loïc Daridon, Eric Delaume, Yann Monerie, Frédéric Perales

► To cite this version:

Loïc Daridon, Eric Delaume, Yann Monerie, Frédéric Perales. Local adaptive refinement method applied to solid mechanics. Applied and Computational Mechanics, 2020, 14 (2), 10.24132/acm.2020.570. hal-03053926

HAL Id: hal-03053926

https://hal.science/hal-03053926

Submitted on 11 Dec 2020

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Local *h*-adaptive refinement method applied to solid mechanics

L. Daridon<sup>*a,c,*</sup>, E. Delaume<sup>*b,c*</sup>, Y. Monerie<sup>*a,c*</sup>, F. Perales<sup>*b,c*</sup>

<sup>*a*</sup>LMGC, University of Montpellier, CNRS, Montpellier, France

<sup>*b*</sup>Institut de Radioprotection et de Sûreté Nucléaire (IRSN), Cadarache

<sup>*c*</sup>MIST University of Montpellier, CNRS, IRSN, France

Received 13 November 2019; accepted 9 November 2020

## Abstract

A good spatial discretization is of prime interest in the accuracy of the Finite Element Method. This paper presents a new refinement criterion dedicated to an *h*-type refinement method called Conforming Hierarchical Adaptive Refinement MethodS (CHARMS) and applied to solid mechanics. This method produces conformally refined meshes and deals with refinement from a basis function point of view. The proposed refinement criterion allows adaptive refinement where the mesh is still too coarse and where a strain or a stress field has a large value or a large gradient. The sensitivity of the criterion to the value or to the gradient can be adjusted. The method and the criteria are validated through 2-D test cases. One limitation of the *h*-adaptive refinement method is highlighted: the discretization of boundary curves.

© 2020 University of West Bohemia. All rights reserved.

*Keywords:* conformity, hierarchy, adaptivity, refinement method

## 1. Introduction

Higher accuracy in Finite Element Method (FEM) can be obtained by adjusting spatial resolution of domain discretization. A solution is to perform a uniform mesh refinement until finest singularities are captured. This approach can lead to a cumbersome calculation if the studied structure exhibits large dimensions compared to the size of some local heterogeneity that drives the overall toughness, e.g. a crack across the thickness of a thin rod or a plate. An alternative is to perform local adaptive refinement. For applications we have in mind, we are looking to divide by about 10 the CPU time with a local adaptive refinement compared to a uniform mesh refinement at same finite element accuracy. Four main refinement methods can be listed:

* the *r*-adaptive method consists in moving nodes of an initial mesh in the region of interest [12],

* the *s*-adaptive method, where the initial mesh is overlaid with finer local meshes [11],

* the *p*-adaptive method consists in increasing basis function degree where it is needed [4, 24],

* the *h*-adaptive method, where elements are split in some regions to decrease their size and locally increase the mesh density [4].

<footnote>*Corresponding author. Tel.: +33 467 149 651, e-mail: loic.daridon@umontpellier.fr
https://doi.org/10.24132/acm.2020.570</footnote>
The major disadvantage of the *p*-method is its difficult use in case of geometrical singularities. Splitting elements with a *h*-type method generally leads to non-conformities. These methods have been also combined to increase their efficiency: *hp*-method [2, 4] or *hr*-method [14, 15] among others. A convenient mesh refinement method has both to keep good element shape quality in order to avoid inaccurate numerical predictions due to elongated/flat/degenerated elements and to lead to conformal refined meshes as required by FEM. Several methods have been developed to overcome this issue of non-conforming *h*- and *s*-methods. For example, to restore continuity of finite element space, hanging (constrained) nodes, which are equipped with the Degrees Of Freedom (DOFs) expressed with DOFs of the active nodes, are added in suitable region. Another way is to use Lagrangian multipliers or penalty methods. An alternative is also to split extra elements in order to conveniently connect refined zone to the coarse one. A method based on this idea has been developed and extended to hexahedron in [19] to keep the number of degrees of freedom constant (*r*-method).

In this paper, a *h*-type refinement method called Conforming Hierarchical Adaptive Refinement MethodS (CHARMS) is studied. This method, initially developed in the framework of fluid mechanics [6, 16] and more recently applied to solid mechanics [10, 13], allows both (i) to preserve element shape quality while producing a conforming mesh and (ii) to perform local adaptive refinement and unrefinement. This method deals with refinement from a basis function point of view and not from a geometrical point of view. The general idea is to construct a hierarchy of nested meshes, where at a given level of refinement each basis functions can be expressed as a linear combination of basis functions associated to the immediate finer or coarser mesh. This refinement method is here used with new refinement criterion and is validated through 2-D test cases. Notice that CHARMS method is also valid in three dimensions and preserves conformal finite elements at all stages of refinement even with non-conformal geometric nodes present.

In Section 2, notations, local adaptive refinement procedure and refinement pattern are introduced. Although unrefinement can be easily performed, only mesh refinement is explained and validated in this study. A new refinement criterion is proposed in Section 3. Contrary to usual refinement criteria based on error estimations, this criterion is based on a coupling between intensity of a mechanical field, its gradient and local characteristic size of the mesh. Validation of the method applied to FEM in solid mechanics is exposed through 2-D test cases and compared to fine uniform meshes. In the last section, a limitation of the method is pointed out: the discretization of boundary curves.

## 2. Conforming, Hierarchical, Adaptive Refinement Methods

CHARMS is here applied to conforming Lagrangian elements and based on the refinement/unrefinement of basis functions. This approach relies on the nested relation between approximation spaces $X_0 \subset \dots \subset X_J$, $J \ge 1$, generated by a set of basis functions $B_j$, $j \in [0, J]$, where the coarsest space is $X_0$ and the finer spaces have increasing indexes up to $X_J$. In this setting, adding/removing basis functions in $B_j$ allows both to refine/unrefine the approximation space $X_j$ and to ensure a linear independency between the set of basis functions in $B_j$ [6, 13, 16, 17].

Given that $X_j \subset X_{j+1}$, all basis functions in $B_j$ can be expressed as linear combination of basis functions in $B_{j+1}$. These linear combinations define children-parents relationships between two consecutive refinement levels. To unrefine a child basis function (level $j+1$), this function is replaced by a linear combination of parents basis functions (level $j$).
In order to perform refinement using CHARMS, a refinement pattern satisfying compatibility requirements has to be defined to generate conformal meshes from a finite element point of view. This pattern must be constituted of the same type of elements than the coarse one. The edges of element have thus to be split in several equal edges, here in two equal edges:

* $\mathbb{P}_1$-triangular elements are split in 4 identical triangles, Fig. 12,

* $\mathbb{Q}_1$-quadrangle elements are split in 4 identical quadrangles, Fig. 13,

* $\mathbb{Q}_1$-hexahedron elements are split in 8 identical hexahedra, Fig. 14,

* $\mathbb{P}_1$-tetrahedron elements are based on the strategy given in [10] which ensure that the shape quality of refined tetrahedra does not degenerate.

The key element of the CHARMS method is the refinement equation that links the coarse basis functions to the refined basis functions. This equation depends on the refinement pattern. It can be generalized and written as follows for a basis function $\hat{\varphi}_i$ related to a reference element $\hat{K}$:

$$ \forall i \in [0, \hat{N} - 1], \quad \hat{\varphi}_i^{(0)} = \sum_{k=0}^{\hat{N}^{(1)} - 1} \hat{\beta}_{ik} \hat{\varphi}_k^{(1)} \quad \text{with} \quad \hat{\beta}_{ik} = \hat{\varphi}_i^{(0)} (\hat{a}_k^{(1)}), \quad (1) $$

where superscripts (0) and (1) denote the reference element and the refinement pattern respectively, $\hat{N}$ is the number of nodes, $\hat{\varphi}_\bullet$ are the basis functions and $\hat{a}_k^{(l)}$ is the node linked to the basis function $\hat{\varphi}_k^{(l)}$.

## Refinement algorithm

Let $B^* = \{\varphi_k^{(j)}, \text{ where } j \in [0, J]\}$ be a multi-level basis. A refining step of $B^*$ consists in creating a new multi-level basis $B$ by replacing basis function $\varphi_k^{(j)}$ by all children basis functions $\varphi_l^{(j+1)}$. In Fig. 1 the red and blue zones correspond to an area approximated respectively by $B^*$ and $B$. At the first step of refinement, in the case of the refinement of the square reference element depicted in Fig. 13 the relations between basis function $\varphi_k^{(0)}$ and their children basis functions $\varphi_k^{(1)}$ are given by:

$$ \left\{ \begin{aligned} \hat{\varphi}_0^{(0)} &= \hat{\varphi}_0^{(1)} + \frac{1}{2} (\hat{\varphi}_1^{(1)} + \hat{\varphi}_3^{(1)}) + \frac{1}{4} \hat{\varphi}_4^{(1)}, \\ \hat{\varphi}_1^{(0)} &= \hat{\varphi}_2^{(1)} + \frac{1}{2} (\hat{\varphi}_1^{(1)} + \hat{\varphi}_5^{(1)}) + \frac{1}{4} \hat{\varphi}_4^{(1)}, \\ \hat{\varphi}_2^{(0)} &= \hat{\varphi}_8^{(1)} + \frac{1}{2} (\hat{\varphi}_5^{(1)} + \hat{\varphi}_7^{(1)}) + \frac{1}{4} \hat{\varphi}_4^{(1)}, \\ \hat{\varphi}_3^{(0)} &= \hat{\varphi}_6^{(1)} + \frac{1}{2} (\hat{\varphi}_3^{(1)} + \hat{\varphi}_7^{(1)}) + \frac{1}{4} \hat{\varphi}_4^{(1)}. \end{aligned} \right. \quad (2) $$

For the initial step, the basis functions are thus given by $B = \{\varphi_k^{(0)}, k \in \{1, 2, 3, 4\}\}$ and after the first refinement step there are given by $B^* = \{\varphi_k^{(1)}, k \in \{1, 2, 3, 4, 5, 6, 7, 8\}\}$.
<figure><img src="image_3.png" /><figcaption>Fig. 1. (Left) Initial mesh. The red zone represents the area approximated by the basis functions $B^*$. (Right) Refined mesh. The blue zone represents the area approximated by the basis functions $B$</figcaption></figure>

### 3. Refinement criteria

In order to perform local adaptive refinement, refinement criteria have to be implemented in order to select basis functions that need to be refined. Classically refinement criteria used in the literature are based on error estimation [3, 5, 18, 26, 27]. In the case of error-based adaptivity, two approaches are possible: (i) one uses two meshes [5, 9, 17] for error estimation; (ii) error is obtained via solution of a global problem [17, 18] or some local problems in elements [1, 24]. Additionally, the error-controlled adaptivity can be driven either by the interpolation error [7–9] or by the approximation, modeling and total errors [1, 21, 24, 25]. The estimators can either be recovery [16, 18] or residual ones [1, 21, 24, 25]. The former estimators were developed by Zienkiewicz [26]. In this case, the basic idea is to build a new stress field $\sigma^*$ from the stress field $\sigma^h$ obtained with the FEM, in interpolating $\sigma^*$ either by the same basis function used in the FEM [26] or by a polynomial in the so-called 'super-convergent points' in element patches of the domain [27]. Then, an error estimator $|\sigma_q^* - \sigma_q^h|$ is calculated for each element $q$. Whereas the second construction of $\sigma^*$ is more time consuming, it shows better results [27] and does not underestimate the error. The interpolation error control is usually associated with the iterative adaptation, while the approximation error control may lead to adaptation with fixed number of steps (three- or four-step for hp-approximation, for example) [20, 21, 25]. The latter approach needs availability of the convergence theory in order to relate the discretization parameters to the error level. The estimators can either be recovery [17, 18] or residual ones [1, 21, 24, 25]. There also exists another approach which is called the goal-oriented adaptivity [22], where the error estimation for the quantity of interest [1, 23] can be performed on the chosen part of the mesh.

All the above-mentioned criteria need additional calculations. In order to avoid these additional numerical costs, a new criterion, which belong to the iterative adaptation approach, have been developed. It is based both on the intensity of an ad hoc field $F$ and on its gradient taking also into account the refinement level. The general form of the criterion reads

$$ \forall i \in [1, N_{\text{dof}}], \quad C^{\text{crit}} > C = \left( \left| \frac{F_i}{F^{\text{max}}} \right| \right)^{\alpha} \left( \frac{\| \nabla F_i \|_2}{\| \nabla F \|_2^{\text{max}}} \right)^{1-\alpha} \left( \frac{1}{n} \right)^r, \quad (3) $$

where $C^{\text{crit}} \in [0, 1]$ is a threshold ($C^{\text{crit}} = 0$ means that all basis functions are refined), $F_i$ and $\nabla F_i$ represent the value of a field $F$ and of its gradient respectively for a given node $i$, $F^{\text{max}}$ and $\nabla F^{\text{max}}$ represent the largest value of $F$ and of its gradient respectively over the structure. $N_{\text{dof}}$ is the number of degrees of freedom, $|| \cdot ||_2$ represents the $L^2$-norm and the power $r$ stands for
the refinement level. For applications we have in mind, the term $|F_i/F^{\max}|$ allows to refine areas submitted to high stress or strain where crack can be initiated, while the term $||\nabla F_i||_2/||\nabla F||_2^{\max}$ allows to refine areas where stress or strain fields evolve rapidly in space, i.e. around a crack tip where already existing crack can propagate. The term $(1/n)^r$ allows to take into account the mesh size by favoring refinement on non refined elements where $n$ is the subdivision number of an edge, namely 2 in our study. The power $\alpha \in [0, 1]$ allows to weight the influence of the field intensity or its gradient. The studied mechanical problem induces the choice of the field $F$ and the power $\alpha$. In this work, we focus on the von Mises stress field $F = \sigma^{VM}$ and the von Mises strain field $F = \varepsilon^{VM}$:

$$ 
\begin{aligned}
\sigma^{VM} &= \frac{1}{2} \sqrt{(\sigma_I - \sigma_{II})^2 + (\sigma_I - \sigma_{III})^2 + (\sigma_{II} - \sigma_{III})^2}, \\
\varepsilon^{VM} &= \frac{1}{2} \sqrt{(\varepsilon_I - \varepsilon_{II})^2 + (\varepsilon_I - \varepsilon_{III})^2 + (\varepsilon_{II} - \varepsilon_{III})^2},
\end{aligned} 
\quad (4) $$

where $\sigma_i$ and $\varepsilon_i$ with $i \in \{I, II, III\}$ are, respectively, principal stress and principal strains.

### 3.1. Role of power $\alpha$

Although in the parametric studies presented, the mechanical behavior considered is independent of time, a fictitious time step is introduced in all numerical studies. This fictitious time step allows an arbitrary evolution of the value of the threshold $C^{\text{crit}}$ while the applied load remains constant. By this way, successive refinements can be observed and convergence analysis can be performed.

The first test case deals with a square hole in a square plate, Fig. 2. The material properties and numerical parameters are shown in Table 1 and 2. A constant vertical displacement is applied

<figure><img src="image_4.png" /><figcaption>Fig. 2. Geometry of the square plate with a square hole, dimensions in mm. The blue rectangle corresponds to zooms showed in Fig. 3</figcaption></figure>

<table><caption>Table 1. Square plate with square hole: material properties</caption><tr><td>Young modulus</td><td>$E = 300$ MPa</td></tr><tr><td>Poisson's ratio</td><td>$\nu = 0.4$</td></tr></table>

<table><caption>Table 2. Square plate with square hole: numerical parameters</caption><tr><td>Increment step</td><td>$\Delta t = 1$ s</td></tr><tr><td>Total number of increment</td><td>51</td></tr><tr><td>Element pattern</td><td>square</td></tr><tr><td>Mesh size</td><td>$h_{\text{ref}} = 0.416$ mm<br>$h_{\text{coarse}} = 4h_{\text{ref}}$</td></tr></table>
<figure><img src="image_5.png" /><figcaption>Fig. 3. Refined meshes for $C^{\text{crit}} = 0.75$ and (a) $\alpha = 0$, (b) $\alpha = 0.2$, (c) $\alpha = 0.4$, (d) $\alpha = 0.6$, (e) $\alpha = 0.8$, (f) $\alpha = 1$. The color code stands for the von Mises stress</figcaption></figure>

on the top of the specimen while the bottom remains fixed along the vertical direction. Two refinement steps are applied and six values of $\alpha$ are considered. The field of interest $F$ is the von Mises stress field $\sigma^{\text{VM}}$. This parametric study allows to estimate the role of the power $\alpha$. Fig. 3 shows the refinements obtained for different values of $\alpha$. These refinements are observed in an area of interest depicted by a blue rectangle in Fig. 2 and located in the vicinity of the corner of the square hole. As expected, when $\alpha$ is close to one, the refinement criterion (3) gives more weight to the intensity of the von Mises stress field than to its gradient and the mesh is refined
where large stress values are observed, i.e. within the ligament of the structural problem with respect to the load direction. In Fig. 3(f) mesh is refined everywhere but above the square hole where the stresses are screened. On the contrary, when $\alpha$ is close to zero the refinement criterion only concerns the gradient of the von Mises stress field and the refinement takes place where some stress concentration occurs, here at the corner of the square hole, Fig. 3(a). Extrapolating this situation to crack problems, the role of the power $\alpha$ could be summarized as follows: when $\alpha$ is close to zero, the criterion (3) induces refinement where new cracks can occur, when $\alpha$ is close to one, the mesh is refined at the crack tip of existing cracks.

## 3.2. CHARMS refinement vs homogeneous refinement

At a given number of elements, a mesh refined in convenient zones can give more precise results than a homogeneous mesh. To quantify this enhancement induced by the refinement criterion (3), a normalized error on the overall strain energy is calculated at each node of a mesh

$$e = \frac{|E^{\text{refi}} - E^{\text{ref}}|}{E^{\text{ref}}}, \qquad (5)$$

where $E^{\text{refi}}$ and $E^{\text{ref}}$ are the global strain energy on a refined mesh and on a reference mesh respectively. The studied case is the square plate with a square hole above detailed. The reference mesh is a homogeneous fine mesh of the structure with about 120 000 quadrangles whose mesh size $h_{\text{ref}}$ is given in Table 2.

As previously mentioned, the adaptive mesh refinement is obtained by decreasing the value of the threshold $C^{\text{crit}}$ (by step of 0.02) whilst the load remains constant.

Fig. 4 shows the evolution of this energy error (5) with respect to the number of elements for different values of $\alpha$. The thick black curve shows the energy errors for homogeneous refinement while the colored curves correspond to CHARMS refinement. Both the CHARMS

<figure><img src="image_6.png" /><figcaption>Fig. 4. Square plate with square hole: evolution of the global energy error with respect to the total number of elements, comparison between the CHARMS refinement and the homogeneous refinement for different values of $\alpha$. Stages of refinement are illustrated in Fig. 5</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>Fig. 5. Sketch of various meshes involved from stage 1 (left) to stage 5 (right) in Fig. 4. CHARMS refinements (top row), homogeneous refinements (bottom row). Stages 1 and 5 correspond to the same state</figcaption></figure>

refinement and the homogeneous refinement start from the same coarse mesh with about 7 000 elements (stage 1 in Fig. 4 and 5) and finish at the same fine mesh with about 120 000 elements (stage 5 in Fig. 4 and 5). At these stages 1 and 5, the energy error $e$ is thus the same for the CHARMS approach and the homogeneous one, with a larger error at stage 1 since it corresponds to the coarser mesh. During stages from 2 to 4, the meshes for the CHARMS refinement do not correspond to the meshes of the homogeneous refinement: the refinement takes place where the criterion (3) is reached. At the first, second and third refinements (stages 2, 3 and 4) the enhancement proposed by the CHARMS approach is significant and the energy error is about two times smaller than for the homogeneous refinement. However, no theoretical estimate of this energy error has been carried out and the estimates in Fig. 4 is pragmatic. In this example, only two CHARMS refinement steps are allowed. Then, when the mesh size is 4 times smaller than in the initial coarse mesh, the refinement process takes place somewhere else and the interest of the CHARMS refinement disappears.

Moreover, we can notice in Fig. 4 that the energy error decreases with $\alpha$. In this example, the high local stress concentration due to the square shape of the hole induces a predominant role of the gradient field in the criteria ($\alpha < 0.5$).

Two main conclusions can be drawn on this example: (i) since successive refinements are allowed, the proposed method can significantly reduce the finite element error (here in an energy sense) compared to homogeneous refinements; (ii) the power $\alpha$ has to be adapted to each problem and in practical calculations, when the parametric studies on $\alpha$ are not available, the value of $\alpha$ can be set to the mid value 0.5.

### 3.3. Role of the field $F$

This test case deals with bimaterial square plate (Fig. 6) constituted with stiff inclusions in softer matrix. Material properties and numerical parameters are shown in Table 3 and 4. A constant vertical displacement is applied on the top of the specimen while the bottom remains fixed along the vertical direction. In this study again two refinement steps are applied, the power $\alpha$ is fixed to 0.4 and two fields of interest $F$ are tested namely the von Mises stress Fig. 7(a) and the von Mises strain Fig. 7(b). As previously, the normalized error $e$ (5) is calculated at each node of the refined mesh to evaluate the accuracy of the method. Again, the value of the threshold $C^{\text{crit}}$ is also reduced by 0.02 at each fictitious time step.
<figure><img src="image_8.png" /><figcaption>Fig. 6. Geometry of the bimaterial plate subjected to tensile loading, dimensions are in mm. The inclusions are in grey and the areas of interest where relative energy errors are calculated are in blue</figcaption></figure>

<table><caption>Table 3. Bimaterial: mechanical properties</caption><thead><tr><th></th><th>Matrix</th><th>Inclusions</th></tr></thead><tbody><tr><th>Young modulus</th><td><i>E</i> = 1.5 GPa</td><td><i>E</i> = 15 GPa</td></tr><tr><th>Poisson's ratio</th><td><i>ν</i> = 0.2</td><td><i>ν</i> = 0.4</td></tr></tbody></table>

<table><caption>Table 4. Bimaterial: numerical properties</caption><tr><td>Increment step</td><td>Δ<i>t</i> = 1 s</td></tr><tr><td>Total number of increment</td><td>51</td></tr><tr><td>Element pattern</td><td>triangle</td></tr><tr><td>Mesh size</td><td><i>h</i><sub>ref</sub> = 0.5 mm<br><i>h</i><sub>coarse</sub> = 4<i>h</i><sub>ref</sub></td></tr></table>

Fig. 7 shows that the obtained refined mesh clearly depends on the choice of the field of interest used in the criterion (3). When the field of interest is set to the von Mises strain, the refinement is located around the inclusion where maximum equivalent strain is reached. When von Mises strain is used, the refinement is located inside the inclusions where maximum

<figure><img src="image_9.png" /><figcaption>Fig. 7. Refined meshes for: (left) Criteria using von Mises stress as *F*. The color code stands for the von Mises stress, (right) Criteria using von Mises strain as *F*. The color code stands for the von Mises strain</figcaption></figure>
<figure><img src="image_10.png" /><figcaption>Fig. 8. Bimaterial: evolution of the global energy errors for the refined mesh for von Mises stress and von Mises Strain criteria and homogeneous meshes depending on the total number of elements</figcaption></figure>

equivalent stress is reached. This ability of the criterion to induce the refinement zone of the mesh in relation to the quantity of interest chosen is very useful according to the studies to be carried out. For example, if we are interested in the fragmentation of inclusions or their decohesion, the choice of the quantity of interest is decisive. Fig. 8 shows the evolutions of energy normalized errors for the two criteria. For the two cases, the error rapidly decreases with respect to the number of elements and since the calculations are performed in linear elasticity the field of interest (stress or strain) has no influence on the FEM quality of the refined meshes.

## 4. Limitation: discretization of curved edge

In this section, a limitation of the CHARMS method is pointed out. Indeed, this method is based on the refinement of the basis functions (h-method) and not on remeshing. Geometric approximation of boundary curves are thus not improved during the refinement steps. We will highlight this limitation on the classical test case of an infinitely long cylinder loaded with an internal imposed displacement, with internal radius a and external radius b. The closed-form solution of this problem reads

$$ \forall r \in [a, b], \quad u_r(r) = Ar + \frac{B}{r} \quad \text{with} \quad A = \frac{u_0}{a(1-\xi)}, \quad B = \frac{a\xi u_0}{\xi - 1}, \quad \xi = \left(\frac{b}{a}\right)^2, \quad (6) $$

where the imposed internal displacement is $u_r(a) = u_0 \neq 0$ and the external displacement $u_r(b) = 0$. For this test case, a and b have been chosen so that $\xi = 2$.

The field F is the von Mises stress and the power $\alpha = 1$ because there is no stress concentration in this example. The refinement pattern used is the $P_1$-triangle element (Appendix A) and the maximum refinement level is 2. For the simulation, we consider a quarter of disc in plane strain assumption (Fig. 9). The material and numerical properties are presented in Table 5 and 6. The value of the threshold $C^{crit}$ is decreased by 0.1 at each time step so that at final time all basis functions are refined.

Knowing the closed-form solution, an exact relative error is calculated at the red dot exhibited in Fig. 9:

$$ \tilde{e} = \frac{|u_r^{\text{sim}} - u_r^{\text{exact}}|}{|u_r^{\text{exact}}|}, \quad (7) $$
<figure><img src="image_11.png" /><figcaption>Fig. 9. Geometry of the internal displacement imposed cylinder: $a = 100$ mm is the internal radius and $b = 200$ mm is the external radius. The red point is the point of interest</figcaption></figure>

<table><caption>Table 5. Hollow cylinder: material parameters</caption><tr><td>Young modulus</td><td>$E = 300$ MPa</td></tr><tr><td>Poisson's ratio</td><td>$\nu = 0.4$</td></tr></table>

<table><caption>Table 6. Hollow cylinder: numerical parameters</caption><tr><td>Increment step</td><td>$\Delta t = 1$ s</td></tr><tr><td>Total number of increments</td><td>10</td></tr><tr><td>Element pattern</td><td>Triangle</td></tr><tr><td>Mesh size</td><td>$h_{\text{ref}} = 0.625$ mm<br>$h_{\text{coarse}} = 4h_{\text{ref}}$</td></tr></table>

where $u_r^{\text{sim}}$ is the displacement calculated at the red dot with adaptive mesh refinement and $u_r^{\text{exact}}$ is the analytical displacement obtained with (6).

In this particular case of an initial coarse mesh of a curved boundary, the curves in Fig. 10 exhibits that the error is more important for a refined mesh using CHARMS method than for a mesh geometrically refined. Effectively, unlike what is classically exhibited using geometric

<figure><img src="image_12.png" /><figcaption>Fig. 10. Hollow cylinder: evolution of the exact relative error (7) calculated at the red node obtained with mesh refinement (blue curve) and with homogeneous meshes (black thick curve) with respect to the number of elements</figcaption></figure>
<figure><img src="image_13.png" /><figcaption>Fig. 11. Mesh obtained with the spatial discretization of a curved geometry. The real curvature is in black dots, the discretization obtained with CHARMS is in blue dots (the coarse mesh is the black square) and the discretization using a homogeneous mesh is in green dots. $e_d$ represents the spatial discretization error between the homogeneous mesh and the refined mesh</figcaption></figure>

mesh refinement and previous results (Fig. 4 and 8), we do not observe a monotonous decrease in error during successive refinement steps using the CHARMS method.

Fig. 11 illustrates the fact that geometric remeshing naturally induces a better edge interpolation. In fact using the CHARMS method coarser the initial mesh is, higher the geometrical error is while refining. The refinement method CHARMS cannot increase geometrical interpolation of an initial coarse mesh. Just as it is done in fluid mechanics for meshes to capture boundary layer behavior, it is important in solid mechanics for structures with large radii of curvature to use an initial mesh with good geometry interpolation before making refinement steps.

## 5. Conclusions

A local adaptive method using a basis function point of view is presented and validated in solid mechanics with elastic behavior. This method enables to refine elements without degenerating their shape quality and to implicitly handle the non-conformities. During a step of refinement the current basis functions are replaced by a linear combination of these functions based on a refinement pattern. A refinement criterion based on the value of an ad hoc field and its gradient has been developed. Compared to classical criteria, this approach is not time-consuming since it does not require any additional calculation but interpolated errors cannot be quantified *a priori*. The presented test cases show good results and put in evidence that the choice of the field of interest influences the region of refinement. At the end, one limitation of the *h*-adaptive refinement method is pointed out. If the specimen has a severe curvature, the initial mesh has to be fine enough because the CHARMS method, based on the interpolation of the shape functions, does not allow an improvement of the geometrical interpolation.

In future work, CHARMS will be applied in case of crack propagation where refinement and unrefinement can be performed adaptively.

## References

[1] Ainsworth, M., Oden, J. T., A posteriori error estimation in finite element analysis, New York, John Wiley & Sons, Inc., 2000. https://doi.org/10.1002/9781118032824
[2] Ainsworth, M., Senior, B., Aspects of an adaptive *hp*-finite element method: Adaptive strategy, conforming approximation and efficient solvers, International Journal for Numerical Methods in Engineering 150 (1997) 65–87. https://doi.org/10.1016/S0045-7825(97)00101-1

[3] Babuška, I., Rheinboldt, W. C., A-posteriori error estimates for the finite element method, International Journal for Numerical Methods in Engineering 12 (1978) 1597–1615.
https://doi.org/10.1002/nme.1620121010

[4] Babuška, I., Suri, M., The *p*- and *h*-p versions of the finite element method, an overview, Computer Methods in Applied Mechanics and Engineering 80 (1990) 5–26.
https://doi.org/10.1016/0045-7825(90)90011-A

[5] Biotteau, E., Combescure, A., Gravouil, A., Lubrecht, A., Stratégie multigrille et méthode de raffinement automatique à précision contrôlée pour la modélisation des structures élastoplastiques 3D en dynamique transitoire non-linéaire, 10e Colloque National en Calcul des Structures, Giens, France, 2011, pp. al-00731 723.

[6] Boyer, F., Lapuerta, C., Minjeaud, S., Piar, B., A multigrid method with local adaptive refinement, application to a ternary Cahn-Hilliard model, ESAIM: Proceedings (2009) hal-00307 186v1.

[7] Demkowicz, L., Computing with *hp*-adaptive finite elements. Vol. 1. One- and two-dimensional elliptic and Maxwell problems, Boca Raton, Chapman & Hall/CRC, 2007.
https://doi.org/10.1201/9781420011685

[8] Demkowicz, L., Devloo, P., Oden, J. T., On an *h*-type mesh-refinement strategy based on minimization of interpolation errors, Computer Methods in Applied Mechanics and Engineering 53 (1985) 67–89. https://doi.org/10.1016/0045-7825(85)90076-3

[9] Demkowicz, L., Rachowicz, W., Devloo, P., Fully automatic *hp*-adaptivity, Journal of Scientific Computing 17 (2002) 117–142. https://doi.org/10.1023/A:1015192312705

[10] Endres, L., Krysl, P., Octasection-based refinement of finite element approximations on tetrahedral meshes that guarantees shape quality, International Journal for Numerical Methods in Engineering 59 (2004) 69–82. https://doi.org/10.1002/nme.863

[11] Fish, J., The *s*-version of the finite element method, Computers and Structures 43 (1992) 539–547.
https://doi.org/10.1016/0045-7949(92)90287-A

[12] Ghosh, S., Manna, S. K., *R*-adapted arbitrary Lagrangian-Eulerian finite-element method in metal-forming simulation, Journal of Materials Engineering and Performance 2 (1993) 271–282.
https://doi.org/10.1007/BF02660296

[13] Grinsprun, E., Krysl, P., Schröder, P., CHARMS: A simple framework for adaptive simulation, ACM Transactions on Graphics (TOG) – Proceedings of ACM SIGGRAPH 21 (2002) 281–290.
https://doi.org/10.1145/566654.566578

[14] Hénap, G., Szabó, L., On numerical solution of elasto-plastic problems by using configurational force driven adaptive methods, Finite Elements in Analysis and Design 92 (2014) 50–59.
https://doi.org/10.1016/j.finel.2014.08.002

[15] Kita, E., Higuchi, K., Kamiya, N., Application of *r*- and *hr*-adaptive BEM to two-dimensional elastic problem, Engineering Analysis with Boundary Elements 24(4) (2000) 317–324.
https://doi.org/10.1016/S0955-7997(00)00011-4

[16] Krysl, P., Grinsprun, E., Schröder, P., Natural hierarchical refinement for finite element methods, International Journal for Numerical Methods in Engineering 56 (2003) 1109–1124.
https://doi.org/10.1002/nme.601

[17] Krysl, P., Trivedi, A., Zhu, B., Object-oriented hierarchical mesh refinement with charms, International Journal for Numerical Methods in Engineering 60 (2004) 1401–1424.
https://doi.org/10.1002/nme.1008

[18] Ladevèze, P., Maunder, E., A general method for recovering equilibrating element traction, Computer Methods in Applied Mechanics and Engineering 137 (1996) 111–151.
https://doi.org/10.1016/S0045-7825(96)01067-5
[19] Nicolas, G., Fouquet, T., Geniaut, S., Cuvilliez, S., Improved adaptive mesh refinement for conformal hexahedral meshes, Advances in Engineering Software 102 (2016) 14–28.
https://doi.org/10.1016/j.advengsoft.2016.07.014

[20] Oden, J. T., Error estimation and control in computational fluid dynamics, The O. C. Zienkiewicz Lecture. Proc. Math. of Finite Elements – MAFELAP VIII, Brunnel University, Uxbridge, 1993, pp. 1–36.

[21] Oden, J. T., Cho, J. R., Adaptive *hpq* finite element methods of hierarchical models for plate- and shell-like structures, Computer Methods in Applied Mechanics and Engineering 136(3–4) (1996) 317–345.
https://doi.org/10.1016/0045-7825(95)00986-8

[22] Oden, J. T., Prudhomme, S., Goal-oriented error estimation and adaptivity for the finite element method, Computers and Mathematics with Applications 41 (2001) 735–756.
https://doi.org/10.1016/S0898-1221(00)00317-5

[23] Prudhomme, S., A posteriori error estimates of quantities of interest, In: Encyclopedia of Applied and Computational Mathematics, Berlin, Heidelberg, Springer Berlin Heidelberg, 2015, pp. 1–5.
https://doi.org/10.1007/978-3-540-70529-1_549

[24] Salagame, R., Belegunda, A., A simple *p*-adaptive refinement procedure for structural shape optimization, Finite Elements in Analysis and Design 24 (1997) 133–155.
https://doi.org/10.1016/S0168-874X(96)00036-4

[25] Zboiński, G., Adaptive *hpq* finite element methods for the analysis of 3D-based models of complex structures. Part 2 – A posteriori error estimation, Computer Methods in Applied Mechanics and Engineering 267 (2013) 531–565.
https://doi.org/10.1016/j.cma.2013.08.018

[26] Zienkiewiecz, O. C., Zhu, J. Z., A simple error estimator and adaptive procedure for practical engineering analysis, Computer Methods in Applied Mechanics and Engineering 24 (1987) 337–357.
https://doi.org/10.1002/nme.1620240206

[27] Zienkiewiecz, O. C., Zhu, J. Z., The superconvergent patch recovery (SPR) and adaptive finite element refinement, Computer Methods in Applied Mechanics and Engineering 101 (1992) 207–224.
https://doi.org/10.1016/0045-7825(92)90023-D
# Appendix

## A. $\mathbb{P}_1$-triangle refinement pattern

The refinement pattern of a $\mathbb{P}_1$-triangle and its equations are specified in Fig. 12.

<figure><img src="image_14.png" /><figcaption>Fig. 12. Refinement pattern of a $\mathbb{P}_1$-triangle element and associated refinement equations</figcaption></figure>

$$ \hat{\varphi}_0^{(0)} = \hat{\varphi}_0^{(1)} + \frac{1}{2}\hat{\varphi}_1^{(1)} + \frac{1}{2}\hat{\varphi}_3^{(1)} $$

$$ \hat{\varphi}_1^{(0)} = \hat{\varphi}_2^{(1)} + \frac{1}{2}\hat{\varphi}_1^{(1)} + \frac{1}{2}\hat{\varphi}_4^{(1)} $$

$$ \hat{\varphi}_2^{(0)} = \hat{\varphi}_5^{(1)} + \frac{1}{2}\hat{\varphi}_3^{(1)} + \frac{1}{2}\hat{\varphi}_5^{(1)} $$

## B. $\mathbb{Q}_1$-quadrangle refinement pattern

The refinement pattern is detailed in Fig. 13 and the corresponding equations are given in (2).

<figure><img src="image_15.png" /><figcaption>Fig. 13. (Left) Square reference element and (Right) associate refinement pattern with basis functions</figcaption></figure>
C. $\mathbb{Q}_1$-cuboid refinement pattern

The refinement pattern and its equations are detailed in Fig. 14.

<figure><img src="image_16.png" /></figure>

<table><thead><tr><th>Refined function</th><th>Coordinates</th><th>Refined function</th><th>Coordinates</th></tr></thead><tbody><tr><td>φ<sub>0</sub><sup>(1)</sup></td><td>(0, 0, 0)</td><td>φ<sub>1</sub><sup>(1)</sup></td><td>(0.5, 0, 0)</td></tr><tr><td>φ<sub>2</sub><sup>(1)</sup></td><td>(1, 0, 0)</td><td>φ<sub>3</sub><sup>(1)</sup></td><td>(0, 0.5, 0)</td></tr><tr><td>φ<sub>4</sub><sup>(1)</sup></td><td>(0.5, 0.5, 0)</td><td>φ<sub>5</sub><sup>(1)</sup></td><td>(1, 0.5, 0)</td></tr><tr><td>φ<sub>6</sub><sup>(1)</sup></td><td>(0, 1, 0)</td><td>φ<sub>7</sub><sup>(1)</sup></td><td>(0.5, 1, 0)</td></tr><tr><td>φ<sub>8</sub><sup>(1)</sup></td><td>(1, 1, 0)</td><td>φ<sub>9</sub><sup>(1)</sup></td><td>(0, 0, 0.5)</td></tr><tr><td>φ<sub>10</sub><sup>(1)</sup></td><td>(0.5, 0, 0.5)</td><td>φ<sub>11</sub><sup>(1)</sup></td><td>(1, 0, 0.5)</td></tr><tr><td>φ<sub>12</sub><sup>(1)</sup></td><td>(0, 0.5, 0.5)</td><td>φ<sub>13</sub><sup>(1)</sup></td><td>(0.5, 0.5, 0.5)</td></tr><tr><td>φ<sub>14</sub><sup>(1)</sup></td><td>(1, 0.5, 0.5)</td><td>φ<sub>15</sub><sup>(1)</sup></td><td>(0, 1, 0.5)</td></tr><tr><td>φ<sub>16</sub><sup>(1)</sup></td><td>(0.5, 1, 0.5)</td><td>φ<sub>17</sub><sup>(1)</sup></td><td>(1, 1, 0.5)</td></tr><tr><td>φ<sub>18</sub><sup>(1)</sup></td><td>(0, 0, 1)</td><td>φ<sub>19</sub><sup>(1)</sup></td><td>(0.5, 0, 1)</td></tr><tr><td>φ<sub>20</sub><sup>(1)</sup></td><td>(1, 0, 1)</td><td>φ<sub>21</sub><sup>(1)</sup></td><td>(0, 0.5, 1)</td></tr><tr><td>φ<sub>22</sub><sup>(1)</sup></td><td>(0.5, 0.5, 1)</td><td>φ<sub>23</sub><sup>(1)</sup></td><td>(1, 0.5, 1)</td></tr><tr><td>φ<sub>24</sub><sup>(1)</sup></td><td>(0, 1, 1)</td><td>φ<sub>25</sub><sup>(1)</sup></td><td>(0.5, 1, 1)</td></tr><tr><td>φ<sub>26</sub><sup>(1)</sup></td><td>(1, 1, 1)</td><td></td><td></td></tr></tbody></table>

$$ \hat{\varphi}_{0}^{(0)} = \hat{\varphi}_{0}^{(1)} + \frac{1}{2} (\hat{\varphi}_{1}^{(1)} + \hat{\varphi}_{3}^{(1)} + \hat{\varphi}_{9}^{(1)}) + \frac{1}{4} (\hat{\varphi}_{4}^{(1)} + \hat{\varphi}_{10}^{(1)} + \hat{\varphi}_{12}^{(1)}) + \frac{1}{8} \hat{\varphi}_{13}^{(1)} $$

$$ \hat{\varphi}_{1}^{(0)} = \hat{\varphi}_{2}^{(1)} + \frac{1}{2} (\hat{\varphi}_{1}^{(1)} + \hat{\varphi}_{11}^{(1)} + \hat{\varphi}_{5}^{(1)}) + \frac{1}{4} (\hat{\varphi}_{4}^{(1)} + \hat{\varphi}_{10}^{(1)} + \hat{\varphi}_{14}^{(1)}) + \frac{1}{8} \hat{\varphi}_{13}^{(1)} $$

$$ \hat{\varphi}_{2}^{(0)} = \hat{\varphi}_{8}^{(1)} + \frac{1}{2} (\hat{\varphi}_{7}^{(1)} + \hat{\varphi}_{5}^{(1)} + \hat{\varphi}_{17}^{(1)}) + \frac{1}{4} (\hat{\varphi}_{16}^{(1)} + \hat{\varphi}_{4}^{(1)} + \hat{\varphi}_{14}^{(1)}) + \frac{1}{8} \hat{\varphi}_{13}^{(1)} $$

$$ \hat{\varphi}_{3}^{(0)} = \hat{\varphi}_{6}^{(1)} + \frac{1}{2} (\hat{\varphi}_{7}^{(1)} + \hat{\varphi}_{3}^{(1)} + \hat{\varphi}_{15}^{(1)}) + \frac{1}{4} (\hat{\varphi}_{16}^{(1)} + \hat{\varphi}_{4}^{(1)} + \hat{\varphi}_{12}^{(1)}) + \frac{1}{8} \hat{\varphi}_{13}^{(1)} $$

$$ \hat{\varphi}_{4}^{(0)} = \hat{\varphi}_{18}^{(1)} + \frac{1}{2} (\hat{\varphi}_{19}^{(1)} + \hat{\varphi}_{21}^{(1)} + \hat{\varphi}_{9}^{(1)}) + \frac{1}{4} (\hat{\varphi}_{22}^{(1)} + \hat{\varphi}_{10}^{(1)} + \hat{\varphi}_{12}^{(1)}) + \frac{1}{8} \hat{\varphi}_{13}^{(1)} $$

$$ \hat{\varphi}_{5}^{(0)} = \hat{\varphi}_{20}^{(1)} + \frac{1}{2} (\hat{\varphi}_{19}^{(1)} + \hat{\varphi}_{23}^{(1)} + \hat{\varphi}_{11}^{(1)}) + \frac{1}{4} (\hat{\varphi}_{22}^{(1)} + \hat{\varphi}_{14}^{(1)} + \hat{\varphi}_{10}^{(1)}) + \frac{1}{8} \hat{\varphi}_{13}^{(1)} $$

$$ \hat{\varphi}_{6}^{(0)} = \hat{\varphi}_{26}^{(1)} + \frac{1}{2} (\hat{\varphi}_{25}^{(1)} + \hat{\varphi}_{23}^{(1)} + \hat{\varphi}_{17}^{(1)}) + \frac{1}{4} (\hat{\varphi}_{14}^{(1)} + \hat{\varphi}_{22}^{(1)} + \hat{\varphi}_{16}^{(1)}) + \frac{1}{8} \hat{\varphi}_{13}^{(1)} $$

$$ \hat{\varphi}_{7}^{(0)} = \hat{\varphi}_{24}^{(1)} + \frac{1}{2} (\hat{\varphi}_{21}^{(1)} + \hat{\varphi}_{25}^{(1)} + \hat{\varphi}_{15}^{(1)}) + \frac{1}{4} (\hat{\varphi}_{12}^{(1)} + \hat{\varphi}_{16}^{(1)} + \hat{\varphi}_{22}^{(1)}) + \frac{1}{8} \hat{\varphi}_{13}^{(1)} $$

<caption>Fig. 14. $\mathbb{Q}_1$-cuboid refinement pattern, the refined basis functions coordinates and refinement equations</caption>