# Error Indicator to Assess the Quality of Finite Element Frame Analysis

Shahrokh Ghavamian, Gilles Pijaudier-Cabot, Jacky Mazars

► To cite this version:

Shahrokh Ghavamian, Gilles Pijaudier-Cabot, Jacky Mazars. Error Indicator to Assess the Quality of Finite Element Frame Analysis. Journal of Structural Engineering, ASCE, 2002, 128 (1), pp.129-132. 10.1061/(ASCE)0733-9445(2002)128:1(129) . hal-01006805

HAL Id: hal-01006805

https://hal.science/hal-01006805

Submitted on 28 Jan 2017

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
# Error Indicator to Assess Quality of Finite-Element Frame Analyses

Shahrokh Ghavamian<fnref n="1" />, Gilles Pijaudier-Cabot<fnref n="2" /> and Jacky Mazars<fnref n="1" />

**Abstract:** The purpose of this article is to describe a tool developed in order to equip a finite-element software based on nonlinear beam analysis with an error indicator aimed at measuring discretization errors. The technique is an extension of the error estimate devised by Ladeve`ze et al. in 1991 in which the finite-element solution is compared to a statically admissible distribution of the generalized stresses. This last distribution is built element by element which is fast and easy to implement. A nonlinear analysis of a single bay frame in static's is presented as an example.

**keywords:** Finite elements; Frames; Analysis; Computer software; Errors.

## Introduction

Advances achieved in construction engineering and more restrictive safety requirements, especially with respect to severe and exceptional loads such as earthquakes, have induced major needs in understanding and predicting the response of civil engineering structures. In most instances, the validity of computational results and their efficiency are of great concern. Quality is, however, quite intricate to evaluate in a broad sense because several factors are important and interact. A finite-element model is based on a series of assumptions, with respect to its geometry, to the boundary conditions, to the representation of the applied loads, to the construction process, and also to the constitutive relations of the different materials. The purpose of this technical note is to present a simple tool for evaluating one of these components: the quality of the discretisation of the finite-element (FE) model.

Many techniques can be found in the literature, related to adaptive meshing and mesh optimization. A posteriori error indicators are mostly elaborated from four different approaches. The first one relies on the analysis of the stress field in the discretised structure (see, e.g., Zienckiewicz and Zhu 1987). Two sets of results are compared: one being discrete, the other being smoothed. The method devised by Babuska and Rheinholdt (1982) is based on the measure of the residual forces in the equations of equilibrium which are never totally balanced in nonlinear FE analyses. Huerta and Diez (2000) use a comparison between two finite-element models. A very fine discretization is considered to be the reference (i.e., quasi-exact) solution, against which the user's mesh is compared. Finally, Ladeve`ze et al. (1991) define an

error indicator as the distance between the finite-element solution and a statically admissible solution.

Most of these indicators provide an upper bound of the true error in elastic analyses. It is not rigorously demonstrated that they provide an upper bound of the error in general nonlinear analyses but they are convergent, i.e., equal to zero when the exact solution is found.

