2. Principle of the method#
2.1. Equations to be solved#
We consider solution \((u,s)\) of a linear elastic problem verifying:
the equilibrium equations:
\(\mathrm{\{}\begin{array}{c}\text{Lu}\mathrm{=}q\text{dans}\Omega \\ {\sigma }_{\text{ij}}{n}_{j}\mathrm{=}{\stackrel{ˉ}{t}}_{i}\text{sur}{\Gamma }_{t}\end{array}\)
with \(L={}^{t}\mathrm{BDB}\) elasticity operator
compatibility equations:
\(\mathrm{\{}\begin{array}{c}\varepsilon \mathrm{=}\text{Bu}\\ u\mathrm{=}\stackrel{ˉ}{u}\text{sur}{\Gamma }_{u}\end{array}\)
with \(\Gamma ={\Gamma }_{u}\cup {\Gamma }_{t}\)
the law of behavior:
\(\sigma =D\varepsilon\)
The problem discretized by finite elements consists in finding the \(({u}_{h},{\sigma }_{h})\) solution of:
\({u}_{h}=N{\stackrel{ˉ}{u}}_{h}\) eq 2.1-1
checking \(K{\stackrel{ˉ}{u}}_{h}=\text{f}\)
with \(K=\underset{\Omega }{\int }{}^{t}\text{}(\mathrm{BN})D(\mathrm{BN})d\Omega\)
\(\text{f}=\underset{\Omega }{\int }{}^{t}\text{}\mathrm{Nq}d\Omega +\underset{{G}_{t}}{\int }{}^{t}\text{}N\stackrel{ˉ}{t}\mathrm{dG}\)
where:
\({\stackrel{ˉ}{\text{u}}}_{h}\) represents the nodal unknowns of displacement
\(N\) the associated form functions
The constraints are calculated from the displacements by the relationship:
\({\sigma }_{h}={\mathrm{DBu}}_{h}\) eq 2.1-2
2.2. Error estimator and effectiveness index#
We note |
\(e=u-{u}_{h}\) |
the error on the trips |
\({e}_{s}=s-{s}_{h}\) |
the constraint error |
The energy norm for error \(e\) is written as:
\(\parallel \text{e}\parallel ={(\underset{\Omega }{\int }{}^{t}\text{}\mathrm{eLe}d\Omega )}^{1/2}\)
in the case of elasticity
\(={(\underset{\Omega }{\int }{}^{t}\text{}{e}_{\sigma }{D}^{\text{-}1}{e}_{\sigma }d\Omega )}^{1/2}\) eq 2.2-1
The overall error above is broken down into the following sum of elementary errors:
\({\parallel e\parallel }^{2}=\sum _{i=1}^{N}{\parallel e\parallel }_{i}^{2}\)
where |
\(N\) is the total number of items. |
\({\parallel e\parallel }_{i}\) represents the local error indicator on the \(i\) element. |
The aim is to estimate the exact error using equation [éq 2.2-1] formulated in constraints. The basic idea of the method is to build a new constraint field noted \({\sigma }^{\text{*}}\) from \({\sigma }_{h}\) and such as:
\({e}_{\sigma }\approx {e}_{\sigma }^{\text{*}}={\sigma }^{\text{*}}-{\sigma }_{h}\)
The error estimator will then be written as:
\({}^{0}\text{}\parallel e\parallel ={(\underset{\Omega }{\int }{}^{t}\text{}{e}_{\sigma }^{\text{*}}{D}^{\text{-}1}{e}_{\sigma }^{\text{*}}d\Omega )}^{1/2}\)
The quality of the estimator is measured by the quantity \(\mathrm{\theta }\), called the estimator effectiveness index:
\(\mathrm{\theta }=\frac{{}^{0}\text{}\parallel e\parallel }{\parallel e\parallel }\)
An error estimator is said to be asymptotically accurate if \(\mathrm{\theta }\to 1\) when \(\parallel e\parallel \to 0\) (or when \(h\to 0\)), which means that the estimated error will always converge to the exact error when the error decreases.
Obviously, the reliability of \({}^{0}\text{}\parallel e\parallel\) depends on the « quality » of \({\sigma }^{\text{*}}\).
The two versions of the ZHU - ZIENKIEWICZ estimator differ at this level (see [§3]).
2.3. Construction of an asymptotically accurate estimator#
The characterization of such an estimator is provided by the following theorem (see [bib 2]).
Theorem
If \(\parallel {e}^{\text{*}}\parallel =\parallel {u-u}^{\text{*}}\parallel\) is the error norm of the reconstructed solution, then the error estimator \({}^{0}\text{}\parallel e\parallel\) defined earlier is asymptotically accurate
If \(\frac{\parallel {e}^{\text{*}}\parallel }{\parallel e\parallel }\to 0\) when \(\parallel e\parallel \to 0\)
This condition is met if the convergence rate with \(h\) of \(\parallel {e}^{\text{*}}\parallel\) is greater than that of \(\parallel e\parallel\). Typically, if we assume that the exact error of the finite element approximation converges to \(\parallel e\parallel =0({h}^{p})\)
and the error of the reconstructed solution in
\(\parallel {e}^{\text{*}}\parallel =0({h}^{p\text{+}\alpha })\) with \(\alpha >0\)
then a simple calculation gives:
\(1-0({h}^{\alpha })\le \mathrm{\theta }\le 1+0({h}^{\alpha })\)
and so \(\mathrm{\theta }\to 1\) when \(h\to 0\)