Zdzisław Wi˛eckowski Department of Mechanics of Materials Technical University of Łód´z, Łód´z, Poland
e–mail: zwi@ck-sg.p.lodz.pl
ABSTRACT
Summary A technique of stress recovery for a finite element solution is investigated. A smoothed, stati-cally admissible stress field is computed by solving the local problem of minimization of complementary energy on patches of elements. The energy error of the approximate solution is estimated by the use of the recovered stress field. The method is illustrated by the example of the Poisson equation; numerical results are related to the problem of torsion of a prismatic bar.
The superconvergent patch recovery (SPR) method is a very popular technique of finding a stress field which is smoother and more accurate than the stress field related directly to the displacement field obtained by the finite element analysis. In SPR method, the improved stress field is usually considered in the form of polynomial (the degree of which is higher than the polynomial repre-senting stresses in the original FE solution) which fits stress values at superconvergent points of elements belonging to a patch surrounding a node of the element mesh. Coefficients of the poly-nomial are found by the use of the least square method [1]. SPR approach can be enriched by a requirement that the recovered stress field should satisfy the equilibrium equations in an approxi-mate (e.g. [2]) or exact (e.g. [3]) way. There are other recovery approaches where the equilibrium equations are satisfied on the element patch by recovered stresses. The equilibrium equations can be satisfied exactly by constructing the self-equilibrated element boundary tractions equivalent to those obtained from the kinematically admissible finite element solution and then finding the statically admissible stress field, related to the tractions, inside the element domain [4]. The equi-librium equations can also be satisfied in the weak sense like in the displacement-based finite element analysis—this procedure is known as the recovery by equilibration of patches (REP) [5].
In all the approaches mentioned above, the stress field obtained by the finite element analysis is used as input data in the recovery procedure. In contrast to these methods, the technique described in the present paper does not utilize stresses as input data but displacements which are of higher accuracy than stresses. The recovered stress field is found by minimizing the functional of com-plementary energy on the set of statically admissible stress fields where the displacements on the patch boundary are taken from the kinematically admissible finite element solution.
Let us consider two-dimensional physical problems described by Poisson’s equation, u; =f,
=1;2, where stresses (fluxes) satisfy the equilibrium (balance) equation of the form:q ; =f withq
= u
; . Letp denote the region occupied by an element patch where the stresses are to be recovered. For any patch p the following local minimization problem is solved: To find
q s
p
2Y such thatqsp minimizes the functional of complementary energy
(q)=
on the following set of statically admissible stress (flux) fields:
whereuhdenotes the approximation of the displacement (potential) field andn is the unit vector outwardly normal to the patch boundary. The statically admissible fields of stresses are constructed by the use of the Prandtl stress function, :
q
whereq is a particular solution of equation q ; = f. In the case of triangular or quadrilateral elements with linear or bilinear interpolation functions, respectively, function is approximated by the following polynomial:
=
The above representation of function means that both the components of vector q are linear functions of coordinatesx =x1 andy =x2while the finite element approximations gives these components as piecewise constant functions. After substituting equations (1) and (2) to the expres-sion for functional and calculating its first variation, we obtain the linear system of algebraic equations with the symmetric matrix. After calculation of the vector of coefficientsa1
;a
2
;:::;a
5, the components of stressesqps at the “central” node of the patch can be calculated and used as the nodal values so that the final recovered field of stresses,q, can be written in the form
q
whereNare the same interpolation functions as used in the approximation of the displacement field,qxandqy denote vectors of nodal values of recovered stress components.
When the improved field of stresses is known, the error of the approximate solution can be esti-mated. The estimated error measured by means of energy norm can be expressed as follows:
ke k=
As an example, the problem of torsion of the prismatic bar with triangular cross-section having equal sides is considered. All the above derivations are valid if we putq1
=
31,q2
=
32,=2G where G denotes the shear modulus, and uh, which appears in expression for complementary energy, is the warp of the cross-section of the bar. Calculation have been made assuming that
G=1, the length of the cross-section edge is equal to 1, and the torsion angle for the unit length of the bar is equal to 1. Four uniform meshes generated automatically with element diameters
h=0:25,0:1,0:04and0:016are utilized in the computations.
The relative error of the approximate solution, kek=kqk, has been estimated and compared with its exact value (the exact solution is known for the considered problem) for all the element meshes used in the calculations. The comparison is shown in Fig. 1. In the same figure, the diagram for
0 0.05 0.1 0.15 0.2 0.25 0
10 20 30 40
Element size
Relative error [%]
Estimated error Exact error
0 0.05 0.1 0.15 0.2 0.25
0.92 0.94 0.96 0.98 1
Element size
Effectivity index
Figure 1: Energy error and effectivity index
the effectivity index, ke k=ke k, which measures the accuracy of the error estimator, is also shown. In the last definition, ke kdenotes the exact value of the error. The recovered stress field obtained in the case of rather coarse mesh forh =0:1is shown in Fig. 2 and compared with the exact solution.
As can be seen from the presented results, the recovery procedure described in the paper gives the smoothed stress field very close to the exact one and is characterized by the effectivity index for the error estimator which slightly differs from unity (0.920 .. 0.996 depending on mesh parameterh).
REFERENCES
[1] O.C. Zienkiewicz and J.Z. Zhu. The superconvergent patch recovery and a posteriori error estimates. Part 1: The recovery technique. Int. J. Num. Meth. Eng., 33, 1331–1364, (1992).
[2] N.-E. Wiberg, F. Abdulwahab and S. Ziukas. Enhanced superconvergent patch recovery in-corporating equilibrium and boundary conditions. Int. J. Num. Meth. Eng., 37, 3417–3440, (1994).
[3] T. Kvamsdal, K.M. Okstad and K.M. Mathisen. Error estimation based on statically admissible stress fields. In J.-A. Désidéri et al., editors, Numerical Methods in Engineering ‘96, pages 14–20, Chichester. John Wiley & Sons, (1996).
[4] P. Ladeveze, G. Coffignal and J.P. Pelle. Accuracy of elastoplastic and dynamic analysis. In I. Babuška at al., editors, Accuracy Estimates and Adaptive Refinements in Finite Element Computations, pages 181–203, Chichester. John Wiley & Sons, (1986).
[5] B. Boroomand and O.C. Zienkiewicz. Recovery by equilibrium patches (REP). Int. J. Num.
Meth. Eng., 40, 137–154, (1997).
4
Figure 2: Recovered (left) and exact (right) stresses:31,32and
p