In this short technical note, we show that the technique due to Ladeve`ze et al. can be easily implemented in finite-element analyses of frames.

## Finite-Element Model

The finite-element model used here is a nonlinear static and dynamic computer program based on a layered finite-element description of reinforced concrete beams and frames (Bazant et al. 1987). This program uses Euler-Bernoulli beam elements. The two-dimensional beam elements possess two nodes. The horizontal displacement $u_x$ (parallel to the neutral axis of the beam) is interpolated with a linear function over the element. The vertical displacement $u_y$ (perpendicular to the neutral axis of the beam) is interpolated with a classical third order hermitian function. These displacements are computed at the centroid of the cross section:

$$ \begin{bmatrix} u_x(x) \\ u_y(x) \end{bmatrix} = N \begin{bmatrix} u_i \\ u_j \\ v_i \\ v_j \\ \theta_i \\ \theta_j \end{bmatrix} \qquad (1) $$

where N denotes the matrix of shape functions, expressed in term of the element neutral axis coordinate x. $u_i$, $v_i$, and $\theta_i$ are the horizontal displacement, the vertical displacement and the rotation at node i, respectively. According to the Euler-Bernoulli theory, the horizontal displacement at any point of the cross section of coordinates (x,y) is $u(x,y)=u_x(x)-y\partial_y u_y(x)/\partial x$ and the nonzero component of the strain (in direction x) is $\varepsilon(x,y) = \partial u(x,y)/\partial x$.

The beam cross section is decomposed into layers made of material with nonlinear uniaxial stress-strain laws. The constitutive model used for concrete is a rate independent damage model. A complete description of the constitutive equations can be found

<footnote n="1">¹LMT-Cachan, ENS Cachan, 61 avenue du Président Wilson F-94235 Cachan cedex, France</footnote>

<footnote n="2">²R&DO, Laboratoire de Génie Civil de Nantes—Saint Nazaire, Ecole Centrale de Nantes, 1 rue de la Noe, BP 92101, F-44321 Nantes cedex 3, France</footnote>
<figure><img src="image_2.png" /><figcaption>Fig. 1. Response of concrete constitutive law to tension-compression loading</figcaption></figure>

in Dubé et al. (1996). Fig. 1 shows the concrete response for tension-compression cycles. For the reinforcements, a one-dimensional elastoplastic model with linear hardening is implemented. It should be pointed out that any other constitutive relations could be implemented as well, without any modification of the calculation of the error indicator described in the next section.

## Error Indicator

The finite-element solution provides the values of the generalized stresses, i.e., the normal and shear forces ($N^k, V^k$) and the bending moment $M^k$ at each nodal point. The upper-script $k$ indicates that these expressions are obtained from a displacement-based (finite element) approximation.

The error estimation consists in measuring the difference between two generalized stress fields. The first one is the finite-element (FE) solution. It is kinematically admissible, which means that it satisfies the displacement boundary conditions, compatibility conditions, but that the generalized stresses field derived with the help of the constitutive relations do not satisfy the equilibrium equations pointwise. The second one is statically admissible. It satisfies the equilibrium equations pointwise, and the force boundary conditions as well. With the help of the constitutive relations, strains can be derived from this field but they are, in general, not compatible with a displacement field and with the displacement boundary conditions. Theoretically, if those two solutions are identical, all the governing equations of the problem are satisfied and the exact solution has been found. The distance between the kinematically (FE) and statically admissible fields is an indicator of the error due to mesh discretisation. It does not mean that the FE solution is better than the statically admissible solution or conversely, any solution can be used for design purposes, but both are approximates.

## Statically admissible generalized stress field

The statically admissible solution (generalized stress) denoted ($N^s, V^s, M^s$) should verify exactly the differential equations of equilibrium/motion:

$$ \frac{dN^s(x,t)}{dx} + f_x = \rho \frac{d^2[u_x(x,t)]}{dt^2} \quad (2a) $$

$$ \frac{dV^s(x,t)}{dx} + f_y = \rho \frac{d^2[u_y(x,t)]}{dt^2} \quad (2b) $$

$$ \frac{dM^s(x,t)}{dx} + V^s(x,t) = 0 \quad (2c) $$

where ($f_x, f_y$) are, for most applications, constant distributed forces, $\rho$=mass per unit length of beam, and rotational inertia terms have been neglected.

For the sake of simplicity, let us omit the inertia terms and focus on statics. Equation (2) may be integrated over each beam element (recall that distributed forces are constant):

$$ N^s(x) = a + f_x x \quad (3a) $$

$$ M^s(x) = \frac{(x)^2}{2} + bx + c \quad (3b) $$

Constants ($a, b, c$) need to be computed. Since each finite element is in equilibrium, according to the FE solution, nodal values in Eq. (3) are equated with the FE solution. Constant $a$ is computed from Eq. (3a) by expressing this relation at either nodes of the element with the same result. Constants ($b, c$) are computed from Eq. (3b) expressed at each node:

$$ \begin{aligned} N^s(0) &= a = N^k(0) \\ M^s(0) &= c = M^k(0) \end{aligned} \quad (4) $$

$$ M^s(1) = \frac{L^2}{2} + bL + c = M^k(1) $$

where $L$=element length. Within each beam element, Eq. (4) provides a statically admissible solution. It is entirely local, i.e., expressed at the element level. Hence its resolution is a post processing operation which is performed element per element successively. The same method applies to more sophisticated load distributions without any modifications.

It is important also to satisfy equilibrium across each beam element. As opposed to the case of two-dimensional finite elements where stress vectors are not continuous across each element, this is automatically satisfied in the beam formulation (because of the order of continuity of the interpolation functions). This is the major reason why the statically admissible solution ($N^s, V^s, M^s$) is so easy to construct. Note that another technique due to Carol and Murcia (1989) could have been used for obtaining statically admissible fields. It would have been more computer time consuming because a global set of equations of size equal to the number of degree of freedom of the problem needs to be solved.

In dynamics, inertial forces must be taken into account in the derivation of the «statically» admissible field of generalized forces. With a classical Newmark integration scheme, which assumes a constant acceleration over each time step, these inertia forces can be viewed as additional distributed loads which are constant during each time step and follow the interpolation of the displacements in space. Hence the «statically» admissible field of normal force becomes a polynomial of order 2 and the «statically» admissible distribution of bending moment becomes a polynomial of order 5. The derivation of the statically admissible solution is exactly the same as in static's, except that the external loads are slightly more complex.

## Definition of the error indicator

We have now two solutions, which can be compared. The distance $E$ between the two solutions is defined as an $L_2$ norm according to Ladevèze et al. (1991). $E$ is further denoted as the global error indicator. It is an approximation of the true error
<figure><img src="image_3.png" /><figcaption>Fig. 2. Single bay case study</figcaption></figure>

(which cannot be known without the knowledge of the exact solution). It can be also decomposed according to the discretization:

$$E = \sum_{\text{elements}} e_p$$

with

$$e_p = \frac{\int_{\text{elements}} P \sqrt{(N^s(x) - N^k(x))^2} dx}{\int_{\text{structure}} \sqrt{(N^k(x) + N^s(x))^2} dx} + \frac{\int_{\text{element p}} \sqrt{(M^s(x) - M^k(x))^2} dx}{\int_{\text{structure}} \sqrt{(M^k(x) + M^s(x))^2} dx} \quad (5)$$

where $e_p$ is the local error indicator in element $p$ (distance between the two solutions in the element). Numerically, the local error indicator is computed over the element according to the quadrature of the finite-element solution. We may retain, in each element, the maximum value of the indicator over the entire calculation. Consequently, the local and global error indicators are either constant or monotonically increasing in the course of the calculation. The above definitions can be used in statics and in dynamics.

## Application

As an illustrative example, consider the single bay plane structure shown in Fig. 2. It is made of reinforced concrete. The reinforcement is constant over the frame; it is incorporated in the beam element assuming a perfect bond between steel and concrete. The loading is displacement controlled: a horizontal displacement is applied to the left corner node. The maximum displacement is 0.5 m, decomposed into a loading history of 100 increments of equal size. Four different meshes shown in Fig. 3 are considered. The number of finite elements ranges from 3 to 25. The distributions of tensile damage, of the plastic strain in the reinforcement, and of the local error indicator over the structure at the end of the loading process are also shown in this figure.

Fig. 4 shows the load-deflection curves (horizontal force versus horizontal displacement at the top left corner) for each finite-

<figure><img src="image_4.png" /><figcaption>Fig. 3. Local results (tensile damage, rebar yielding, and local error distributions) for different meshes</figcaption></figure>
<figure><img src="image_5.png" /><figcaption>**Fig. 4.** Force displacement behavior for single bay example</figcaption></figure>

element mesh. It can be observed that they are quite similar for the two finest meshes. In fact the FE solution seems to have converged toward the exact one, in terms of mesh refinement. In terms of local quantities (damage, plastic strain, and local errors), there are still differences as shown in Fig. 3. The minimum and maximum local errors keep on diminishing upon refinement. The maximum local error, however, concentrates in the joints where material nonlinearity is more intense, as expected. According to this analysis, the global response can be obtained quite accurately with a moderately dense mesh, local quantities such as damage or plastic strain are much more demanding from the viewpoint of the finite-element discretization.

## Conclusion

An extension of the error indicator devised by Ladevèze et al. (1991) to the simple case of beam analysis has been presented.

The method relies on the comparison between the finite-element solution and a statically admissible field of generalized stresses. It is quite simple to implement because statically admissible fields of generalized stresses are constructed element per element in a post-processing course of the FE calculation. It can be easily extended to a three-dimensional configuration, or to cases where distributed forces are arbitrary (known) functions, such as in dynamics.

Two forms of error estimate can be considered: the global error which is an overall measure of the quality of the finite-element discretization and the distribution of the local errors which pin-points locations where the discretization is too coarse or could be coarser. This latter information could serve as a point of entry in the design of optimized meshes.

## References

Babuska, I., and Rheinholdt, W. C. (1982). "Computational error estimates and adaptive processes for some nonlinear structural problems." *Comput. Methods Appl. Mech. Eng.*, 34, 895-937.

Bazant, Z. P., Pan, J., and Pijaudier-Cabot, G. (1987). "Softening in reinforced concrete beams and frames." *J. Struct. Eng.* 113(12), 2333-2347.

Carol, I., and Murcia, J. (1989). "Nonlinear time-dependent analysis of planar frames using an "exact" formulation—I. Theory." *Comput. Struct.*, 33, 79-87.

Dubé, J. F., Pijaudier-Cabot, G., and La Borderie, C. (1996). "Rate dependent damage model for concrete in dynamics." *J. Eng. Mech.*, 122(10), pp. 939-947.

Huerta, A., and Diez, P. (2000). "Error estimation including pollution assessment for nonlinear finite element analysis." *Comput. Methods Appl. Mech. Eng.*, 181, 21-411.

Ladevèze, P., Pelle, J. P., and Rougeot, P. (1991). "Error estimation and mesh optimization for classical finite elements." *Eng. Comput.*, 8, 69-80.

Zienkiewicz, O. C., and Zhu, J. Z. (1987). "A simple error estimator and adaptive procedure for practical engineering analysis." *Int. J. Numer. Methods Eng.*, 24, 337-357.