You do not have permission to edit this page, for the following reason:

You are not allowed to execute the action you have requested.


You can view and copy the source of this page.

x
 
1
<!-- metadata commented in wiki content
2
==ERROR ESTIMATION AND MESH ADAPTIVITY IN INCOMPRESSIBLE VISCOUS FLOWS USING  A RESIDUAL POWER APPROACH ==
3
4
'''Eugenio Oñate, Joaquin Arteaga, Julio García and Roberto Flores
5
6
~&nbsp;
7
8
International Center for Numerical Methods in Engineering (CIMNE)
9
10
Universidad Politécnica de Cataluña
11
12
Campus Norte UPC, 08034 Barcelona, Spain
13
14
E-mail: cimne@cimne.upc.es, web page: http://www.cimne.upc.es/'''
15
-->
16
17
==Abstract==
18
19
A methodology for error estimation and mesh adaptation for finite element (FE) analysis of incompressible viscous flow is presented. The error estimation method is based on the evaluation of the energy rate (the power) of the FE residuals of the momentum and incompressibility equations. The residuals are computed using recovered values of the derivatives of the velocity and pressure variable obtained via a nodal derivative recovery technique. Two mesh adaptation procedures based on: a)  the equi-distribution of the residual power  among the elements in the mesh and b) the equi-distribution of the density of the total residual power are presented. The stabilized form of the Navier-Stokes equations using a finite calculus (FIC) formulation is solved with  a fractional step FE scheme. This allows the use of linear triangles and tetrahedra with an equal order interpolation for the velocity and pressure variables. A nodal-based approach is used for computing  the residual power integrals. The examples show the ability of the error estimation  and the mesh adaptation process to capture the high gradients of the solution in the vecinity of boundary layers and in zones of high vorticity of the fluid for both steady state and transient flow situations.
20
21
==1 INTRODUCTION==
22
23
The estimation of the error in the finite element solution of the Navier-Stokes equations for an incompressible fluid is a difficult problem <span id='citeF-1'></span>[[#cite-1|1]]. The presence of convective terms in the momentum equations prevents the use of standard error estimators based on the energy norm of the error developed for elliptic type problems typical of solid mechanics and heat transfer analysis <span id='citeF-2'></span>[[#cite-2|2]]. Extension of these methods to fluid mechanics problems have been reported in <span id='citeF-5'></span>[[#cite-5|5]]. A different class of mesh refinement procedures for fluid flow problems based on error indicators using approximations of the Hessian operator has been reported in <span id='citeF-6'></span>[[#cite-6|6]]. Alternative error estimation approaches for fluid flow problems  can be found in <span id='citeF-11'></span>[[#cite-11|11]].
24
25
The error estimation method presented in this paper is based on the evaluation of a particular form of the energy rate (the power) of the finite element residuals of the momentum and mass balance (incompressibility) equations. The residuals are computed using recovered values of the derivatives of the velocity and pressure variable obtained via a nodal derivative recovery technique. Two mesh adaptation processes based on: a)  the equi-distribution of the residual power  among the elements in the mesh and b) the equi-distribution of the density of the total residual power are presented. The stabilized form of the Navier-Stokes equations using a finite calculus (FIC) procedure developed by the authors <span id='citeF-17'></span>[[#cite-17|17]] is used for the FE solution of the examples presented in the paper using a fractional step scheme. This allows the use of linear triangles and tetrahedra with an equal order interpolation for the velocity and pressure variables. A nodal-based approach is used for computation of the residual power error and the total power generated by the viscous stresses and the pressure terms in the mesh. The residual power form derived from the FIC approach includes power terms emanating from the velocity field and from the gradient of the velocity field. It is shown that accounting for the terms involving the velocity gradient provides an effective  error estimation and mesh adaptation procedure for CFD problems. Also the FIC form of the residual power expression can be effectively extended for error estimation and mesh adaptivity of other problems in computational mechanics <span id='citeF-22'></span>[[#cite-22|22]].
26
27
The examples show the ability of the error estimation  and the mesh adaptation process to capture the high gradients of the solution in the vecinity of boundary layers and in zones of high vorticity of the fluid for both steady state and transient flow situations.
28
29
==2 STABILIZED FINITE ELEMENT SOLUTION OF THE INCOMPRESSIBLE NAVIER-STOKES EQUATIONS==
30
31
We consider here the solution of the Navier-Stokes equations for an incompressible fluid written as
32
33
'''Momentum'''
34
35
{| class="formulaSCP" style="width: 100%; text-align: left;" 
36
|-
37
| 
38
{| style="text-align: left; margin:auto;width: 100%;" 
39
|-
40
| style="text-align: center;" | <math>r_{m_i}=0 \qquad \hbox{in }\Omega  </math>
41
|}
42
| style="width: 5px;text-align: right;white-space: nowrap;" | (1)
43
|}
44
45
'''Mass balance'''
46
47
{| class="formulaSCP" style="width: 100%; text-align: left;" 
48
|-
49
| 
50
{| style="text-align: left; margin:auto;width: 100%;" 
51
|-
52
| style="text-align: center;" | <math>r_v=0 \qquad \hbox{in }\Omega  </math>
53
|}
54
| style="width: 5px;text-align: right;white-space: nowrap;" | (2)
55
|}
56
57
For the transient case
58
59
{| class="formulaSCP" style="width: 100%; text-align: left;" 
60
|-
61
| 
62
{| style="text-align: left; margin:auto;width: 100%;" 
63
|-
64
| style="text-align: center;" | <math>r_{m_i}=\rho \left[{\partial u_i \over \partial t} + u_j {\partial u_i \over \partial x_j}\right]+{\partial p\over \partial x_i}-{\partial \tau _{ij}\over \partial x_j}-b_i </math>
65
|}
66
| style="width: 5px;text-align: right;white-space: nowrap;" | (3)
67
|}
68
69
{| class="formulaSCP" style="width: 100%; text-align: left;" 
70
|-
71
| 
72
{| style="text-align: left; margin:auto;width: 100%;" 
73
|-
74
| style="text-align: center;" | <math>r_v = {\partial u_i\over \partial x_i} </math>
75
|}
76
| style="width: 5px;text-align: right;white-space: nowrap;" | (4)
77
|}
78
79
with <math display="inline">i,j=1,n_d</math> when <math display="inline">n_d</math> is the number of space dimensions (<math display="inline">n_d =2</math> for 2D flow problems).
80
81
In eq.(3) <math display="inline">\rho </math> is the fluid density (here assumed to be constant), <math display="inline">u_i</math> is the velocity component in the <math display="inline">i</math>-th direction, <math display="inline">p</math> the pressure, <math display="inline">b_i</math> the body forces and <math display="inline">\tau _{ij}</math> the viscous stress tensor components related to the velocity gradients through the fluid viscosity <math display="inline">\mu </math> by
82
83
{| class="formulaSCP" style="width: 100%; text-align: left;" 
84
|-
85
| 
86
{| style="text-align: left; margin:auto;width: 100%;" 
87
|-
88
| style="text-align: center;" | <math>\tau _{ij}=2\mu \left(\varepsilon _{ij}-{1\over 3} {\partial u_k\over \partial x_k}\delta _{ij}\right)</math>
89
|}
90
| style="width: 5px;text-align: right;white-space: nowrap;" |  (5a)
91
|}
92
93
with
94
95
{| class="formulaSCP" style="width: 100%; text-align: left;" 
96
|-
97
| 
98
{| style="text-align: left; margin:auto;width: 100%;" 
99
|-
100
| style="text-align: center;" | <math>\varepsilon _{ij}={1\over 2}\left({\partial u_i\over \partial x_j}+{\partial u_j\over \partial x_i}\right)</math>
101
|}
102
| style="width: 5px;text-align: right;white-space: nowrap;" |  (5b)
103
|}
104
105
Summation convention for repeated indexes in products and derivatives is used unless otherwise specified.
106
107
The field equations are completed by the standard boundary conditions
108
109
{| class="formulaSCP" style="width: 100%; text-align: left;" 
110
|-
111
| 
112
{| style="text-align: left; margin:auto;width: 100%;" 
113
|-
114
| style="text-align: center;" | <math>\sigma _{ij}n_j - t_i^p = 0 \quad \hbox{ on }  \Gamma _t</math>
115
| style="width: 5px;text-align: right;white-space: nowrap;" | (6)
116
|-
117
| style="text-align: center;" | <math> u_i-u_i^p  = 0 \quad \hbox{ on }  \Gamma _u </math>
118
| style="width: 5px;text-align: right;white-space: nowrap;" | (7)
119
|}
120
|}
121
122
where <math display="inline">t_i^p</math> and <math display="inline">u_i^p</math> are prescribed tractions and the prescribed displacement at the Neuman and Dirichlet boundaries <math display="inline">\Gamma _t</math> and <math display="inline">\Gamma _u</math>, respectively.
123
124
We will choose an equal order finite element (FE) approximation for the velocity and pressure
125
126
{| class="formulaSCP" style="width: 100%; text-align: left;" 
127
|-
128
| 
129
{| style="text-align: left; margin:auto;width: 100%;" 
130
|-
131
| style="text-align: center;" | <math>u_i\simeq \hat u_i =  \sum \limits _{j=1}^n N_j (\hat u_i)_j = {N}\hat{u}_i </math>
132
|}
133
| style="width: 5px;text-align: right;white-space: nowrap;" |  (8a)
134
|}
135
136
{| class="formulaSCP" style="width: 100%; text-align: left;" 
137
|-
138
| 
139
{| style="text-align: left; margin:auto;width: 100%;" 
140
|-
141
| style="text-align: center;" | <math>p\simeq \hat p =   \sum \limits _{j=1}^n N_j \hat p_j = {N}\hat{p}</math>
142
|}
143
| style="width: 5px;text-align: right;white-space: nowrap;" |  (8b)
144
|}
145
146
where <math display="inline">\hat {(\cdot )}</math> denotes approximate values, <math display="inline">N_j</math> are the shape functions, <math display="inline">(\cdot )_j</math> are nodal values and <math display="inline">n</math> is the number of nodes in the mesh <span id='citeF-1'></span>[[#cite-1|1]]. Standard  three node linear triangles are used in the examples given in the paper.
147
148
It is well known that stabilized finite element schemes are needed for overcoming the numerical instabilities induced by high values of the convective terms and by the incompressibility constraint when equal order interpolation for the velocities and the presure are used. A number of stabilized FE methods are available in the literature and most of them can be casted in the following form
149
150
{| class="formulaSCP" style="width: 100%; text-align: left;" 
151
|-
152
| 
153
{| style="text-align: left; margin:auto;width: 100%;" 
154
|-
155
| style="text-align: center;" | <math>\int _\Omega N_k \hat r_{m_i} d\Omega +\int _{\Gamma _t} N_k (\hat \sigma _{ij} n_j - t_i^p)d\Gamma +\sum _e \int _{\Omega ^e}\tau _{m_i}^e P_{m_i} (N_k)\hat r_{m_i} d\Omega =0 </math>
156
|}
157
| style="width: 5px;text-align: right;white-space: nowrap;" | (9)
158
|}
159
160
{| class="formulaSCP" style="width: 100%; text-align: left;" 
161
|-
162
| 
163
{| style="text-align: left; margin:auto;width: 100%;" 
164
|-
165
| style="text-align: center;" | <math>\int _\Omega N_k \hat r_v d\Omega + \sum _e \int _{\Omega ^e}\tau _v ^e P_v (N_k)\hat r_v d\Omega =0 </math>
166
|}
167
| style="width: 5px;text-align: right;white-space: nowrap;" | (10)
168
|}
169
170
In above <math display="inline">\hat r_{m_i}\equiv r_{m_i}(\hat u_j, \hat p)</math> and <math display="inline">\hat r_v =r_v(\hat u_j, \hat p)</math> denote the approximate finite element residuals for the momentum and mass balance equations, respectively, <math display="inline"> P_{m_i} (N_k)</math> and <math display="inline">P_v (N_k)</math> are appropriate functions of the shape functions <math display="inline">N_k</math> and <math display="inline">\tau _{m_i}^e</math> and <math display="inline">\tau _v ^e</math> are stabilization parameters which depend on the element dimensions and the physical variables of the problem.
171
172
Eqs.(9) and (10) can be interpreted as  extended  Galerkin expressions where the last integrals over the element interiors has been added to the standard Galerkin forms of the momentum and mass balance equations in order to get a physically meanful (stabilized) numerical solution. Many popular stabilized methods (SUPG, GLS, SGS, etc.) can be written in the form of Eqs.(9) and (10) and they just differ in the choice of the weighting functions <math display="inline">P_{m_i}</math> and <math display="inline">P_v</math> and in the values of the stabilization parameters <span id='citeF-1'></span>[[#cite-1|1]].
173
174
A somehow different class of stabilized method developed by the authors is based on a finite calculus (FIC) formulation <span id='citeF-17'></span>[[#cite-17|17]]. Here the standard momentum and mass balance equations are modified by enforcing the balance laws in a domain of finite size. The modified governing equations contain additional terms which play the role of stabilizing terms in the numerical solution. The Galerkin FE form of the FIC equations of an incompressible fluid can also be written by integral forms analogous to Eqs.(9) and (10) and hence we will retain these equations as general expressions applicable to any stabilized method. Details on the FIC formulation used for the solution of the examples presented in the paper are given in a later section.
175
176
==3 RESIDUAL POWER FORMS BASED ON THE RECOVERED BALANCE EQUATIONS==
177
178
We will assume that a finite element solution has been found for the velocity and pressure variables ''using any stabilized FE scheme''. The next step in the error estimation and mesh adaptivity method proposed is to compute the ''enhanced recovered values'' of the derivatives of the velocity and the pressure at the element nodes. The derivative recovery procedure used here is explained in the next section.
179
180
The recovered nodal derivative fields of the velocities and the pressure are used now to compute the finite element residuals of the discretized momentum and mass balance equations (1) and (2) and the Neumann boundary conditions (6) as:
181
182
{| class="formulaSCP" style="width: 100%; text-align: left;" 
183
|-
184
| 
185
{| style="text-align: left; margin:auto;width: 100%;" 
186
|-
187
| style="text-align: center;" | <math>\begin{array}{c}\displaystyle R_{m_i}= \bar{\hat r}_{m_i} \\  ~~\\ \displaystyle R_v=\bar{\hat r}_v \\  ~~\\ \displaystyle R_{\Gamma _i}=n_j \bar{\hat \sigma }_{ij}-t_i^p \end{array} </math>
188
|}
189
| style="width: 5px;text-align: right;white-space: nowrap;" | (11)
190
|}
191
192
where <math display="inline">\bar{\hat{(\cdot )}}</math> denote values computed using the recovered derivative fields. In Eqs.(11) <math display="inline">R_{m_i}</math> can be interpreted as the ith component of an out-of-balance body force, whereas <math display="inline">R_v</math> represents an spurious non-zero volumetric strain rate. Similarly, the recovered residuals <math display="inline">R_{\Gamma _i}</math> are out-of-balance tractions at the Neuman boundary.
193
194
The ''power'' due to the recovered residuals can be written as
195
196
{| class="formulaSCP" style="width: 100%; text-align: left;" 
197
|-
198
| 
199
{| style="text-align: left; margin:auto;width: 100%;" 
200
|-
201
| style="text-align: center;" | <math>\displaystyle P_{m_i}= \int _\Omega W_i^u  R_{m_i} d\Omega + \int _{\Gamma _t} \bar W_i^u  R_{\Gamma _i} d\Gamma \qquad \hbox{no sum in }i</math>
202
|}
203
| style="width: 5px;text-align: right;white-space: nowrap;" |  (12a)
204
|}
205
206
{| class="formulaSCP" style="width: 100%; text-align: left;" 
207
|-
208
| 
209
{| style="text-align: left; margin:auto;width: 100%;" 
210
|-
211
| style="text-align: center;" | <math>\displaystyle P_v= \int _\Omega \bar W^p  R_v d\Omega </math>
212
|}
213
| style="width: 5px;text-align: right;white-space: nowrap;" |  (12b)
214
|}
215
216
where the weights <math display="inline">(W_i^u,\bar W_i^u)</math> and <math display="inline">W^p</math> are functions of the approximate velocity vector field and the pressure, respectively. Indeed the simplest choice is <math display="inline">W_i^u=\bar W_i^u =\hat u_i</math> and <math display="inline">W^p =\hat p</math>. However, many other options are possible. In this work we have chosen
217
218
{| class="formulaSCP" style="width: 100%; text-align: left;" 
219
|-
220
| 
221
{| style="text-align: left; margin:auto;width: 100%;" 
222
|-
223
| style="text-align: center;" | <math>W_i^u=\hat u_i + {h_j\over 2} {\partial \hat u_i\over \partial x_j}\quad ,\quad \bar W_i^u  = 0\quad \hbox{and} \quad  W^p =\hat p + {h_j\over 2} {\partial \hat p\over \partial x_j} </math>
224
|}
225
| style="width: 5px;text-align: right;white-space: nowrap;" | (13)
226
|}
227
228
where <math display="inline">h_j</math> is a characteristic length of the order of the element size. This choice is a result of the ''finite calculus'' formulation described in a next section.
229
230
Note that at the solid boundaries <math display="inline">\hat {u}=0</math>  and  at the exterior boundary <math display="inline">\hat \sigma _{ij}n_j\simeq t_i</math> on <math display="inline">\Gamma _t</math>.  This justifies neglecting the  boundary integral in Eq.(12a) by choosing <math display="inline">\bar W_i^u  = 0</math>.
231
232
The weights <math display="inline">W_i^u</math> and <math display="inline">W^p</math> as defined in Eq.(13) are discontinuous for standar finite element interpolations. A continuous discrete form of these weights can be found by using the recovered values of the derivatives of <math display="inline">\hat u_i</math> and <math display="inline">\hat p</math> as explained in Section 4.
233
234
The integrals over the domain can be computed by summing up the individual element contributions in the usual manner. However in the examples presented in the paper the volume integrals in Eqs.(12)  are computed using a ''nodal integration scheme'' as
235
236
{| class="formulaSCP" style="width: 100%; text-align: left;" 
237
|-
238
| 
239
{| style="text-align: left; margin:auto;width: 100%;" 
240
|-
241
| style="text-align: center;" | <math>\begin{array}{l}\displaystyle P_{m_i} = \sum \limits _{k=1}^{\bar N} [P_{m_i}]_k= \sum \limits _{k=1}^{\bar N} \left|\left[W_i^u R_{m_i}\right]_k \right|\Omega _k \quad \hbox{no sum in } i\\ \\ \displaystyle P_v= \sum \limits _{k=1}^{\bar N} [P_v]_k= \sum \limits _{k=1}^{\bar N} \left|\left[W_i^p  R_v\right]_k\right|\Omega _k  \end{array} </math>
242
|}
243
| style="width: 5px;text-align: right;white-space: nowrap;" | (14)
244
|}
245
246
where <math display="inline">[\cdot ]_k</math> denotes values computed at node <math display="inline">k</math>, <math display="inline">\bar N</math> is the number of ''internal'' nodes in the mesh and <math display="inline">\Omega _k </math> is the tributary area of node <math display="inline">k</math> (Figure 1). In Eqs.(14) absolute values are used to avoid random cancelation of the residual power terms. Indeed Eqs.(14) imply nodal continuity of the weighting functions obtained via the derivative recovery process explained in the next Section. The global error for a given mesh is estimated as
247
248
{| class="formulaSCP" style="width: 100%; text-align: left;" 
249
|-
250
| 
251
{| style="text-align: left; margin:auto;width: 100%;" 
252
|-
253
| style="text-align: center;" | <math>\Vert P\Vert _\Omega = \sum \limits _{k=1}^{\bar N} \Vert P\Vert _{\Omega ^k}  </math>
254
|}
255
| style="width: 5px;text-align: right;white-space: nowrap;" | (15)
256
|}
257
258
where <math display="inline">\Vert P\Vert _{\Omega ^k} </math> is the estimated error value for an individual node defined as
259
260
{| class="formulaSCP" style="width: 100%; text-align: left;" 
261
|-
262
| 
263
{| style="text-align: left; margin:auto;width: 100%;" 
264
|-
265
| style="text-align: center;" | <math>\Vert P\Vert _{\Omega ^k} = \sum \limits _{i=1}^{n_d} [P_{m_i}]_k +[P_v]_k </math>
266
|}
267
| style="width: 5px;text-align: right;white-space: nowrap;" | (16)
268
|}
269
270
===Remark===
271
272
We note that the residual power expressions chosen can be interpreted as the projection of the recovered residuals of the momentum and incompressibility equations against the weighting functions <math display="inline">P_{m_i}</math> and <math display="inline">P_v</math>, respectively. In this sense, Eqs.(12) can be viewed as a particular class of residual projection method. A review of these procedures can be found in <span id='citeF-15'></span>[[#cite-15|15]] and references herein.
273
274
==4 RECOVERY OF THE NODAL DERIVATIVES OF THE VELOCITIES AND THE PRESSURE==
275
276
The recovered derivative field <math display="inline">\left(\displaystyle \overline{\partial \hat u_i\over \partial x_j}\right)</math> is defined so as to satisfy the following variational form
277
278
{| class="formulaSCP" style="width: 100%; text-align: left;" 
279
|-
280
| 
281
{| style="text-align: left; margin:auto;width: 100%;" 
282
|-
283
| style="text-align: center;" | <math>\int _\Omega N_k \left[\left(\overline{\partial \hat u_i\over \partial x_j}\right)- {\partial \hat u_i\over \partial x_j}\right]d\Omega =0\quad \hbox{no sum in }i </math>
284
|}
285
| style="width: 5px;text-align: right;white-space: nowrap;" | (17)
286
|}
287
288
where <math display="inline">N_k</math> are the same shape functions used to approximate the velocity field and <math display="inline">\displaystyle{\partial \hat u_i\over \partial x_j}</math> are known derivatives values directly obtained from the finite element solution.
289
290
A standard finite element approximation is chosen for the recovered derivative field; i.e.
291
292
{| class="formulaSCP" style="width: 100%; text-align: left;" 
293
|-
294
| 
295
{| style="text-align: left; margin:auto;width: 100%;" 
296
|-
297
| style="text-align: center;" | <math>\left(\overline{\partial \hat u_i\over \partial x_j}\right)=\sum \limits _k N_k \left(\overline{\partial \hat u_i\over \partial x_j}\right)_k </math>
298
|}
299
| style="width: 5px;text-align: right;white-space: nowrap;" | (18)
300
|}
301
302
Substituting Eq.(18) into (17) gives the following system of equations from which the recovered nodal derivative values can be obtained
303
304
{| class="formulaSCP" style="width: 100%; text-align: left;" 
305
|-
306
| 
307
{| style="text-align: left; margin:auto;width: 100%;" 
308
|-
309
| style="text-align: center;" | <math>{M} \bar{\hat {d}}_{ij}= {f}_{ij} </math>
310
|}
311
| style="width: 5px;text-align: right;white-space: nowrap;" | (19)
312
|}
313
314
with
315
316
{| class="formulaSCP" style="width: 100%; text-align: left;" 
317
|-
318
| 
319
{| style="text-align: left; margin:auto;width: 100%;" 
320
|-
321
| style="text-align: center;" | <math>\bar{\hat {d}}_{ij}= \left[\left(\overline{\partial \hat u_i\over \partial x_j}\right)_1,\left(\overline{\partial \hat u_i\over \partial x_j}\right)_2,\cdots , \left(\overline{\partial \hat u_i\over \partial x_j}\right)_N\right]^T </math>
322
|}
323
| style="width: 5px;text-align: right;white-space: nowrap;" | (20)
324
|}
325
326
where <math display="inline">N</math> is the number of nodes in the mesh.
327
328
The element contributions to <math display="inline">{M}</math> and '''f''' are
329
330
{| class="formulaSCP" style="width: 100%; text-align: left;" 
331
|-
332
| 
333
{| style="text-align: left; margin:auto;width: 100%;" 
334
|-
335
| style="text-align: center;" | <math>M_{kl}^e =\int _{\Omega ^e} N_k N_l d\Omega \quad , \quad (f_{ij}^e)_k = \int _{\Omega ^e} N_k {\partial \hat u_i\over \partial x_j} d\Omega  </math>
336
|}
337
| style="width: 5px;text-align: right;white-space: nowrap;" | (21)
338
|}
339
340
We recall that in above <math display="inline">\displaystyle{\partial \hat u_i\over \partial x_j}</math> is the (discontinuous) derivative field obtained from the finite element solution and <math display="inline">\left(\displaystyle \overline{\partial \hat u_i\over \partial x_j}\right)_k</math> are the recovered derivative values at node <math display="inline">k</math>. The process is repeated for <math display="inline">i,j=1,n_d</math>.
341
342
The full form of matrix <math display="inline">{M}</math> is used for the solution of Eq.(19) using a standard Jacobian iterative scheme.
343
344
The same procedure is followed in order to recover the nodal derivatives of the pressure and the nodal second derivatives of the velocity field <math display="inline">\left(\displaystyle{\partial ^2 \hat u_k\over \partial x_i\partial x_j}\right)</math> needed for computation of all the terms in the recovered momentum equations <math display="inline">R_{m_i}</math>.
345
346
==5 ERROR INDICATORS AND MESH ADAPTATION PROCEDURES==
347
348
A ''global error parameter'' is defined as
349
350
{| class="formulaSCP" style="width: 100%; text-align: left;" 
351
|-
352
| 
353
{| style="text-align: left; margin:auto;width: 100%;" 
354
|-
355
| style="text-align: center;" | <math>\xi _g = {\Vert P\Vert _\Omega \over \mu U} </math>
356
|}
357
| style="width: 5px;text-align: right;white-space: nowrap;" | (22)
358
|}
359
360
where <math display="inline">\mu </math> is the prescribed global error value and <math display="inline">U</math> is the total power in the mesh defined here as
361
362
{| class="formulaSCP" style="width: 100%; text-align: left;" 
363
|-
364
| 
365
{| style="text-align: left; margin:auto;width: 100%;" 
366
|-
367
| style="text-align: center;" | <math>U= \int _\Omega \left[\tau _{ij} \varepsilon _{ij} + {\partial p\over \partial x_i}u_i\right]d\Omega  </math>
368
|}
369
| style="width: 5px;text-align: right;white-space: nowrap;" | (23)
370
|}
371
372
The first term in Eq.(23) is the power of the viscous stresses, whereas the second one is the power of the pressure gradient terms. The integral in Eq.(23) is computed using nodal information by an expression analogous to Eq.(15).
373
374
The mesh adaption process aims to reduce the global error (leading to a value of <math display="inline">\xi _g =1</math>) while satisfying a local error condition for each individual element. The two options considered in this work for the local error condition are presented next.
375
376
===5.1 Mesh adaptation based on the equidistribution of the global error===
377
378
A ''local error parameter'' <math display="inline"> \xi _g</math> is defined for each element as
379
380
{| class="formulaSCP" style="width: 100%; text-align: left;" 
381
|-
382
| 
383
{| style="text-align: left; margin:auto;width: 100%;" 
384
|-
385
| style="text-align: center;" | <math>\xi ^e= {\Vert P\Vert _{\Omega ^e} N\over \Vert P\Vert _\Omega } </math>
386
|}
387
| style="width: 5px;text-align: right;white-space: nowrap;" | (24)
388
|}
389
390
where <math display="inline">\Vert P\Vert _{\Omega ^e}</math> is the residual power error for element <math display="inline">e</math> and <math display="inline">N</math> is the total number of elements in the mesh. Eq.(24) is based in the so called equi-distribution of the global error, i.e. if <math display="inline">\xi ^e= 1</math> for all elements then the total residual power error is uniformly distributed in the mesh <span id='citeF-4'></span>[[#cite-4|4]].
391
392
The error indicator for the element is defined as
393
394
{| class="formulaSCP" style="width: 100%; text-align: left;" 
395
|-
396
| 
397
{| style="text-align: left; margin:auto;width: 100%;" 
398
|-
399
| style="text-align: center;" | <math>\gamma ^e = (\xi ^e)^\alpha (\xi _g)^\beta  </math>
400
|}
401
| style="width: 5px;text-align: right;white-space: nowrap;" | (25)
402
|}
403
404
Parameters <math display="inline">\alpha </math> and <math display="inline">\beta </math> are found by analyzing the convergence rate of the global and local error parameters as described in <span id='citeF-4'></span>[[#cite-4|4]]. In the examples presented in the paper the following values of <math display="inline">\alpha </math> and <math display="inline">\beta </math> have been used
405
406
{| class="formulaSCP" style="width: 100%; text-align: left;" 
407
|-
408
| 
409
{| style="text-align: left; margin:auto;width: 100%;" 
410
|-
411
| style="text-align: center;" | <math>\alpha = {1\over p+n_d}\quad ,\quad \beta = {1\over p+1} </math>
412
|}
413
| style="width: 5px;text-align: right;white-space: nowrap;" | (26)
414
|}
415
416
where <math display="inline">p</math> is the order of the finite element interpolation (i.e. <math display="inline">p=1</math> for linear elements).
417
418
The new element size <math display="inline">h_{new}</math> is computed in terms of the error indicator <math display="inline">\gamma ^e</math> and the current element size as
419
420
{| class="formulaSCP" style="width: 100%; text-align: left;" 
421
|-
422
| 
423
{| style="text-align: left; margin:auto;width: 100%;" 
424
|-
425
| style="text-align: center;" | <math>h_{new}= h_{old} [\gamma ^e]^{-1} </math>
426
|}
427
| style="width: 5px;text-align: right;white-space: nowrap;" | (27)
428
|}
429
430
===5.2 Mesh adaptation based on the equidistribution of the error density===
431
432
An alternative local error parameter based on the equidistribution of the density of the residual power can be defined as
433
434
{| class="formulaSCP" style="width: 100%; text-align: left;" 
435
|-
436
| 
437
{| style="text-align: left; margin:auto;width: 100%;" 
438
|-
439
| style="text-align: center;" | <math>\xi ^e= {\Omega \Vert P\Vert _{\Omega ^e} \over \Omega ^e \Vert P\Vert _\Omega } </math>
440
|}
441
| style="width: 5px;text-align: right;white-space: nowrap;" | (28)
442
|}
443
444
where <math display="inline">\Omega </math> and <math display="inline">\Omega ^e</math> are the area of the analysis domain and the area of an individual element, respectively. Now <math display="inline">\xi ^e = 1</math> if the density of the error (error per unit area) is uniformly distributed over the mesh.
445
446
The expression for the error indicator for the element is again given by Eq.(25). The following values of <math display="inline">\alpha </math> and <math display="inline">\beta </math> are now chosen as deduced following the procedure described in <span id='citeF-4'></span>[[#cite-4|4]]
447
448
{| class="formulaSCP" style="width: 100%; text-align: left;" 
449
|-
450
| 
451
{| style="text-align: left; margin:auto;width: 100%;" 
452
|-
453
| style="text-align: center;" | <math>\alpha = {1\over p} \quad ,\quad \beta = {1\over p} </math>
454
|}
455
| style="width: 5px;text-align: right;white-space: nowrap;" | (29)
456
|}
457
458
The new element sizes are computed  from Eq.(27).
459
460
This mesh adaptation criteria leads to meshes where the smaller elements are concentrated in zones where the higher gradients of the velocity and pressure occur, whereas much larger elemens are placed in other zones of the domain <span id='citeF-4'></span>[[#cite-4|4]]. This is adequate for capturing accurate solutions in the vecinity of boundary layers or singularity zones. The drawback of this procedure is that it tends to refine progressively (and very fast) the mesh in these zones until <math display="inline">\xi ^e=1</math>. This  leads to regions with very small elements and hence to an extremely large number of elements in the mesh. A limiting value in the size of the smallest element in the new mesh is mandatory in these cases.
461
462
Alternatively the criterium based on the equidistribution of the global error described in the previous section leads to meshes with a “smoother” distribution of the element sizes between the regions where high gradients occur and the other zones in the domain.
463
464
A comparison of these two mesh adaptions procedures for linear structural analysis can be found in <span id='citeF-4'></span>[[#cite-4|4]].
465
466
===5.3 Mesh refinement strategy===
467
468
The mesh refinement strategy is as follows.
469
470
'''Step 1''' Compute the nodal recovered values of the derivatives of the velocity and pressure fields using the algorithm described in  Section 4. '''Step 2''' Compute the residual power error for each element in the mesh and the total residual power error using Eqs.(15) and (16). '''Step 3''' Compute the global and local (element) error parameters <math display="inline">\xi _g</math> and <math display="inline">\xi ^e</math>. '''Step 4''' Compute the new element sizes using Eq.(27). Perform a smoothing of the new element sizes to avoid problems during the regeneration of the mesh. '''Step 5''' Regenerate a new mesh. Minimum and maximum sizes of the elements in the mesh are prescribed in order to have a better control of the mesh generation process. The mesh generation  is carried out using the GiD pre/postprocessing system <span id='citeF-23'></span>[[#cite-23|23]].
471
472
Above mesh refinement strategy is equally applicable for the two mesh adaption procedures explained in the previous sections.
473
474
For steady-state flow problems the error estimation and the mesh refinement process are performed after a steady state solution has been found. For transient problems the mesh refinement  is performed after a prescribed number of time steps.
475
476
==6 STABILIZED FEM FOR THE NAVIER-STOKES EQUATIONS USING FINITE CALCULUS==
477
478
===6.1 Finite calculus form of the Navier-Stokes equations===
479
480
The FIC governing equations for incompressible viscous flows are obtained by applying the standard conservation laws expressing balance of momentum and mass over a control domain. Assuming that the control domain has finite dimensions and expressing the variation of mass and momentum over the domain using Taylor series expansions of one order higher than those used in the standard infinitesimal  theory, the following expressions are found <span id='citeF-17'></span>[[#cite-17|17]]:<br/>
481
482
''Momentum''
483
484
{| class="formulaSCP" style="width: 100%; text-align: left;" 
485
|-
486
| 
487
{| style="text-align: left; margin:auto;width: 100%;" 
488
|-
489
| style="text-align: center;" | <math>r_{m_i} - \underline{{1\over 2} h_{j} {\partial r_{m_i}\over \partial x_ j}}{1\over 2} h_{j} {\partial r_{m_i}\over \partial x_ j} =0 \qquad \hbox{in } \Omega  </math>
490
|}
491
| style="width: 5px;text-align: right;white-space: nowrap;" | (30)
492
|}
493
494
''Mass balance''
495
496
{| class="formulaSCP" style="width: 100%; text-align: left;" 
497
|-
498
| 
499
{| style="text-align: left; margin:auto;width: 100%;" 
500
|-
501
| style="text-align: center;" | <math>r_v - \underline{{1\over 2} h_{j}{\partial r_v \over \partial x_j}}{1\over 2} h_{j}{\partial r_v \over \partial x_j} =0 \qquad \hbox{in } \Omega  </math>
502
|}
503
| style="width: 5px;text-align: right;white-space: nowrap;" | (31)
504
|}
505
506
where <math display="inline">r_{m_i}</math> and <math display="inline">r_v</math> were defined in Eqs.(3) and (4), respectively.
507
508
Eqs.(30) and (31) are the ''FIC forms'' of the governing differential equations for an incompressible flow. The terms underlined in Eqs.(30) and (31) introduce naturally the necessary stabilization at the discretization level. The so called ''characteristic length'' vector <math display="inline">{h}</math> is defined as (for 2D problems)
509
510
{| class="formulaSCP" style="width: 100%; text-align: left;" 
511
|-
512
| 
513
{| style="text-align: left; margin:auto;width: 100%;" 
514
|-
515
| style="text-align: center;" | <math>{h}=\left\{\begin{matrix} h_{1}\\ h_{2}\\\end{matrix}\right\} </math>
516
|}
517
| style="width: 5px;text-align: right;white-space: nowrap;" | (32)
518
|}
519
520
where <math display="inline">h_{1}</math> and <math display="inline">h_{2}</math> are characteristic dimensions of the finite control domain where balance of momentum and mass conservation is enforced. The components of vector <math display="inline">{h}</math> introduce the necessary stabilization along the streamline and transverse directions to the flow in the discrete problem <span id='citeF-17'></span>[[#cite-17|17]].
521
522
The method to derive the modified differential equations (30) and (31) incorporating the stabilization terms is termed ''finite calculus'' (FIC) as a reference to the standard infinitesimal calculus  where the size of the domain where balance of mass and momentum is enforced is assumed to be negligible. Note that for <math display="inline">{h}\to 0</math> the standard infinitesimal form of the momentum and mass balance equations is recovered <span id='citeF-1'></span>[[#cite-1|1]].
523
524
Eqs.(30) and (31) are completed by the following FIC boundary conditions <span id='citeF-17'></span>[[#cite-17|17]].
525
526
<br/>
527
528
''Equilibrium of surface tractions  at the boundary'' <math display="inline">\Gamma _t</math>
529
530
{| class="formulaSCP" style="width: 100%; text-align: left;" 
531
|-
532
| 
533
{| style="text-align: left; margin:auto;width: 100%;" 
534
|-
535
| style="text-align: center;" | <math>n_j \sigma _{ij}-t_i^p + \underline{{1\over 2} h_{j}n_j r_{m_i}}{1\over 2} h_{j}n_j r_{m_i} = 0\qquad \hbox{on}~~\Gamma _t  </math>
536
|}
537
| style="width: 5px;text-align: right;white-space: nowrap;" | (33)
538
|}
539
540
''Prescribed velocity at  the boundaries ''
541
542
{| class="formulaSCP" style="width: 100%; text-align: left;" 
543
|-
544
| 
545
{| style="text-align: left; margin:auto;width: 100%;" 
546
|-
547
| style="text-align: center;" | <math>u_t =u_t^p \qquad \hbox{on}~~\Gamma _{u_t} </math>
548
|}
549
| style="width: 5px;text-align: right;white-space: nowrap;" | (34)
550
|}
551
552
{| class="formulaSCP" style="width: 100%; text-align: left;" 
553
|-
554
| 
555
{| style="text-align: left; margin:auto;width: 100%;" 
556
|-
557
| style="text-align: center;" | <math>u_n - \underline{{1\over 2}h_{i}n_ir_v}{1\over 2}h_{i}n_ir_v=u_n^p \qquad \hbox{on} ~~\Gamma _{u_n} </math>
558
|}
559
| style="width: 5px;text-align: right;white-space: nowrap;" | (35)
560
|}
561
562
In eq.(34) <math display="inline">u_t</math> and <math display="inline">u_t^p</math> denote the tangential velocity to the boundary and its prescribed value, respectively.
563
564
Eq.(35) expresses the balance of mass on an arbitrary domain next to the boundary and <math display="inline">u_n</math> and <math display="inline">u_n^p</math> denote the velocity normal to the boundary and its prescribed value, respectively. The value of <math display="inline">u_n^p</math> is  zero on solid walls and stationary free surfaces.
565
566
Also in eqs.(34) and (35) <math display="inline">\Gamma _{u_t}</math> and <math display="inline">\Gamma _{u_n}</math> are the parts of the boundary <math display="inline">\Gamma </math> of <math display="inline">\Omega </math> where the tangential and normal velocities are prescribed, respectively. The Dirichlet boundary is defined as <math display="inline">\Gamma _u=\Gamma _{u_t}\cup \Gamma _{u_n}</math>.
567
568
The underlined terms in eqs.(33) and (35) introduce the necessary stabilization at the boundaries in a form consistent with that of eqs.(30) and (31). Details of the derivation of eqs.(30&#8211;35) can be found in <span id='citeF-17'></span>[[#cite-17|17]] whereas the derivation of eq.(9) is shown in <span id='citeF-18'></span>[[#cite-18|18]].
569
570
===6.2 Residual power forms of FIC equations===
571
572
We will choose again an equal order finite element approximation for the velocity and pressure variables given by Eqs.(8)
573
574
Again we will assume that a finite element solution has been found for the velocity and the pressure variables. The finite element algorithm used for the examples solved in this paper is based on a ''fractional step scheme'' <span id='citeF-19'></span>[[#cite-19|19]]. The next step is to obtain enhanced recovered values of the derivatives of the velocity and the pressure at the element nodes as explained in previouss section.
575
576
The recovered nodal derivative fields of the velocities and the pressure are used now to compute the FIC residuals of the momentum and mass balance equations (1) and (2):
577
578
{| class="formulaSCP" style="width: 100%; text-align: left;" 
579
|-
580
| 
581
{| style="text-align: left; margin:auto;width: 100%;" 
582
|-
583
| style="text-align: center;" | <math>\begin{array}{c}\displaystyle R_{m_i}= \bar{\hat r}_{m_i} - {h_{j} \over 2} {\partial \bar{\hat r}_{m_i}\over \partial x_j}\\  \\ \displaystyle R_v=\bar{\hat r}_v - {h_{j}\over 2} {\partial \bar{\hat r}_v \over \partial x_j}  \end{array} </math>
584
|}
585
| style="width: 5px;text-align: right;white-space: nowrap;" | (36)
586
|}
587
588
where once more <math display="inline">\bar{\hat{(\cdot )}}</math> denote values computed using the recovered derivative fields.
589
590
The total  power due to the FIC residuals can be written as
591
592
{| class="formulaSCP" style="width: 100%; text-align: left;" 
593
|-
594
| 
595
{| style="text-align: left; margin:auto;width: 100%;" 
596
|-
597
| style="text-align: center;" | <math>\begin{array}{c}\displaystyle P_{m_i}= \int _\Omega \hat u_i R_{m_i} d\Omega + \int _{\Gamma _t} \hat u_i  \left(n_j \bar{\hat \sigma }_{ij} - t_i + {h_{j}\over 2} n_j \bar{\hat r}_{m_i}\right)d\Gamma \\ \\ \displaystyle P_v= \int _\Omega \hat p R_v d\Omega - \int _{\Gamma _{u_n}} \hat p \left(\hat u_n -u_n^p - {h_{j}\over 2}n_j  \bar{\hat r}_v\right)d\Gamma  \end{array} </math>
598
|}
599
| style="width: 5px;text-align: right;white-space: nowrap;" | (37)
600
|}
601
602
Substituting the expressions of <math display="inline">R_{m_i}</math> and <math display="inline">R_v</math> of Eqs.(36) into (37) and integrating by parts some derivative terms gives
603
604
{| class="formulaSCP" style="width: 100%; text-align: left;" 
605
|-
606
| 
607
{| style="text-align: left; margin:auto;width: 100%;" 
608
|-
609
| style="text-align: center;" | <math>P_{m_i}= \int _\Omega \left[\hat u_i + {h_j\over 2}\left({\partial \hat u_i \over  \partial x_j} \right)\right]\bar{\hat r}_{m_i} d\Omega +  \int _{\Gamma _t} \hat u_i \left(\bar{\hat \sigma }_{ij}n_j  -t_i\right)d\Gamma  - \int _{\Gamma _u}  {h_j\over 2} n_j \hat u_i \bar{\hat r}_{m_i} d\Gamma \quad \hbox{no sum in } i</math>
610
|}
611
| style="width: 5px;text-align: right;white-space: nowrap;" |  (38a)
612
|}
613
614
{| class="formulaSCP" style="width: 100%; text-align: left;" 
615
|-
616
| 
617
{| style="text-align: left; margin:auto;width: 100%;" 
618
|-
619
| style="text-align: center;" | <math>P_v= \int _\Omega \left[\hat p  +{h_j\over 2} \left({\partial \hat p\over \partial x_j}\right)\right] \bar{\hat r}_v d\Omega - \int _{\Gamma _{u_n}} \hat p (\hat u_n - u_n^p)d\Gamma - \int _{\Gamma -\Gamma _{u_n}}\hat p {h_j\over 2}n_j \bar{\hat r}_v d\Gamma </math>
620
|}
621
| style="width: 5px;text-align: right;white-space: nowrap;" |  (38b)
622
|}
623
624
Numerical experiments have  shown that the contribution of the boundary integrals in Eqs.(38) is negligible versus the value of the integrals over the domain. Hence the contribution of the surface integrals in Eqs.(38) will not be taken into account from now onwards. An explanation for this is that at the solid boundaries  in a viscous flow <math display="inline">\hat u_j=0</math> and hence the line integrals in Eq.(38a) vanish there. Note also that at the exterior boundary <math display="inline">\hat \sigma _{ij}n_j\simeq t_i</math> on <math display="inline">\Gamma _t</math> and <math display="inline">r_{m_i}\simeq 0</math> on <math display="inline">\Gamma _u</math> and this justifies neglecting the line integral at those boundaries. In Eq.(38b) the first line integral vanishes if <math display="inline">\hat u_n</math> is made equal to the prescribed normal velocity <math display="inline">u_n^p</math> when solving the system of discretized equations.  The second integral in Eq.(38b) can be neglected by accepting that either <math display="inline">p=0</math> (external boundary) or <math display="inline">\bar{\hat r}_v \simeq 0</math> at the <math display="inline">\Gamma - \Gamma _{u_n}</math> boundary.
625
626
In summary, the relevant residual power terms for the FIC momentum and mass balance equations are
627
628
{| class="formulaSCP" style="width: 100%; text-align: left;" 
629
|-
630
| 
631
{| style="text-align: left; margin:auto;width: 100%;" 
632
|-
633
| style="text-align: center;" | <math>\displaystyle P_{m_i}= \int _\Omega \left[\hat u_i + {h_j\over 2} \left({\partial \hat u_i \over  \partial x_j}\right)\right]\bar{\hat r}_{m_i} d\Omega \qquad \hbox{no sum in }i</math>
634
|}
635
| style="width: 5px;text-align: right;white-space: nowrap;" |  (39a)
636
|}
637
638
{| class="formulaSCP" style="width: 100%; text-align: left;" 
639
|-
640
| 
641
{| style="text-align: left; margin:auto;width: 100%;" 
642
|-
643
| style="text-align: center;" | <math>\displaystyle P_v= \int _\Omega \left[\hat p + {h_j\over 2} \left({\partial \hat p\over \partial x_j}\right)\right]\bar{\hat r}_v  d\Omega </math>
644
|}
645
| style="width: 5px;text-align: right;white-space: nowrap;" |  (39b)
646
|}
647
648
In the examples presented in the paper the integrals in Eq.(39) are computed using a nodal integration scheme. This requires computing the recovered nodal value of the derivative of the velocity and pressure field. This is an unexpensive task as the recovered nodal values are needed to evaluate the recovered residual. The nodal computation of the residual power expressions is therefore written as
649
650
{| class="formulaSCP" style="width: 100%; text-align: left;" 
651
|-
652
| 
653
{| style="text-align: left; margin:auto;width: 100%;" 
654
|-
655
| style="text-align: center;" | <math>\begin{array}{l}\displaystyle P_{m_i}= \sum \limits _{k=1}^{\bar N}\left|\left[\left(\hat u_i + {h_j\over 2}\left(\overline{\partial \hat u_i \over  \partial x_j}\right)\right)\bar{\hat r}_{m_i}\right]_k\right|\Omega _k\\ \\ \displaystyle P_v= \sum \limits _{k=1}^{\bar N} \left|\left[\left(\hat p  +{h_j\over 2} \left(\overline{\partial \hat p\over \partial x_j}\right)\right)\bar{\hat r}_v \right]_k \right|\Omega _k  \end{array} </math>
656
|}
657
| style="width: 5px;text-align: right;white-space: nowrap;" | (40)
658
|}
659
660
the different terms have the meaning explained for Eq.(14).
661
662
We note that in the examples presented next the characteristic length parameters have been defined as the characteristic sizes of  each tributary domain for a node (see Figure 1).
663
664
<div id='img-1'></div>
665
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
666
|-
667
|[[Image:Draft_Samper_249588187-figura0.png|300px|Tributary domain for a node k and definition of characteristic lengths.]]
668
|- style="text-align: center; font-size: 75%;"
669
| colspan="1" | '''Figure 1:''' Tributary domain for a node <math>k</math> and definition of characteristic lengths.
670
|}
671
672
The global and element error norm are defined by Eqs.(15) and (16), respectively. The error estimation process and the mesh adaptivity scheme follow precisely the steps explained in Section 5.
673
674
===Remark===
675
676
Eqs.(40) define the weighting functions in the general residual power forms of Eqs.(12) and (13) as
677
678
{| class="formulaSCP" style="width: 100%; text-align: left;" 
679
|-
680
| 
681
{| style="text-align: left; margin:auto;width: 100%;" 
682
|-
683
| style="text-align: center;" | <math>\begin{array}{l}\displaystyle{W_i^u = \hat u_i + {h_j\over 2} \left(\overline{\partial \hat u_i \over  \partial x_j}\right)} \quad , \quad \bar W_i^u  = 0 \\ \\ \displaystyle{W^p  = \hat p +  {h_j\over 2}  \left(\overline{\partial \hat p\over \partial x_j}\right)} \end{array} </math>
684
|}
685
| style="width: 5px;text-align: right;white-space: nowrap;" | (41)
686
|}
687
688
These expressions, obtained using arguments from the FIC formulation, can be used to compute variational residual forms of the error in order to derive effective mesh refinement procedures for other problems of computational mechanics using the general scheme presented in Section 5.3 <span id='citeF-19'></span>[[#cite-19|19]].
689
690
We note again however that Eqs.(41) are not the only option for <math display="inline">W_i^u</math> and <math display="inline">W^p</math> and some other choices are explored in <span id='citeF-22'></span>[[#cite-22|22]] to derive a variety of error estimation and mesh refinement procedures for convection-diffusion problems. Numerical results reported in <span id='citeF-22'></span>[[#cite-22|22]] show that the expressions given by Eqs.(41) have a better performance.
691
692
==7 EXAMPLES==
693
694
Typical parameters used for the  examples presented are: prescribed global  error <math display="inline">\mu = 1%</math>, maximum and minimum sizes of the new elements generated <math display="inline">h_{max}=0.3</math>m and <math display="inline">h_{min}=0.001</math>m.
695
696
===7.1 Flow past a NACA airfoil===
697
698
This a transient problem. The mesh refinement process has been carried out at a time <math display="inline">t=2s</math>. The value of the Reynolds number for this problem is 100. Figure 2 shows the initial mesh of 2574 nodes and 5148 linear elements used.
699
700
<div id='img-2'></div>
701
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
702
|-
703
|[[Image:Draft_Samper_249588187-fig1.png|600px|Flow past a NACA airfoil. Characteristic length =1m, Re=100. Original mesh]]
704
|- style="text-align: center; font-size: 75%;"
705
| colspan="1" | '''Figure 2:''' Flow past a NACA airfoil. Characteristic length =1m, <math>Re=100</math>. Original mesh
706
|}
707
708
The mesh adaptation procedure chosen first is that of  equidistribution of the global error. Figures 3 and 4 show the mesh obtained after five refinement steps. The distribution of the velocity field for that mesh is shown in Figures 5 and 6. Table 1 shows the value of the global error and the evolution of the number of elements and nodes during the mesh refinement process. The percentage errors in the drag and lift values are also given in Table 1. We note that the solution was found using a laminar flow model. Therefore some error in the drag values is to be expected even for fine meshes. The error of 6.4% obtained in this case is quite acceptable.
709
710
<div id='img-3'></div>
711
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
712
|-
713
|[[Image:Draft_Samper_249588187-fig2.png|600px|Mesh obtained after five refinement steps using the equidistribution of the global error.]]
714
|- style="text-align: center; font-size: 75%;"
715
| colspan="1" | '''Figure 3:''' Mesh obtained after five refinement steps using the equidistribution of the global error.
716
|}
717
718
<div id='img-4'></div>
719
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
720
|-
721
|[[Image:Draft_Samper_249588187-fig3.png|468px|]]
722
|[[Image:Draft_Samper_249588187-fig3_b.png|534px|Details of the refined mesh of Figure 3 in the vecinity of the airfoil]]
723
|- style="text-align: center; font-size: 75%;"
724
| colspan="2" | '''Figure 4:''' Details of the refined mesh of Figure 3 in the vecinity of the airfoil
725
|}
726
727
<div id='img-5'></div>
728
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
729
|-
730
|[[Image:Draft_Samper_249588187-naca_line_1.png|600px|Velocity field after five refinement steps]]
731
|- style="text-align: center; font-size: 75%;"
732
| colspan="1" | '''Figure 5:''' Velocity field after five refinement steps
733
|}
734
735
<div id='img-6'></div>
736
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
737
|-
738
|[[Image:Draft_Samper_249588187-naca_line_2.png|510px|Detail of the velocity field]]
739
|- style="text-align: center; font-size: 75%;"
740
| colspan="1" | '''Figure 6:''' Detail of the velocity field
741
|}
742
743
744
{|  class="floating_tableSCP wikitable" style="text-align: center; margin: 1em auto;min-width:50%;"
745
|+ style="font-size: 75%;" |Table. 1 Flow past a NACA airfoil (equidistribution of global error). Convergence of the global error with the mesh refinement. Minimum element size = 0.001m
746
|- style="border-top: 2px solid;"
747
| style="border-left: 2px solid;border-right: 2px solid;" |  Mesh 
748
| style="border-left: 2px solid;border-right: 2px solid;" | Global error <math display="inline">\left({\Vert P\Vert _\Omega \over U}\right)</math>(%)  
749
| style="border-left: 2px solid;border-right: 2px solid;" | Nodes 
750
| style="border-left: 2px solid;border-right: 2px solid;" | Elements
751
| style="border-left: 2px solid;border-right: 2px solid;" | Drag error %
752
| style="border-left: 2px solid;border-right: 2px solid;" | Lift error% 
753
|- style="border-top: 2px solid;"
754
| style="border-left: 2px solid;border-right: 2px solid;" |  1
755
| style="border-left: 2px solid;border-right: 2px solid;" | 5.55
756
| style="border-left: 2px solid;border-right: 2px solid;" | 2574
757
| style="border-left: 2px solid;border-right: 2px solid;" | 4784
758
| style="border-left: 2px solid;border-right: 2px solid;" | 29.1810942 
759
| style="border-left: 2px solid;border-right: 2px solid;" | 2.99902143
760
|- style="border-top: 2px solid;"
761
| style="border-left: 2px solid;border-right: 2px solid;" |  2
762
| style="border-left: 2px solid;border-right: 2px solid;" | 2.93
763
| style="border-left: 2px solid;border-right: 2px solid;" | 5340
764
| style="border-left: 2px solid;border-right: 2px solid;" | 10680
765
| style="border-left: 2px solid;border-right: 2px solid;" | 28.7443059
766
| style="border-left: 2px solid;border-right: 2px solid;" | 2.5500113
767
|- style="border-top: 2px solid;"
768
| style="border-left: 2px solid;border-right: 2px solid;" |  3
769
| style="border-left: 2px solid;border-right: 2px solid;" | 2.87
770
| style="border-left: 2px solid;border-right: 2px solid;" | 7942
771
| style="border-left: 2px solid;border-right: 2px solid;" | 15884
772
| style="border-left: 2px solid;border-right: 2px solid;" | 22.1078024
773
| style="border-left: 2px solid;border-right: 2px solid;" | 2.00851784
774
|- style="border-top: 2px solid;"
775
| style="border-left: 2px solid;border-right: 2px solid;" |  4
776
| style="border-left: 2px solid;border-right: 2px solid;" | 2.77
777
| style="border-left: 2px solid;border-right: 2px solid;" | 11948
778
| style="border-left: 2px solid;border-right: 2px solid;" | 23896
779
| style="border-left: 2px solid;border-right: 2px solid;" | 15.2893198
780
| style="border-left: 2px solid;border-right: 2px solid;" | 1.54898957
781
|- style="border-top: 2px solid;"
782
| style="border-left: 2px solid;border-right: 2px solid;" |  5
783
| style="border-left: 2px solid;border-right: 2px solid;" | 1.81
784
| style="border-left: 2px solid;border-right: 2px solid;" | 17142
785
| style="border-left: 2px solid;border-right: 2px solid;" | 34284
786
| style="border-left: 2px solid;border-right: 2px solid;" | 10.6478917
787
| style="border-left: 2px solid;border-right: 2px solid;" | 1.17033653
788
|- style="border-top: 2px solid;"
789
| style="border-left: 2px solid;border-right: 2px solid;" |  6
790
| style="border-left: 2px solid;border-right: 2px solid;" | 1.81
791
| style="border-left: 2px solid;border-right: 2px solid;" | 17255
792
| style="border-left: 2px solid;border-right: 2px solid;" | 34510
793
| style="border-left: 2px solid;border-right: 2px solid;" | 8.54878424
794
| style="border-left: 2px solid;border-right: 2px solid;" | 1.02567233
795
|- style="border-top: 2px solid;"
796
| style="border-left: 2px solid;border-right: 2px solid;" |  7
797
| style="border-left: 2px solid;border-right: 2px solid;" | 1.71
798
| style="border-left: 2px solid;border-right: 2px solid;" | 17372
799
| style="border-left: 2px solid;border-right: 2px solid;" | 34744
800
| style="border-left: 2px solid;border-right: 2px solid;" | 6.93290243
801
| style="border-left: 2px solid;border-right: 2px solid;" | 0.86514698
802
|- style="border-top: 2px solid;border-bottom: 2px solid;"
803
| style="border-left: 2px solid;border-right: 2px solid;" |  8
804
| style="border-left: 2px solid;border-right: 2px solid;" | 1.72
805
| style="border-left: 2px solid;border-right: 2px solid;" | 16729
806
| style="border-left: 2px solid;border-right: 2px solid;" | 33458
807
| style="border-left: 2px solid;border-right: 2px solid;" | 6.40377809
808
| style="border-left: 2px solid;border-right: 2px solid;" | 0.97782406
809
810
|}
811
812
The same problem was solved using the criterium of the equidistribution of error density. The same maximum and minimum sizes for the generated elements are used (i.e., <math display="inline">h_{max}=0.3</math>m and <math display="inline">h_{min}=0.001</math>m). In this case,  after the third step the mesher  was unable to generate a new grid from the element size distributions given by the method. Figures 7 and 8 show the mesh obtained after the second step. The failure of the mesh generator is due to the very large gradients of the element size field leading to an extremely high density of elements concentrated near the airfoil.
813
814
<div id='img-7'></div>
815
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
816
|-
817
|[[Image:Draft_Samper_249588187-fig6.png|600px|Mesh obtained after two refinement steps using the equidistribution of the error density (min size = 0.001m)]]
818
|- style="text-align: center; font-size: 75%;"
819
| colspan="1" | '''Figure 7:''' Mesh obtained after two refinement steps using the equidistribution of the error density (min size = 0.001m)
820
|}
821
822
<div id='img-8'></div>
823
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
824
|-
825
|[[Image:Draft_Samper_249588187-fig7.png|474px|]]
826
|[[Image:Draft_Samper_249588187-fig7_b.png|534px|Detail of the refined mesh shown in the  previous figure]]
827
|- style="text-align: center; font-size: 75%;"
828
| colspan="2" | '''Figure 8:''' Detail of the refined mesh shown in the  previous figure
829
|}
830
831
832
{|  class="floating_tableSCP wikitable" style="text-align: center; margin: 1em auto;min-width:50%;"
833
|+ style="font-size: 75%;" |Table. 2 Flow past a NACA profile (equidistribution of error density). Convergence of the global error with the mesh refinement (min. element size = 0.001m)
834
|- style="border-top: 2px solid;"
835
| style="border-left: 2px solid;border-right: 2px solid;" |  Mesh 
836
| style="border-left: 2px solid;border-right: 2px solid;" | Global error <math display="inline">\left({\Vert P\Vert _\Omega \over U}\right)</math>(%)  
837
| style="border-left: 2px solid;border-right: 2px solid;" | Nodes 
838
| style="border-left: 2px solid;border-right: 2px solid;" | Elements
839
| style="border-left: 2px solid;border-right: 2px solid;" | Drag error %
840
| style="border-left: 2px solid;border-right: 2px solid;" | Lift error% 
841
|- style="border-top: 2px solid;"
842
| style="border-left: 2px solid;border-right: 2px solid;" |  1
843
| style="border-left: 2px solid;border-right: 2px solid;" | 5.55
844
| style="border-left: 2px solid;border-right: 2px solid;" | 2574
845
| style="border-left: 2px solid;border-right: 2px solid;" | 4784
846
| style="border-left: 2px solid;border-right: 2px solid;" | 29.1810942 
847
| style="border-left: 2px solid;border-right: 2px solid;" | 201.915023
848
|- style="border-top: 2px solid;"
849
| style="border-left: 2px solid;border-right: 2px solid;" |  2
850
| style="border-left: 2px solid;border-right: 2px solid;" | 2.18
851
| style="border-left: 2px solid;border-right: 2px solid;" | 10594
852
| style="border-left: 2px solid;border-right: 2px solid;" | 21188
853
| style="border-left: 2px solid;border-right: 2px solid;" | 10.4796861
854
| style="border-left: 2px solid;border-right: 2px solid;" | 1.16034921
855
|- style="border-top: 2px solid;border-bottom: 2px solid;"
856
| style="border-left: 2px solid;border-right: 2px solid;" |  3
857
| style="border-left: 2px solid;border-right: 2px solid;" | 2.14
858
| style="border-left: 2px solid;border-right: 2px solid;" | 39895
859
| style="border-left: 2px solid;border-right: 2px solid;" | 79790
860
| style="border-left: 2px solid;border-right: 2px solid;" | 3,17339951
861
| style="border-left: 2px solid;border-right: 2px solid;" | 0.5988364
862
863
|}
864
865
<div id='img-9'></div>
866
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
867
|-
868
|[[Image:Draft_Samper_249588187-Fig8.png|600px|Mesh obtained after four refinement steps using the equidistribution of the error density (min size = 0.002m).]]
869
|- style="text-align: center; font-size: 75%;"
870
| colspan="1" | '''Figure 9:''' Mesh obtained after four refinement steps using the equidistribution of the error density (min size = 0.002m).
871
|}
872
873
<div id='img-10'></div>
874
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
875
|-
876
|[[Image:Draft_Samper_249588187-Fig9.png|480px|]]
877
|[[Image:Draft_Samper_249588187-Fig9_b.png|540px|Detail of the refined mesh of Figure 9.]]
878
|- style="text-align: center; font-size: 75%;"
879
| colspan="2" | '''Figure 10:''' Detail of the refined mesh of Figure 9.
880
|}
881
882
Table 2 shows that the criterium of equidistribution of error density leads to a smaller error for the drag and lift  in  fewer time steps. Note however the increase in the number of elements and nodes  with respect to that of Table 1.  In order to avoid this problem, the minimum element size was increased next to 0.002m. This allowed to complete the mesh refinement process and the results are shown in Table 3 and Figures 9 and 10. Note that the drag and lift errors have now slightly increased, showing the dependence of these values with the minimum element size. Again we note that these errors can be reduced using a turbulence model.
883
884
885
{|  class="floating_tableSCP wikitable" style="text-align: center; margin: 1em auto;min-width:50%;"
886
|+ style="font-size: 75%;" |Table. 3 Flow past a NACA airfoil (equidistribution of error density). Convergence of the global error with the mesh refinement. Minimum  element size = 0.002m
887
|- style="border-top: 2px solid;"
888
| style="border-left: 2px solid;border-right: 2px solid;" |  Mesh 
889
| style="border-left: 2px solid;border-right: 2px solid;" | Global error <math display="inline">\left({\Vert P\Vert _\Omega \over U}\right)</math>(%)  
890
| style="text-align: right;border-left: 2px solid;border-right: 2px solid;" | Nodes 
891
| style="text-align: right;border-left: 2px solid;border-right: 2px solid;" | Elements
892
| style="border-left: 2px solid;border-right: 2px solid;" | Drag error %
893
| style="border-left: 2px solid;border-right: 2px solid;" | Lift error% 
894
|- style="border-top: 2px solid;"
895
| style="border-left: 2px solid;border-right: 2px solid;" |  1
896
| style="border-left: 2px solid;border-right: 2px solid;" | 6.38
897
| style="text-align: right;border-left: 2px solid;border-right: 2px solid;" | 2574 
898
| style="text-align: right;border-left: 2px solid;border-right: 2px solid;" | 4784 
899
| style="border-left: 2px solid;border-right: 2px solid;" | 29.1810942
900
| style="border-left: 2px solid;border-right: 2px solid;" | 2.99902143
901
|- style="border-top: 2px solid;"
902
| style="border-left: 2px solid;border-right: 2px solid;" |  2 
903
| style="border-left: 2px solid;border-right: 2px solid;" | 3.42  
904
| style="text-align: right;border-left: 2px solid;border-right: 2px solid;" | 9871 
905
| style="text-align: right;border-left: 2px solid;border-right: 2px solid;" | 19742 
906
| style="border-left: 2px solid;border-right: 2px solid;" | 11.4092413 
907
| style="border-left: 2px solid;border-right: 2px solid;" | 1.21135183
908
|- style="border-top: 2px solid;"
909
| style="border-left: 2px solid;border-right: 2px solid;" |  3 
910
| style="border-left: 2px solid;border-right: 2px solid;" | 2.23  
911
| style="text-align: right;border-left: 2px solid;border-right: 2px solid;" | 29414 
912
| style="text-align: right;border-left: 2px solid;border-right: 2px solid;" | 58828 
913
| style="border-left: 2px solid;border-right: 2px solid;" | 6.34299015 
914
| style="border-left: 2px solid;border-right: 2px solid;" | 0.78983056
915
|- style="border-top: 2px solid;"
916
| style="border-left: 2px solid;border-right: 2px solid;" |  4 
917
| style="border-left: 2px solid;border-right: 2px solid;" | 2.11  
918
| style="text-align: right;border-left: 2px solid;border-right: 2px solid;" | 50062 
919
| style="text-align: right;border-left: 2px solid;border-right: 2px solid;" | 100124 
920
| style="border-left: 2px solid;border-right: 2px solid;" | 6.31490459 
921
| style="border-left: 2px solid;border-right: 2px solid;" | 0.78318179
922
|- style="border-top: 2px solid;border-bottom: 2px solid;"
923
| style="border-left: 2px solid;border-right: 2px solid;" |  5 
924
| style="border-left: 2px solid;border-right: 2px solid;" | 2.15  
925
| style="text-align: right;border-left: 2px solid;border-right: 2px solid;" | 71826 
926
| style="text-align: right;border-left: 2px solid;border-right: 2px solid;" | 143652 
927
| style="border-left: 2px solid;border-right: 2px solid;" | 6.76638966 
928
| style="border-left: 2px solid;border-right: 2px solid;" | 0.79634854
929
930
|}
931
932
Finally Figures 11 and 12 respectively show the evolution of the drag  and lift errors for the three cases considered and Figure 13 shows the convergence of the global error.
933
934
<div id='img-11'></div>
935
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
936
|-
937
|[[Image:Draft_Samper_249588187-graph1.png|600px|Flow past a NACA airfoil. Convergence of the drag error for both refinement methods]]
938
|- style="text-align: center; font-size: 75%;"
939
| colspan="1" | '''Figure 11:''' Flow past a NACA airfoil. Convergence of the drag error for both refinement methods
940
|}
941
942
<div id='img-12'></div>
943
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
944
|-
945
|[[Image:Draft_Samper_249588187-graph2.png|600px|Flow past a NACA airfoil. Convergence of the lift error for both refinement methods. Numbers next to legend denote the minimum element size prescribed]]
946
|- style="text-align: center; font-size: 75%;"
947
| colspan="1" | '''Figure 12:''' Flow past a NACA airfoil. Convergence of the lift error for both refinement methods. Numbers next to legend denote the minimum element size prescribed
948
|}
949
950
<div id='img-13'></div>
951
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
952
|-
953
|[[Image:Draft_Samper_249588187-graf_dof_naca.png|600px|Convergence of the global error for both refinement methods]]
954
|- style="text-align: center; font-size: 75%;"
955
| colspan="1" | '''Figure 13:''' Convergence of the global error for both refinement methods
956
|}
957
958
===7.2 <span id='lb-7.2'></span>Flow past a cylinder Re=100===
959
960
The mesh refinement process for this transient problem has been carried out at a time <math display="inline">t=60s</math>. Figure 14 shows the original mesh of 2135 nodes and 5365 linear elements used.
961
962
Using the equidistribution of the global error, the maximum and minimum size of the elements to be generated  are 0.3m and 0.001m, respectively.  Figures 15 and 16 show the mesh obtained after two refinement steps. Figures 17 and 18 show the velocity field for the third mesh.  Table 4 shows the value of the global error and the evolution of the number of elements and nodes during the mesh refinement process.
963
964
<div id='img-14'></div>
965
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
966
|-
967
|[[Image:Draft_Samper_249588187-fig8.png|600px|Flow past a cylinder. Radius=1m, Re=100. Original mesh]]
968
|- style="text-align: center; font-size: 75%;"
969
| colspan="1" | '''Figure 14:''' Flow past a cylinder. Radius=1m, <math>Re=100</math>. Original mesh
970
|}
971
972
<div id='img-15'></div>
973
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
974
|-
975
|[[Image:Draft_Samper_249588187-fig9.png|528px|Mesh obtained after two refinement steps using the equidistribution of the global error]]
976
|- style="text-align: center; font-size: 75%;"
977
| colspan="1" | '''Figure 15:''' Mesh obtained after two refinement steps using the equidistribution of the global error
978
|}
979
980
<div id='img-16'></div>
981
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
982
|-
983
|[[Image:Draft_Samper_249588187-fig10.png|600px|Detail of the refined mesh in the neighbourhood of the cylinder]]
984
|- style="text-align: center; font-size: 75%;"
985
| colspan="1" | '''Figure 16:''' Detail of the refined mesh in the neighbourhood of the cylinder
986
|}
987
988
<div id='img-17'></div>
989
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
990
|-
991
|[[Image:Draft_Samper_249588187-cilin_line_1.png|600px|Velocity field after two refinement steps. Equidistribution of global error]]
992
|- style="text-align: center; font-size: 75%;"
993
| colspan="1" | '''Figure 17:''' Velocity field after two refinement steps. Equidistribution of global error
994
|}
995
996
<div id='img-18'></div>
997
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
998
|-
999
|[[Image:Draft_Samper_249588187-cilin_line_2.png|486px|Detail of the velocity field]]
1000
|- style="text-align: center; font-size: 75%;"
1001
| colspan="1" | '''Figure 18:''' Detail of the velocity field
1002
|}
1003
1004
<br/>
1005
1006
1007
{|  class="floating_tableSCP wikitable" style="text-align: center; margin: 1em auto;min-width:50%;"
1008
|+ style="font-size: 75%;" |Table. 4 Flow past a cylinder (equidistribution of global error). Convergence of the global error with the mesh refinement
1009
|- style="border-top: 2px solid;"
1010
| style="border-left: 2px solid;border-right: 2px solid;" |  Mesh 
1011
| style="border-left: 2px solid;border-right: 2px solid;" | Global error <math display="inline">\left({\Vert P\Vert _\Omega \over U}\right)</math>(%) 
1012
| style="border-left: 2px solid;border-right: 2px solid;" | Nodes 
1013
| style="border-left: 2px solid;border-right: 2px solid;" | Elements
1014
|- style="border-top: 2px solid;"
1015
| style="border-left: 2px solid;border-right: 2px solid;" |  1
1016
| style="border-left: 2px solid;border-right: 2px solid;" | 8.72 
1017
| style="border-left: 2px solid;border-right: 2px solid;" | 2135
1018
| style="border-left: 2px solid;border-right: 2px solid;" | 5365
1019
|- style="border-top: 2px solid;"
1020
| style="border-left: 2px solid;border-right: 2px solid;" |  2
1021
| style="border-left: 2px solid;border-right: 2px solid;" | 2.67
1022
| style="border-left: 2px solid;border-right: 2px solid;" | 7863
1023
| style="border-left: 2px solid;border-right: 2px solid;" | 15875
1024
|- style="border-top: 2px solid;"
1025
| style="border-left: 2px solid;border-right: 2px solid;" |  3
1026
| style="border-left: 2px solid;border-right: 2px solid;" | 2.05
1027
| style="border-left: 2px solid;border-right: 2px solid;" | 12035
1028
| style="border-left: 2px solid;border-right: 2px solid;" | 24264
1029
|- style="border-top: 2px solid;"
1030
| style="border-left: 2px solid;border-right: 2px solid;" |  4
1031
| style="border-left: 2px solid;border-right: 2px solid;" | 2
1032
| style="border-left: 2px solid;border-right: 2px solid;" | 14448
1033
| style="border-left: 2px solid;border-right: 2px solid;" | 29114
1034
|- style="border-top: 2px solid;border-bottom: 2px solid;"
1035
| style="border-left: 2px solid;border-right: 2px solid;" |  5
1036
| style="border-left: 2px solid;border-right: 2px solid;" | 1.86
1037
| style="border-left: 2px solid;border-right: 2px solid;" | 16661
1038
| style="border-left: 2px solid;border-right: 2px solid;" | 33539
1039
1040
|}
1041
1042
The same problem was solved using the criterium based on the equidistribution of the error density. Now the minimum size of the elements generated is  increased to 0.003m, in order to avoid problems in the regeneration of the  mesh.  The mesh obtained after two  steps is shown in Figures 19 and 20. Table 5 shows the evolution of the global error and the number of elements and nodes during the mesh refinement process.
1043
1044
Figure 21 shows the convergence  of the global error for both cases.
1045
1046
<div id='img-19'></div>
1047
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
1048
|-
1049
|[[Image:Draft_Samper_249588187-fig13.png|600px|Mesh obtained after two refinement steps using the equidistribution of the error density]]
1050
|- style="text-align: center; font-size: 75%;"
1051
| colspan="1" | '''Figure 19:''' Mesh obtained after two refinement steps using the equidistribution of the error density
1052
|}
1053
1054
<div id='img-20'></div>
1055
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
1056
|-
1057
|[[Image:Draft_Samper_249588187-fig14.png|600px|Details of the refined mesh of Figure 19]]
1058
|- style="text-align: center; font-size: 75%;"
1059
| colspan="1" | '''Figure 20:''' Details of the refined mesh of Figure 19
1060
|}
1061
1062
1063
{|  class="floating_tableSCP wikitable" style="text-align: center; margin: 1em auto;min-width:50%;"
1064
|+ style="font-size: 75%;" |Table. 5 Flow past a cylinder (equidistribution of error density). Convergence of the global error with the mesh refinement
1065
|- style="border-top: 2px solid;"
1066
| style="border-left: 2px solid;border-right: 2px solid;" |  Mesh 
1067
| style="border-left: 2px solid;border-right: 2px solid;" | Global error <math display="inline">\left({\Vert P\Vert _\Omega \over U}\right)</math>(%) 
1068
| style="border-left: 2px solid;border-right: 2px solid;" | Nodes 
1069
| style="border-left: 2px solid;border-right: 2px solid;" | Elements
1070
|- style="border-top: 2px solid;"
1071
| style="border-left: 2px solid;border-right: 2px solid;" |  1
1072
| style="border-left: 2px solid;border-right: 2px solid;" | 8.6
1073
| style="border-left: 2px solid;border-right: 2px solid;" | 2653
1074
| style="border-left: 2px solid;border-right: 2px solid;" | 5092
1075
|- style="border-top: 2px solid;"
1076
| style="border-left: 2px solid;border-right: 2px solid;" |  2
1077
| style="border-left: 2px solid;border-right: 2px solid;" | 2.96
1078
| style="border-left: 2px solid;border-right: 2px solid;" | 9198
1079
| style="border-left: 2px solid;border-right: 2px solid;" | 18057
1080
|- style="border-top: 2px solid;border-bottom: 2px solid;"
1081
| style="border-left: 2px solid;border-right: 2px solid;" |  3
1082
| style="border-left: 2px solid;border-right: 2px solid;" | 2.37
1083
| style="border-left: 2px solid;border-right: 2px solid;" | 77064
1084
| style="border-left: 2px solid;border-right: 2px solid;" | 153190
1085
1086
|}
1087
1088
<div id='img-21'></div>
1089
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
1090
|-
1091
|[[Image:Draft_Samper_249588187-graf_dof_cilin.png|600px|The convergence of the global error for both cases.]]
1092
|- style="text-align: center; font-size: 75%;"
1093
| colspan="1" | '''Figure 21:''' The convergence of the global error for both cases.
1094
|}
1095
1096
===7.3 Driven cavity flow===
1097
1098
This is another well known benchmark problem in the CFD literature. The problem is solved for a Reynolds number of 1000. The mesh regeneration was performed at the time <math display="inline">t=2s</math>.
1099
1100
The initial mesh of  477 nodes and 952 linear elements is shown in Figure 22. The mesh obtained after three refinement steps  and the corresponding velocity  field  are shown in Figure 23. The mesh adaptation criterium based on the equidistribution of the global error was used.  Table 6 shows the value of the global error and the evolution of the number of elements and nodes during the mesh refinement process.
1101
1102
Figure 24 shows the convergence of the global error.
1103
1104
<div id='img-22'></div>
1105
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
1106
|-
1107
|[[Image:Draft_Samper_249588187-fig20.png|348px|Driven cavity flow. Side length =1m, Re=1000. Original mesh]]
1108
|- style="text-align: center; font-size: 75%;"
1109
| colspan="1" | '''Figure 22:''' Driven cavity flow. Side length =1m, <math>Re=1000</math>. Original mesh
1110
|}
1111
1112
<div id='img-23'></div>
1113
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
1114
|-
1115
|[[Image:Draft_Samper_249588187-fig21.png|324px|]]
1116
|[[Image:Draft_Samper_249588187-cavity_line_4.png|318px|Mesh obtained after three refinement steps (equidistribution of the global error) and velocity field.]]
1117
|- style="text-align: center; font-size: 75%;"
1118
| colspan="2" | '''Figure 23:''' Mesh obtained after three refinement steps (equidistribution of the global error) and velocity field.
1119
|}
1120
1121
1122
{|  class="floating_tableSCP wikitable" style="text-align: center; margin: 1em auto;min-width:50%;"
1123
|+ style="font-size: 75%;" |Table. 6 Driven cavity flow (equidistribution of global error). Convergence of the global error with the mesh refinement
1124
|- style="border-top: 2px solid;"
1125
| style="border-left: 2px solid;border-right: 2px solid;" |  Mesh 
1126
| style="border-left: 2px solid;border-right: 2px solid;" | Global error <math display="inline">\left({\Vert P\Vert _\Omega \over U}\right)</math>(%) 
1127
| style="border-left: 2px solid;border-right: 2px solid;" | Nodes 
1128
| style="border-left: 2px solid;border-right: 2px solid;" | Elements
1129
|- style="border-top: 2px solid;"
1130
| style="border-left: 2px solid;border-right: 2px solid;" |  1
1131
| style="border-left: 2px solid;border-right: 2px solid;" | 26.88
1132
| style="border-left: 2px solid;border-right: 2px solid;" | 477
1133
| style="border-left: 2px solid;border-right: 2px solid;" | 872
1134
|- style="border-top: 2px solid;"
1135
| style="border-left: 2px solid;border-right: 2px solid;" |  2
1136
| style="border-left: 2px solid;border-right: 2px solid;" | 14.5
1137
| style="border-left: 2px solid;border-right: 2px solid;" | 3089
1138
| style="border-left: 2px solid;border-right: 2px solid;" | 5985
1139
|- style="border-top: 2px solid;"
1140
| style="border-left: 2px solid;border-right: 2px solid;" |  3
1141
| style="border-left: 2px solid;border-right: 2px solid;" | 9.56
1142
| style="border-left: 2px solid;border-right: 2px solid;" | 6475
1143
| style="border-left: 2px solid;border-right: 2px solid;" | 12521
1144
|- style="border-top: 2px solid;"
1145
| style="border-left: 2px solid;border-right: 2px solid;" |  4
1146
| style="border-left: 2px solid;border-right: 2px solid;" | 9.42
1147
| style="border-left: 2px solid;border-right: 2px solid;" | 8103
1148
| style="border-left: 2px solid;border-right: 2px solid;" | 15700
1149
|- style="border-top: 2px solid;"
1150
| style="border-left: 2px solid;border-right: 2px solid;" |  5
1151
| style="border-left: 2px solid;border-right: 2px solid;" | 9.14
1152
| style="border-left: 2px solid;border-right: 2px solid;" | 9411
1153
| style="border-left: 2px solid;border-right: 2px solid;" | 18280
1154
|- style="border-top: 2px solid;border-bottom: 2px solid;"
1155
| style="border-left: 2px solid;border-right: 2px solid;" |  6
1156
| style="border-left: 2px solid;border-right: 2px solid;" | 6.63
1157
| style="border-left: 2px solid;border-right: 2px solid;" | 12997
1158
| style="border-left: 2px solid;border-right: 2px solid;" | 25207
1159
1160
|}
1161
1162
<div id='img-24'></div>
1163
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
1164
|-
1165
|[[Image:Draft_Samper_249588187-graf_dof_cavity.png|600px|Convergence of the global error in the cavity problem.]]
1166
|- style="text-align: center; font-size: 75%;"
1167
| colspan="1" | '''Figure 24:''' Convergence of the global error in the cavity problem.
1168
|}
1169
1170
<div id='img-25'></div>
1171
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
1172
|-
1173
|[[Image:Draft_Samper_249588187-fig23.png|600px|High-lift airfoil configuration. Characteristic length =0.6m, Re=100. Initial mesh]]
1174
|- style="text-align: center; font-size: 75%;"
1175
| colspan="1" | '''Figure 25:''' High-lift airfoil configuration. Characteristic length =0.6m, Re=100. Initial mesh
1176
|}
1177
1178
<div id='img-26'></div>
1179
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
1180
|-
1181
|[[Image:Draft_Samper_249588187-fig24.png|480px|]]
1182
|[[Image:Draft_Samper_249588187-fig24.png|468px|Mesh obtained after two refinement steps using equiditribution of global error (μ=0.01)]]
1183
|- style="text-align: center; font-size: 75%;"
1184
| colspan="2" | '''Figure 26:''' Mesh obtained after two refinement steps using equiditribution of global error (<math>\mu =0.01</math>)
1185
|}
1186
1187
<div id='img-27'></div>
1188
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
1189
|-
1190
|[[Image:Draft_Samper_249588187-fig25.png|480px|]]
1191
|[[Image:Draft_Samper_249588187-fig25.png|480px|Mesh obtained in the fith after changing μ=0.01 to μ= 0.005 in the fourth step]]
1192
|- style="text-align: center; font-size: 75%;"
1193
| colspan="2" | '''Figure 27:''' Mesh obtained in the fith after changing <math>\mu =0.01</math> to <math>\mu = 0.005</math> in the fourth step
1194
|}
1195
1196
===7.4 Flow past a  wing in high-lift configuration===
1197
1198
This last example is a multiple  airfoil configuration including slat, wing and flap. The problem has been solved with a laminar solver and a Reynolds number  of 100.  In this example the maximum size of the  generated elements has been reduced to 0.05m. The initial mesh of 11703 nodes and 23480 elements is shown in Figure 25.   The resulting mesh after two refinement steps is shown in Figure 26. The convergence of the global error with the mesh refinement is shown in Table 7.  The algorithm based on the equidistribution of the global error has been used  in this example.
1199
1200
1201
{|  class="floating_tableSCP wikitable" style="text-align: center; margin: 1em auto;min-width:50%;"
1202
|+ style="font-size: 75%;" |Table. 7 High-lift configuration (equidistribution of global error). Convergence of the global error with the mesh refinement
1203
|- style="border-top: 2px solid;"
1204
| style="border-left: 2px solid;border-right: 2px solid;" |  Mesh 
1205
| style="border-left: 2px solid;border-right: 2px solid;" | Global error <math display="inline">\left({\Vert P\Vert _\Omega \over U}\right)</math>(%) 
1206
| style="border-left: 2px solid;border-right: 2px solid;" | Nodes 
1207
| style="text-align: right;border-left: 2px solid;border-right: 2px solid;" | Elements
1208
|- style="border-top: 2px solid;"
1209
| style="border-left: 2px solid;border-right: 2px solid;" |  1
1210
| style="border-left: 2px solid;border-right: 2px solid;" | 3.76
1211
| style="border-left: 2px solid;border-right: 2px solid;" | 11703
1212
| style="text-align: right;border-left: 2px solid;border-right: 2px solid;" | 22874
1213
|- style="border-top: 2px solid;"
1214
| style="border-left: 2px solid;border-right: 2px solid;" |  2
1215
| style="border-left: 2px solid;border-right: 2px solid;" | 2.51
1216
| style="border-left: 2px solid;border-right: 2px solid;" | 18334
1217
| style="text-align: right;border-left: 2px solid;border-right: 2px solid;" | 35909
1218
|- style="border-top: 2px solid;"
1219
| style="border-left: 2px solid;border-right: 2px solid;" |  3
1220
| style="border-left: 2px solid;border-right: 2px solid;" | 1.95
1221
| style="border-left: 2px solid;border-right: 2px solid;" | 26718
1222
| style="text-align: right;border-left: 2px solid;border-right: 2px solid;" | 52487
1223
|- style="border-top: 2px solid;"
1224
| style="border-left: 2px solid;border-right: 2px solid;" |  4
1225
| style="border-left: 2px solid;border-right: 2px solid;" | 1.79
1226
| style="border-left: 2px solid;border-right: 2px solid;" | 27845
1227
| style="text-align: right;border-left: 2px solid;border-right: 2px solid;" | 54810
1228
|- style="border-top: 2px solid;border-bottom: 2px solid;"
1229
| style="border-left: 2px solid;border-right: 2px solid;" |  5
1230
| style="border-left: 2px solid;border-right: 2px solid;" | 1.58
1231
| style="border-left: 2px solid;border-right: 2px solid;" | 51361
1232
| style="text-align: right;border-left: 2px solid;border-right: 2px solid;" | 101532
1233
1234
|}
1235
1236
After the fourth step the prescribed global error was reduced from <math display="inline">\mu = 0.01</math> to <math display="inline">\mu = 0.005</math>. The result is shown in Figure 27.
1237
1238
The velocity field for the fifth  step is shown in Figure 28.
1239
1240
Figures 29 and 30 show the convergence of the drag and lift error, respectively. Figure 31 shows the evolution of the global error.
1241
1242
<div id='img-28'></div>
1243
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
1244
|-
1245
|[[Image:Draft_Samper_249588187-dasault_line.png|600px|Velocity field obtained after five refinement steps]]
1246
|- style="text-align: center; font-size: 75%;"
1247
| colspan="1" | '''Figure 28:''' Velocity field obtained after five refinement steps
1248
|}
1249
1250
<div id='img-29'></div>
1251
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
1252
|-
1253
|[[Image:Draft_Samper_249588187-graph3.png|600px|The convergence of the drag error]]
1254
|- style="text-align: center; font-size: 75%;"
1255
| colspan="1" | '''Figure 29:''' The convergence of the drag error
1256
|}
1257
1258
<div id='img-30'></div>
1259
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
1260
|-
1261
|[[Image:Draft_Samper_249588187-graph4.png|600px|The convergence of the lift error]]
1262
|- style="text-align: center; font-size: 75%;"
1263
| colspan="1" | '''Figure 30:''' The convergence of the lift error
1264
|}
1265
1266
<div id='img-31'></div>
1267
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
1268
|-
1269
|[[Image:Draft_Samper_249588187-graf_dof_dassault.png|600px|The convergence of the global error ]]
1270
|- style="text-align: center; font-size: 75%;"
1271
| colspan="1" | '''Figure 31:''' The convergence of the global error 
1272
|}
1273
1274
==8 CONCLUDING REMARKS==
1275
1276
An error estimator for the finite element analysis of  incompressible fluid flow problems has been proposed. The estimator is based on the weigthed form (the power) of the residuals of the stabilized momentum and mass balance equations. The particular form of the weighting functions used has been found from a finite calculus formulation of the Navier-Stokes equations. This weighted residual form derived is also applicable for error estimation and mesh adaptivity in other problems in computational mechanics <span id='citeF-22'></span>[[#cite-22|22]].
1277
1278
A simple mesh refinement strategy based on the equal distribution of the total residual power among the elements in the mesh has been presented. The refinement strategy has proven to be more cost-effective than that based on the equal distribution of the density of the error also presented in the paper. The efficiency of the  error estimation and the mesh adaptation process has been verified in the solution of steady state and transient flow problems  using a FIC stabilized formulation for the Navier-Stokes equations and a standard fractional step scheme.
1279
1280
The method yields refined meshes which capture very well the features of the velocity and pressure fields in high gradient zones near the solid boundaries and within the analysis domain. Note  that  in the problems presented here the global error is not reduced below its limited value. The reason is that a minimun  element size has been prescribed in order to facilitate the mesh generation process. Should this limit be reduced then the global error  would be subsequently reduced.
1281
1282
Further work in this field will focus in the relationship of the residual power error values and the error in selected individual variables (velocity, pressure, etc) or other solutions outputs target (drag, lift, etc.).
1283
1284
==ACKNOWLEDGEMENTS==
1285
1286
The numerical results were obtained using the finite element code Tdyn provided by the company COMPASS Ingeniería y Sistemas SA <span id='citeF-31'></span>[[#cite-31|31]]. Thanks are also given to Dassault Aviation for providing the multiple airfoil geometry for the last example.
1287
1288
===BIBLIOGRAPHY===
1289
1290
<div id="cite-1"></div>
1291
'''[[#citeF-1|[1]]]''' O.C. Zienkiewicz and  R.L. Taylor. ''The finite element method''. 5th Edition, 3 Volumes, Butterworth&#8211;Heinemann, 2000.
1292
1293
<div id="cite-2"></div>
1294
'''[[#citeF-2|[2]]]''' O.C. Zienkiewicz and J.Z. Zhu. A simple error estimator and  adaptive procedure for practical engineering analysis. ''Int. J. Num.  Meth. Eng.'', '''24''', 337-57, 1987.
1295
1296
<div id="cite-3"></div>
1297
'''[3]''' J.Z. Zhu and O.C. Zienkiewicz. Adaptive techniques in the finite  element method. ''Comm. App. Num. Math.'', '''4''', 197&#8211;204, 1988.
1298
1299
<div id="cite-4"></div>
1300
'''[[#citeF-4|[4]]]'''  Oñate E. and G. Bugeda. A study of mesh optimality criteria in adaptive finite element analysis. ''Engineering Computations'', Vol. '''3''', 307-321, 1994.
1301
1302
<div id="cite-5"></div>
1303
'''[[#citeF-5|[5]]]'''  J. Wu, Z.J. Zhu, J. Szmelter and O.C. Zienkiewicz. Error estimation and adaptivity in Navier-Stokes incompressible flows. ''Comput Mechancis'', '''6''', 254&#8211;70, 1990.
1304
1305
<div id="cite-6"></div>
1306
'''[[#citeF-6|[6]]]''' R. Löhner, K. Morgan and O.C. Zienkiewicz. An adaptive  finite element procedure for compressible high speed flows. ''Comp. Meth. Appl. Mech. Eng.'', '''51''', 441&#8211;65, 1985.
1307
1308
<div id="cite-7"></div>
1309
'''[7]''' R. Löhner, K Morgan and O.C. Zienkiewicz. Adaptive grid refinement for the Euler and compressible Navier-Stokes equations. in  '' Accuracy Estimates and Adaptive Refinements in Finite Element Computations''  (eds I. Babuska, O.C. Zienkiewicz, J. Gago y  E.R.  de A. Oliveira), cap. 15, pp. 281&#8211;98, Wiley, Chichester, 1986.
1310
1311
<div id="cite-8"></div>
1312
'''[8]''' J.  Peraire, M. Vahdati, K. Morgan and O.C. Zienkiewicz.   Adaptive remeshing for compressible flow computations. ''J. Comp. Phys.'', '''72''', 449&#8211;66, 1987.
1313
1314
<div id="cite-9"></div>
1315
'''[9]''' J.T. Oden, T. Strouboulis and P. Devloo. Adaptive finite element methods for high speed compressible flows. ''Int. J. Num. Meth. Eng.'', '''7''', 1211&#8211;28, 1987.
1316
1317
<div id="cite-10"></div>
1318
'''[10]''' R. Löhner. Adaptive remeshing for transient problems. '' Comput. Meth. Appl. Mech. Eng.'', '''75''', 195&#8211;214, 1989.
1319
1320
<div id="cite-11"></div>
1321
'''[[#citeF-11|[11]]]'''  J.T. Oden, W. Wu and M. Ainsworth. An a posteriori error estimate for finite element approximations of the Navier-Stokes equations, '' Comput. Meth. Appl. Mech. Eng.'', '''111''', 185&#8211;202, 1993.
1322
1323
<div id="cite-12"></div>
1324
'''[12]'''  C. Johnson, R. Rannacher and M. Boman. Numerics and hydrodynamic stability: Towards error control in CFD. ''SIAM J. Numer. Anal.'', '''32''', 1058&#8211;1079, 1995.
1325
1326
<div id="cite-13"></div>
1327
'''[13]'''  R.C. Almeida, R.A. Feijoo, A.C. Galeao, C. Padra and R.S. Silva. Adaptive finite element computational fluid dynamics using an anisotropic error estimator. ''Comput. Meth. Appl. Mech. Eng.'', '''182''', 379&#8211;400, 2000.
1328
1329
<div id="cite-14"></div>
1330
'''[14]'''  J.M. Cascón, G.C. García and R. Rodríguez. A priori and a posteriori error analysis for a large-scale ocean circulation finite element model. ''Comput. Meth. Appl. Mech. Eng.'', '''192''', 5305&#8211;5327, 2003.
1331
1332
<div id="cite-15"></div>
1333
'''[[#citeF-15|[15]]]'''  W. Bangerth and R. Rannacher. ''Adaptive finite elements for differential equations''. Birkhäuser Verlag, Basel, 2003.
1334
1335
<div id="cite-16"></div>
1336
'''[16]'''  D.L. Darmofal, D.A. Venditti. Anisotropic grid adaptation for functional outputs: Application to two-dimensional viscous flows. ''Journal of Computational Physics'', '''187''', 22&#8211;66, 2003.
1337
1338
<div id="cite-17"></div>
1339
'''[[#citeF-17|[17]]]''' E. Oñate. Derivation of stabilized equations for  advective-diffusive transport and fluid flow problems.  ''Comput. Methods Appl. Mech. Engrg.'', '''151:1-2''', 233&#8211;267, 1998.
1340
1341
<div id="cite-18"></div>
1342
'''[[#citeF-18|[18]]]'''  E. Oñate. A stabilized finite element method for incompressible viscous flows using a finite increment calculus formulation. ''Comput. Methods Appl. Mech. Engrg.'', '''182:1&#8211;2''', 355&#8211;370,  2000.
1343
1344
<div id="cite-19"></div>
1345
'''[[#citeF-19|[19]]]'''  E. Oñate and  J. García. A finite element method for fluid-structure interaction with surface waves using a finite calculus formulation. in ''Comput. Methods Appl. Mech. Engrg.'', '''191:6-7''', 635-660, 2001.
1346
1347
<div id="cite-20"></div>
1348
'''[20]'''  E. Oñate,  J. García, G. Bugeda and S.R. Idelsohn. A general stabilized formulation for incompressible fluid flow using finite calculus and the FEM. In ''Towards a new fluid dynamics with its challenges in aeronautics'', J. Periaux, D. Chamption, O. Pironneau and Ph. Thomas (Eds.), CIMNE, Barcelona, Spain 2002.
1349
1350
<div id="cite-21"></div>
1351
'''[21]'''   E. Oñate, J. García, S.R. Idelsohn and F. Del Pin. FIC formulations for finite element analysis of incompressible flows. Eulerian, ALE and Lagrangian approaches, Submitted to ''Comput. Methods Appl. Mech. Engng.'', July 2004.
1352
1353
<div id="cite-22"></div>
1354
'''[[#citeF-22|[22]]]'''   E. Oñate, J. Arteaga and R. Flores. A weigthed recovered residual method for error estimation and mesh adaptivity in convection-diffusion problems. Research Report CIMNE, Barcelona, September 2004.
1355
1356
<div id="cite-23"></div>
1357
'''[[#citeF-23|[23]]]'''  GiD. The personal pre/postprocessing system. CIMNE, Barcelona, www.gidhome.com, 2004.
1358
1359
<div id="cite-24"></div>
1360
'''[24]''' M.A. Cruchaga and E. Oñate. A finite element formulation for incompressible flow problems using a generalized streamline operator. ''Comput. Methods Appl. Mech. Engng.'', '''143''', 49&#8211;67, 1997.
1361
1362
<div id="cite-25"></div>
1363
'''[25]''' T.J.R. Hughes and M. Mallet. A new finite element formulations  for computational fluid dynamics: III. The generalized streamline operator  for multidimensional advective-diffusive systems. ''Comput Methods Appl.  Mech. Engrg.'', '''58''', 305&#8211;328, (1986a).
1364
1365
<div id="cite-26"></div>
1366
'''[26]''' T.J.R. Hughes, L.P. Franca and G.M. Hulbert. A new finite  element formulation for computational fluid dynamics: VIII. The  Galerkin/least-squares method for advective-diffusive equations. '' Comput. Methods Appl. Mech. Engrg.'', '''73''', 173&#8211;189, (1989).
1367
1368
<div id="cite-27"></div>
1369
'''[27]''' T.J.R. Hughes. Multiscale phenomena: Green functions,  subgrid scale models, bubbles and the origins of stabilized methods.  ''Comput. Methods Appl. Mech. Engrg'', '''127''', 387&#8211;401, (1995).
1370
1371
<div id="cite-28"></div>
1372
'''[28]'''  R. Codina. Comparison of some finite element methods for solving the diffusion-convection-reaction equation. ''Comput. Methods Appl. Mech. Engrg.'', '''156''', 185&#8211;210, (1998).
1373
1374
<div id="cite-29"></div>
1375
'''[29]'''  R. Codina. Stabilization of incompressibility and convection through orthogonal sub-scales in finite element method. ''Comput. Methods Appl. Mech. Engrg.'', '''190''', 1579&#8211;1599, (2000).
1376
1377
<div id="cite-30"></div>
1378
'''[30]''' O.C. Zienkiewicz and  R. Codina. A general algorithm for  compressible and incompressible flow. Part I: The split characteristic  based scheme. ''Int. J. Num. Meth. in Fluids'', '''20''', 869-85, (1995).
1379
1380
<div id="cite-31"></div>
1381
'''[[#citeF-31|[31]]]'''  Tdyn. An unstructured finite element code for fluid dynamic analysis, COMPASS Ingeniería y Sistemas SA, www.compassis.com, 2003.
1382

Return to Onate et al 2006a.

Back to Top

Document information

Published on 01/01/2006

DOI: 10.1016/j.cma.2004.07.054
Licence: CC BY-NC-SA license

Document Score

0

Times cited: 11
Views 44
Recommendations 0

Share this document

claim authorship

Are you one of the authors of this document?