Sans titre

5 downloads 0 Views 711KB Size Report
The flow of a viscous incompressible fluid depends on the temperature where this temperature presents high variations, which occurs only in a part of the ...
Automatic coupling and finite element discretization

The flow of a viscous incompressible fluid depends on the temperature where this temperature presents high variations, which occurs only in a part of the domain. In this case, it seems reasonable to couple the Navier-Stokes and heat equations only in this subdomain. A posteriori analysis leads to a numerical determination of the subdomain where all equations must be coupled. M. Braack, A. Ern

of the Navier–Stokes and heat equations • The full model

Christine Bernardi

• The simplified model • The discrete problem

Laboratoire Jacques-Louis Lions C.N.R.S. & Universit´ e Pierre et Marie Curie

• A posteriori analysis

with Fr´ ed´ eric Coquel and Pierre-Arnaud Raviart

• A strategy for the automatic coupling • About a turbulence model

Ω : connected bounded open set in IRd, d = 2 or 3, with a Lipschitz– continuous boundary.

 % &   −div ν(T ) ∇u + (u · ∇)u + grad p =           div u = 0     

f

in Ω,

−α ∆T + (u · ∇)T = g         u=0          T =T 0

The full model

in Ω,

in Ω, on ∂Ω, on ∂Ω.

Unknowns : velocity u, pressure p, temperature T . This model seems realistic in a large number of sitiuations.

Variational formulation 1 Find (u, p, T ) in H01(Ω)d × L2 0 (Ω) × H (Ω) such that

T = T0 Data : f , g and the boundary temperature T0. The parameter α is a positive constant. Here, we work without restriction with a viscosity which is truncated ∀ξ ∈ IR,

and regularized : ν ∈ W 2,∞(IR).

on ∂Ω,

and ∀v ∈ H01(Ω)d,

'



ν(T ) ∇ u : ∇ v dx +

'



(u · ∇)u · v(x) dx −

ν1 ≤ ν(ξ) ≤ ν2,

∀q ∈ L2 0 (Ω), ∀S ∈ H01(Ω),

α

'





'



'



(div v)(x)p(x) dx = (f , v),

(div u)(x)q(x) dx = 0,

grad T · grad S dx +

'



(u · ∇)T · S(x) dx = (g, S).

It is readily checked that this problem is equivalent to the previous system.

Existence result

Uniqueness of the solution if regular enough for sufficiently small data.

Lifting of the trace T0 by a simple version of the Hopf lemma, plus two key theorems : • the fixed point of Brouwer, • the dominated convergence of Lebesgue.

A regularity property

1 H 2 (∂Ω),

Theorem. For any data f in H −1(Ω)d, g in H −1(Ω) and T0 in 1 the problem has a solution (u, p, T ) in H01(Ω)d × L2 0 (Ω) × H (Ω). Moreover, this solution satisfies ν1 *u*H 1(Ω)d ≤ c *f *H −1(Ω)d ,

(

We assume for simplicity that Ω is a polygon or a polyhedron. Proposition. There exist a real number metry of Ω and on the ratio ν2/ν1 and ding on the geometry of Ω such that, 2 ≤ q + ≤ q1, and for any data (f , g) in 2− 1 ,q +

% & ν 1 *p*L2(Ω) ≤ c (1 + 2 ) *f *H −1(Ω)d + 2 *f *2 , H −1 (Ω)d ν1 ν1

α *T *H 1(Ω) ≤ c *f *H −1(Ω)d + *g*H −1(Ω) + α *T0*

N.G. Meyers

1

H 2 (∂Ω)

)

q0 > 2 depending on the geoa real number q1 only depenfor any q, 2 < q ≤ q0, and q +, + the space W −1,q (Ω)d × Lq (Ω)

and T0 in W q+ (∂Ω), any solution (u, p, T ) of the problem belongs to + W 1,q (Ω)d × Lq (Ω) × W 2,q (Ω). .

As well-known, q1 is ≥ 4 3 for a general domain Ω and ≥ 2 when Ω is convex.

We introduce an arbitrary decomposition of the domain without overlap Ω = Ωs ∪ Ωf ,

Ωs ∪ Ωf = ∅.

The key idea is as follows : We have decided to replace ν(T ) by a positive constant ν0 on the part Ωs of the domain.

The simplified model

∀ξ ∈ IR,

ν ∗(x, ξ) =

 ν(ξ) ν

for a.e. x in Ωf , for a.e. x in Ωs,

0

Note that the function : ξ 0→ ν ∗(·, ξ) satisfies exactly the same properties as the function ν(·). For instance, min{ν0, ν1} ≤ ν ∗(x, ξ) ≤ ν2

∀ξ ∈ IR,

for a.e. x in Ω.

Variational formulation 1 Find (u∗, p∗, T ∗) in H01(Ω)d × L2 0 (Ω) × H (Ω) such that

T = T0  % &   −div ν ∗(·, T ∗) ∇u∗ + (u∗ · ∇)u∗ + grad p∗ =           div u∗ = 0      −α ∆T ∗ + (u∗ · ∇)T ∗ = g          u∗ = 0        ∗  T = T0

f

in Ω, in Ω,

on ∂Ω,

and ∀v ∈ H01(Ω)d,

'



ν ∗(x, T ∗) ∇ u∗ : ∇ v dx +

'





in Ω, ∀q ∈ L2 0 (Ω),

on ∂Ω, on ∂Ω.

(u∗ · ∇)u∗ · v(x) dx

∀S ∈ H01(Ω),

α

'





'



'



(div v)(x)p∗(x) dx = (f , v),

(div u∗)(x)q(x) dx = 0,

grad T ∗ · grad S dx +

'



(u∗ · ∇)T ∗ · S(x) dx = (g, S).

The properties which are stated below are derived thanks to exactly the same arguments as for the full problem.

A regularity property N.G. Meyers Existence result We assume for simplicity that Ω is a polygon or a polyhedron.

1

Theorem. For any data f in H −1(Ω)d, g in H −1(Ω) and T0 in H 2 (∂Ω), the 1 problem has a solution (u∗, p∗, T ∗) in H01(Ω)d × L2 0 (Ω) × H (Ω). Moreover, this solution satisfies min{ν0, ν1} *u*H 1(Ω)d ≤ c *f *H −1(Ω)d , % & ν2 1 *p*L2(Ω) ≤ c (1 + ) *f *H −1(Ω)d + *f *2 , H −1 (Ω)d min{ν0, ν1} min{ν0, ν1}2 (

α *T *H 1(Ω) ≤ c *f *H −1(Ω)d + *g*H −1(Ω) + α *T0*

1

H 2 (∂Ω)

Proposition. There exist a real number q0∗ > 2 depending on the geometry of Ω and on the ratio ν2/ min{ν0, ν1} that, for the real number q1introduced above, for any q, 2 < q ≤ q0∗ , and q +, 2 ≤ q + ≤ q1, and for any 2− 1 ,q +

+

)

data (f , g) in the space W −1,q (Ω)d ×Lq (Ω) and T0 in W q+ (∂Ω), any solu+ tion (u, p, T ) of problem (2.4)−(2.5) belongs to W 1,q (Ω)d ×Lq (Ω)×W 2,q (Ω). .

Note however that the real number q0∗ tends to 2 whence the ratio ν2/ min{ν0, ν1} tends to +∞.

A consistency result Proposition. Let (Ωn)n be an increasing sequence of domains with Lipschitz-continuous boundaries converging to Ω. If (un, pn, T n) denotes a solution of the simplified problem with Ωf = Ωn and Ωs = Ω \ Ωf , +

+

+

there exists a subsequence (un , pn , T n )n+ which converges to a solution (u, p, T ) of the full problem.

The discrete problem

This result is of course not sufficient to evaluate the distance *u − u∗*H 1(Ω)d + *p − p∗*L2(Ω) + *T − T ∗*H 1(Ω). Only a posteriori analysis wiil provide a more precise but numerical result.

Find (uh, ph, Th) in Xh × Mh × Yh such that (Th)h : regular family of triangulations of Ω such that each element of Th is contained either in Ωs or in Ωf . Taylor–Hood finite elements *

∀vh ∈ Xh,

+

Xh = vh ∈ H01(Ω)d; ∀K ∈ Th, vh|K ∈ P2(K)d ,

Yh =

*

Th = T0h '



ν ∗(x, Th) ∇ uh : ∇ vh dx +

+ Sh ∈ H 1(Ω); ∀K ∈ Th, Sh|K ∈ P2(K) ,

Y0h =

Yh ∩ H01(Ω).

We also define T0h as the Lagrange interpolate of T0 with values in the trace space of Yh.

'



(uh · ∇)uh · vh(x) dx −

+

Mh = qh ∈ H 1(Ω) ∩ L20(Ω); ∀K ∈ Th, qh|K ∈ P1(K) .

*

on ∂Ω,

and

∀qh ∈ Mh, ∀Sh ∈ Y0 h,

α

'





'



'



(div vh)(x)ph(x) dx = (f , vh),

(div uh)(x)qh(x) dx = 0,

grad Th · grad Sh dx +

'



(uh · ∇)Th · Sh(x) dx = (g, Sh).

The a priori analysis of this problem makes use of the discrete implicit function theorem due to F. Brezzi, J. Rappaz, P.-A. Raviart. Thus, it requires a different formulation of both reduced and discrete problems.

Let S(ξ) denote the Stokes operator which, with a distribution F in H −1(Ω)d, associates the part u of the solution of the Stokes problem  % & ∗    −div ν (·, ξ) ∇u + grad p =

div u = 0    u=0

F

in Ω, in Ω, on ∂Ω.

∀vh ∈ Xh,

Let L denote the (inverse) Laplace operator which, with a distribution 1

G in H −1(Ω and a function T0 in H 2 (∂Ω), associates the solution T of the Laplace equation ,

−α ∆T = G T = T0

We also use the notation U ∗ = written in an abridged way F∗(U ∗) = U ∗+

-

.

S(T ∗) 0 0 L

in Ω, on ∂Ω,

-

. u∗ . Thus the reduced problem can be ∗ T

G∗(U ∗) = 0,

with

-

Let Sh(ξ) denote the operator which, with a distribution F in H −1(Ω)d, associates the part uh in Xh of the solution of the discrete Stokes problem :

.

(u · ∇)u − f & G∗(U ) = % . (u · ∇)T − g, T0

'



νh(ξ) ∇ uh ∀qh ∈ Mh,

*Sh(ξ)F*H 1(Ω)d ≤ c *F*H −1(Ω)d ,

and, for any G in H −1(Ω),

*Lh(G, 0)*H 1(Ω) ≤ c *G*H −1(Ω).

&

% & * S(ξ) − Sh(ξ) F*H 1(Ω)d ≤ c hs *S(ξ)F*H s+1(Ω)d + *F*H s−1(Ω)d , 5 and, with d−1 2 < σ ≤ 2, ( ) *(L −L h)(G, R0)*H 1(Ω) ≤ c hs *LG*H s+1(Ω) + hσ *R0* σ+ 1 . H 2 (∂Ω)



'



'



(div vh)(x)ph(x) dx = (F, vh),

(div uh)(x)qh(x) dx = 0.

1

a continuous function R0 in H 2 (∂Ω), associates the solution Th in Yh, equal to the Lagrange interpolate of R0 on ∂Ω, of the discrete Laplace equation ∀Sh ∈ Y0 h,

When setting Uh =

-

α

'



grad Th · grad Sh dx = (G, Sh).

.

uh , the reduced problem can be written Th

Fh(Uh) = Uh +

-

Sh(kh) 0 0 Lh

.

G(Uh) = 0.

Y = H01(Ω)d × H01(Ω).

Assumption 1. We consider a solution (u∗, p∗, k∗) of the reduced problem • which belongs to H s+1(Ω)d × H s(Ω) × H s+1(Ω) for a real number s, d 2 − 1 < s ≤ 2, • such that DF ∗(U ∗) is an isomorphism of Y.

hmin = min hK .

• Convergence properties For a real number s, 0 ≤ s ≤ 2, with appropriate regularity assumptions, %

∇ vh dx −

Let Lh denote the operator which, with a distribution G in H −1(Ω) and

About the Stokes and Laplace operators : • Stability properties For any F in H −1(Ω)d,

:

K∈Th

Assumption 2. The following property holds lim λh hs = 0,

h→0

with λh equal to | log hmin| for d = 2, to h−1 min for d = 3.

By combining all this, we prove the next lemmas. Lemma. There exists an h0 > 0 such that, for all h ≤ h0, DFh(U ∗) is an isomorphism of Y. Moreover, the norm of its inverse is bounded independently of h. Lemma. The mapping Fh satisfies the Lipschitz property, for all V1 and V2 in a bounded subset of Y, *DFh(V1) − DFh(V2)*E(Y) ≤ c λh *V1 − V2*Y .

F. Brezzi, J. Rappaz, P.-A. Raviart Theorem. If Assumptions 1 and 2 hold, there exist positive real numbers κ et h0 such that, for all h ≤ h0, • the discrete problem has a unique solution (uh, ph, kh) such that (uh, Th) belongs to the ball of Y with centre (u∗, T ∗) and radius κ λ−1 h , • this solution satisfies

Lemma. The next estimate is satisfied for a constant c(u∗, p∗, T ∗) only depending on the solution (u∗, p∗, T ∗) *Fh(U ∗)*Y ≤ c(u∗, p∗, T ∗) hs.

*u∗ − uh*H 1(Ω)d + *p∗ − ph*L2(Ω) + *T ∗ − Th*H 1(Ω)∩L∞(Ω) ≤ c(u∗, p∗, T ∗) hs,

where the constant c(u∗, p∗, T ∗) only depends on the solution (u∗, p∗, T ∗).

Conclusions • The previous estimate is fully optimal and yields a convergence of order h2 for a smooth solution. • The regularity assumption on the solution is likely in dimension d = 2 but not in dimension d = 3.

A posteriori analysis

However the convergence of the discretization can be proved when only the following assumptions hold lim λh h2 = 0.

h→0

This condition is much more likely. • Similar results hold for any triple (Xh, Mh, Yh) (of order at least 2 in dimension 3) such that an optimal inf-sup condition holds between Xh and Mh.

Some notation As usual for multi-step discretizations, the idea is to evaluate separately the two parts of the error, namely

(f )

Th

(s)

• the error due to the simplification of the model, • the error due to the discretization. The a posteriori analysis of this problem makes use of the non discrete implicit function theorem due to J. Pousin, J. Rappaz for bounding both errors from above.

Th

: set of elements of Th which are contained in Ωf ,

: set of elements of Th which are contained in Ωs.

For each K in Th, • EK : set of of edges (d = 2) or faces (d = 3) of K which are not contained in ∂Ω, • ωK : union of elements of Th that share at least an edge (d = 2) otr a face (d = 3) with K.

fh and gh : piecewise constant approximations of f and g, respectively.

Evaluating the error due to the simplification of the model Remark. In the implementation of the discrete problem, ν(·) is usually replaced by its Lagrange interpolate. For any continuous function ξ on Ω, we denote by νh(ξ) the function such that • its restriction to any K in Th belongs to P1(K), • which is equal to ν(ξ) at all vertices of K.

∀ξ ∈ IR,

 ν (ξ) h νh∗(x, ξ) = ν 0

for x in Ωf , for x in Ωs,

Let S0 denote the Stokes operator which, with a distribution F in H −1(Ω)d, associates the part u of the solution of the Stokes problem    −ν0 ∆u + grad p =

in Ω, in Ω, on ∂Ω.

F

div u = 0 u=0

 

The residual equation reads U − U∗ +

-

. (

S0 0 0 L

where the residual R is given by R=

-

%

)

G(U ) − G(U ∗) +

-

.

S0 0 0 L

div (ν(T ∗) − ν ∗(·, T ∗)) ∇u∗ (0, 0)

&.

.

R = 0,

• Indicator for the error due to the reduction of the model

J. Pousin, J. Rappaz 1,ρ

If Assumption 1(w) holds, there exists a neighbourhood of U in W0 (Ω)d× W 1,ρ(Ω) and a constant c only depending on U such that the following estimate holds for any solution (u∗, p∗, T ∗) of the reduced problem such that U ∗ = (u∗, T ∗) belongs to this neighbourhood &

%

*u − u∗*H 1(Ω)d + *T − T ∗*H 1(Ω) ≤ c * ν ∗(·, T ∗) − ν(T ∗) ∇u∗*L2(Ω)d×d .

≤c

(

/

(s) K∈Th

%

(s)

(ηK )2 + ε2 K

&) 1 2

%

&

ηK = * ν0 − νh(Th) ∇uh*L2(K)d×d .

A first technical difficulty. We have to evaluate the error due to the interpolation of the function ν, namely the quantity %

&

εK = * ν(Th) − νh(Th) ∇uh*L2(K)d×d .

Lemma. The following estimate holds εK ≤ c hK *grad Th*L∞(K)d *∇uh*L2(K)d×d .

Theorem. If Assumption 1(w) is satisfied and the solution U ∗ = (u∗, T ∗) 1,ρ exhibited above belongs to W0 (Ω)d × W 1,ρ(Ω), ρ >d , the following a posteriori estimate holds *u − u∗*H 1(Ω)d + *p − p∗*L2(Ω) + *T − T ∗*H 1(Ω)

%

(s)

Assumption 1(w). We consider a solution (u, p, T ) of the full problem 1,ρ • such that U = (u, T ) belongs to W0 (Ω)d × W 1,ρ(Ω), ρ >d , • DF (U ) is an isomorphism of H01(Ω)d × H 1(Ω).

&

+ c+ *u∗ − uh*H 1(Ωs)d + *T ∗ − Th*H 1(Ωs) ,

Proposition. Let (u∗, p∗, T ∗) be a solution of the simplified problem such 1,ρ that U ∗ = (u∗, T ∗) belongs to W0 (Ω)d × W 1,ρ(Ω), ρ > d. Thus, the follo(s)

wing estimate holds for the correponding error indicators ηK (

/

(s) K∈Th

(s)

(ηK )2

)1 2

(

≤ c *u − u∗*H 1(Ω)d + *p − p∗*L2(Ω) + *T − T ∗*H 1(Ω) + *u∗ − uh*H 1(Ωs)d + *T ∗ − Th*H 1(Ωs) +

where the constants c and c+ only depend on the norms of u∗ and T ∗. Proof. The estimate above, plus triangle inequalities, plus an inf-sup condition for the pressure term.

(d)1

(d)2

+ ηK

,

with & % (d)1 ηK = hK *fh + div νh∗(·, Th) ∇uh − (uh · ∇)uh − grad ph*L2(K)d 1 / + he2 *νh∗(·, Th) [∂nuh]e*L2(e)d + *div uh*L2(K), e∈EK

and / (d)2 ηK = hK *gh + α ∆Th − (uh · ∇)Th*L2(K) + e∈EK

ε2 K

(s)

K∈Th

&1 ) 2

.

J. Pousin, J. Rappaz

• Indicators linked to the finite element discretization (d)

/

Proof. We use the same triangle inequalities& as in the previous proof. % Next, bounding the quantity * ν ∗(·, T ∗) − ν(T ∗) ∇u∗*L2(Ω)d×d follows from the variational form of the residual equation.

Evaluating the error due to the discretization

ηK = ηK

%

1 h2

e *α [∂nTh ]e*L2 (e) .

Theorem. If (u∗, p∗, T ∗) is a solution of the reduced problem such that 1,ρ U ∗ = (u∗, T ∗) belongs to W0 (Ω)d × W 1,ρ(Ω), ρ > d, and DF ∗(U ∗) is an isomorphism of H01(Ω)d × H01(Ω), we have the following estimate *u∗ − uh*H 1(Ω)d + *p∗ − ph*L2(Ω) + *T ∗ − Th*H 1(Ω) ≤c + c++

( /

K∈Th

%

( /

K∈Th

(d)

(ηK )2

)1 2

+ c+

(

/

ε2 K

(f ) K∈Th

2 2 h2 K *f − fh *L2 (K)d + *g − gh *L2 (K)

&) 12

)1

2

+ c+++ *T0 − T0h*

1

H 2 (∂Ω)

,

for any solution (uh, ph, Th) of the reduced problem in a fixed neighbourhood of (u∗, p∗, T ∗).

Fully standard arguments (relying on appropriate choices of the function V in the residual equation) yield the following results. Proposition. Let (u∗, p∗, T ∗) be a solution of problem (3.5) − (3.6) such 1,ρ that U ∗ = (u∗, T ∗) belongs to W0 (Ω)d × W 1,ρ(Ω), ρ > d. Thus, the follo(d)2 (d)1 wing estimate holds for each error indicator ηK and ηK , K ∈ Th, ( (d)1 ηK ≤ c *u∗ − uh*H 1(ω )d + *p∗ − ph*L2(ω ) + *T ∗ − Th*H 1(ω ) K K K &) / % + hκ *f − fh*L2(κ)d + εκ . κ⊂ωK (d)2

ηK

(

≤ c *u∗ − uh*H 1(ω )d + *T ∗ − Th*H 1(ω ) + K K

/

κ⊂ωK

)

hκ *g − gh*L2(κ) .

Conclusions • The previous estimates are fully optimal ! Up to the terms ( /

K∈Th

%

2 2 h2 K *f − fh *L2 (K)d + *g − gh *L2 (K)

&) 12

and

*T0 − T0h*

1

H 2 (∂Ω)

,

which only depend on the data, the full error

E = *u − u∗*H 1(Ω)d + *p − p∗*L2(Ω) + *T − T ∗*H 1(Ω)

+ *u∗ − uh*H 1(Ω)d + *p∗ − ph*L2(Ω) + *T ∗ − Th*H 1(Ω)

satisfies the following equivalence property c

(

/

(s) K∈Th

(s)

(ηK )2 +

&) 21 / % (d) ηK )2 − ε2 ≤E K

K∈Th

E ≤ c+

(

/

(s)

(s)

(ηK )2 +

K∈Th

&) 12 / % (d) ηK )2 + ε2 , K K∈Th

with equivalence constants independent of the subdomains Ωs and Ωf and h.

A strategy for the automatic coupling

(d)

• The upper bound for the ηK is local. So it can be thought that these indicators provide an efficient tool for mesh adaptivity. (s)

• Even if no local estimate is established for the ηK , it can be hoped that they give a good representation of the simplification error.

The key point is to determine a deomposition of the domain in two subdomains Ωs and Ωf . • On Ωs, the Navier-Stokes and heat equations are uncoupled. • On Ωf , we solve the Navier-Stokes coupled with the heat equation. The decomposition is optimal if the errors due to the simplification of the model and to the discretization are of the same order. Initialization step : We choose a triangulation Th0 of the domain Ω that leads to a right approximation of the data. W. D¨ orfler We take Ω0 f =∅

et

Ω0 s = Ω,

and solve the discrete Navier–Stokes equations on Ω, next the heat equation.

Adaptation step : Assume that we are given a solution which is comn n puted with the decomposition into Ωn s and Ωf and a triangulation Th . (s)

ηh

=

1 n(s)

,(Th

)

/

(s)

ηK ,

n(s)

K∈Th

(d)

ηh

=

/ 1 (d) η . ,(Thn) K∈T n K h

Three substeps are needed for the adaptation :

• First m-substep : All elements K of Thn such that (s)

(s)

(s) (s) (d) ηK ≥ min{η h , η h } n+1 0 are inserted in a new subdomain Ω . f

ηK ≥ η h

or

0n+1 is regularized in order to obtain a new • Second m-substep : This Ω t partition of Ω into Ωn+1 and Ωn+1 . s f

Plus two algorithms



(d)

d-substep : By comparing the ηK

(d)

to their mean value η h , a new

triangulation Thn+1 is built such that the following property holds : The diameter of a new triangle or tetrahedron containing K or contained (d) (d) in K is of the same order as hK times η h /ηK .

• Uncoupling the unknowns m Find (um h , ph ) in Xh × Mh such that

∀vh ∈ Xh,

'



νh∗(x, Thm−1) ∇ um h : ∇ vh dx +

'

− Proving the convergence of mesh adaptation relies on now standard arguments. W. D¨ orfler P. Binev, W. Dahmen, R. DeVore But it seems that the convergence of the decomposition could only be checked by numerical experiments !

∀qh ∈ Mh, Find khm in ∀Sh ∈ Yh,



Yh such that α

'



'



'Ω

m (um h · ∇)uh · vh (x) dx



(div vh)(x)pm h (x) dx = (f , vh ),

(div um h )(x)qh (x) dx = 0,

grad Thm · grad Sh dx +

'



m (um h · ∇)Th · Sh (x) dx = (g, Sh ).

The convergence of this algorithm towards a solution of the discrete problem is proved for small enough data. T. Chac´ on Rebollo, S. Del Pino, D. Yakoubi

• Handling the nonlinear term in Navier–Stokes equations by the characteristic method O. Pironneau

About a turbulence model Joint work with

Four iterative algorithms are involved, without taking into account the final conjugate gradient.

Tom´ as Chac´ on Rebollo, Fr´ ed´ eric Hecht and Roger Lewandowski Of course, all iterations can be mixed.

The full model

α is a positive constant. B. Mohammadi, O. Pironneau

Ω : connected bounded open set in IRd, d = 2 or 3, with a Lipschitz– continuous boundary.  % &   −div ν(k) ∇u + (u · ∇)u + grad p =           div u = 0      −α ∆k = ν(k) |∇u|2          u=0         k=0

The viscosity ν is usually given by ∀ξ ≥ 0,

f

in Ω, in Ω, in Ω, on ∂Ω, on ∂Ω.

Unknowns : velocity u, pressure p, turbulent kinetic energy k. This is only a first draft of the real model !

ν(ξ) = ν0 + ν1

1

ξ.

Here, we work without restriction with a viscosity which is truncated ∀ξ ∈ IR,

and regularized : ν ∈ W 2,∞(IR).

ν0 ≤ ν(ξ) ≤ ν2,

Variational formulation

1,r +

Find (u, p, k) in H01(Ω)d × L2 0 (Ω) × W0 ∀v ∈ H01(Ω)d,

'



1 1 + + = 1. r r

and

r>d

ν(k) ∇ u : ∇ v dx +

'

(Ω) such that



∀q ∈ L2 0 (Ω), 1,r

∀χ ∈ W0 (Ω),

α

'







T. Gallou¨ et, R. Herbin R. Lewandowski C.B., T. Chac´ on Rebollo, R.L., F. Murat J. Lederer, R.L.

(u · ∇)u · v(x) dx −

'

• Existence result :

'



(div v)(x)p(x) dx = (f , v),

(div u)(x)q(x) dx = 0,

grad k · grad χ dx =

'



ν(k) |∇u|2 χ(x) dx.

• Uniqueness of the solution in dimension d = 2 for more regular and sufficiently small data. • A weaker regularity property as above.

It is readily checked that this problem is equivalent to the previous system.

The reduced model We introduce an arbitrary decomposition of the domain without overlap Ω = Ωt ∪ Ω.,

Ωt ∪ Ω. = ∅.  % &   −div ν ∗(·, k∗) ∇u∗ + (u∗ · ∇)u∗ + grad p∗ =           div u∗ = 0      −α ∆k∗ = ν ∗(·, k∗) |∇u∗|2        ∗   u =0        ∗  k =0

The key idea is as follows : • On Ω., we only consider and solve Navier–Stokes equations. • On Ωt, we solve the full problem.  ν(ξ) 4 ν + ν √ξ for a.e. x in Ωt, 0 1 ∀ξ ∈ IR, ν ∗(x, ξ) = ν for a.e. x in Ω., 0 Note that the function : ξ 0→ ν ∗(·, ξ) satisfies exactly the same properties as the function ν(·). For instance, ν0 ≤ ν ∗(x, ξ) ≤ ν2

∀ξ ∈ IR,

f

in Ω, in Ω, in Ω, on ∂Ω, on ∂Ω,

for a.e. x in Ω.

Variational formulation

r>d

and

1 1 + + = 1. r r

1,r +

Find (u∗, p∗, k∗) in H01(Ω)d × L2 0 (Ω) × W0 ∀v ∈ H01(Ω)d,

'



ν ∗(x, k∗) ∇ u∗ : ∇ v dx +

(Ω) such that

'



(u∗ · ∇)u∗ · v(x) dx '

− ∀q ∈ L2 0 (Ω), 1,r

∀χ ∈ W0 (Ω),

α

'





'





(div v)(x)p∗(x) dx = (f , v),

(div u∗)(x)q(x) dx = 0,

grad k∗ · grad χ dx =

'



ν ∗(x, k∗) |∇u∗|2 χ(x) dx.

The existence result is derived thanks to exactly the same arguments as for the full problem.

The discrete problem (Th)h : regular family of triangulations of Ω such that each element of Th is contained either in Ωt or in Ω.. Taylor–Hood finite elements *

+

Xh = vh ∈ H01(Ω)d; ∀K ∈ Th, vh|K ∈ P2(K)d , *

+

Mh = qh ∈ H 1(Ω) ∩ L20(Ω); ∀K ∈ Th, qh|K ∈ P1(K) .

*

+

Yh = χh ∈ W01,r (Ω); ∀K ∈ Th, χh|K ∈ P2(K) .

Find (uh, ph, kh) in Xh × Mh × Yh such that '

∀vh ∈ Xh,



ν ∗(x, kh) ∇ uh : ∇ vh dx +

'



Assumption 3. The following stability property holds

(uh · ∇)uh · vh(x) dx '

− ∀qh ∈ Mh, ∀χh ∈ Yh,

α

'





'





(div uh)(x)qh(x) dx = 0,

grad kh · grad χh dx =

'



*LhG*L∞(Ω) ≤ c *LG*L∞(Ω).

(div vh)(x)ph(x) dx = (f , vh),

νh(kh) |∇uh|2 χh(x) dx.

The a priori analysis of this problem again makes use of the discrete implicit function theorem due to F. Brezzi, J. Rappaz, P.-A. Raviart. However, a further assumption is needed.

A.H. Schatz, L.B. Wahlbin Theorem. If Assumptions 1, 2 and 3 hold, there exist positive real numbers κ et h0 such that, for all h ≤ h0, • the discrete problem has a unique solution (uh, ph, kh) such that (uh, kh) belongs to the ball of Y with centre (u∗, k∗) and radius κ λ−2 h , • this solution satisfies *u∗ − uh*H 1(Ω)d + *p∗ − ph*L2(Ω) + *k∗ − kh*H 1(Ω)∩L∞(Ω) ≤ c(u∗, p∗, k∗) hs,

where the constant c(u∗, p∗, k∗) only depends on the solution (u∗, p∗, k∗).

A posteriori analysis Conclusions

Let ρ and ρ+ be such that

• The previous estimate is fully optimal and yields a convergence of order h2 for a smooth solution. • The regularity assumption on the solution is likely in dimension d = 2 but not in dimension d = 3. However the convergence of the discretization can be proved when only the following assumptions hold

2 < ρ 3 such that, for all q, 2 ≤ q ≤ q˜0, and for any data f in W −1,q (Ω)d, any solution (u, p, k) of the full or reduced problem belongs to W 1,q (Ω)d × Lq (Ω) × W 1,q (Ω).

*u − u∗*W 1,ρ(Ω)d + *p − p∗*Lρ(Ω) + *k − k∗*W 1,ρ(Ω) %

&

≤ c η m + c(f ) *u∗ − uh*W 1,ρ(Ω )d + *k∗ − kh*W 1,ρ(Ω ) , . .

for any solution (u∗, p∗, k∗) of the reduced problem in a fixed neighbourhood of (u, p, k). It is also easy to bound from above the indicator η m.

Evaluating the error due to the discretization The same technical difficulty, due to the replacement of ν ∗(·, kh) by νh(kh). R. Verf¨ urth

εS K = hK *div

%%

& & ν ∗(·, kh) − νh(kh) ∇uh *Lρ(K)d 1 % & / + heρ *[ ν ∗(·, kh) − νh(kh) ∂nuh]e*Lρ(e)d , e∈EK % & ∗ 2 εL K = hK * ν (·, kh ) − νh (kh ) |∇uh | *Lρ (K) .

• Indicators for the error due to the discretization For each K in Th,

S + ηL , η K = ηK K

with %

&

S = h *f + div ν (k ) ∇u ηK K h h h h − (uh · ∇)uh − grad ph *Lρ (K)d

+

/

e∈EK

1

heρ *[νh(kh) ∂nuh]e*Lρ(e)d + *div uh*Lρ(K),

L = h *ν (k ) |∇u |2 + α ∆k * ηK K h h h h Lρ (K) +

Lemma. The following estimates hold εS K ≤ c hK *grad kh *L∞(K)d *∇uh *Lρ (ωK )d×d , 2 2 εL K ≤ c hK *grad kh *L∞(K)d *∇uh *L2ρ (K)d×d .

1

/

e∈EK

heρ *[α ∂nkh]e*Lρ(e).

Conclusions : The previous estimates are fully optimal : J. Pousin, J. Rappaz Theorem. If (u∗, p∗, k∗) is a solution of the full problem such that DF∗(U ∗) is an isomorphism of X , we have the following estimate

However, in view of the automatical insertion of the model (i.e. of the choice of Ωt and Ω.), we prefer to replace ηm by local quantities.

*u∗ − uh*W 1,ρ(Ω)d + *p∗ − ph*Lρ(Ω) + *k∗ − kh*W 1,ρ(Ω)

%

for any solution (uh, ph, kh) of the reduced problem in a fixed neighbourhood of (u∗, p∗, k∗). S and η L , K ∈ T , It is also easy to bound from above the indicators ηK h K by the local error.

&

η m = *div µ∗(·, kh) ∇uh *W −1,ρ(Ω)d + *µ∗(·, kh) |∇uh|2*W −1,ρ(Ω).

( / % &) 1ρ ρ ρ ρ ρ L ρ ≤c ηK + (εS , K ) + (εK ) + hK *f − fh *Lρ (K)d K∈Th

Th. : subset of elements of Th which are contained in Ω.. For each K in Th., we set 1

m = *µ∗ (·, k ) ∇u * ∗ 2 ηK h h Lρ (K)d×d + *µ (·, kh ) 2 ∇uh *Lρ (K)d×d .

Omega-T IsoValue 0 1 2

A first numerical experiment on the code FreeFem++ by F. Hecht and O. Pironneau

Lemma. The following upper bound holds ηm ≤

% /

K∈Th.

m )ρ (ηK

&1 ρ

.

The strategy for the automatic insertion of the turbulence model is exactly the same as previously.

Conclusions The convergence of the algorithm highly depends on the values of ν0 and ν1. But it leads to good results in some realistic situations. In view of the first experiments, we hope that the method will be efficient in non academic three-dimensional cases and, among them, for unsteady flows.

References [1] C. Bernardi, T. Chac´ on Rebollo, F. Hecht, R. Lewandowski — Automatic insertion of a turbulence model in the finite element discretization of the Navier–Stokes equations, Math. Models and Methods in Applied Sciences 19 (2009), 1139–1183. [2] C. Bernardi, F. Coquel, P.-A. Raviart — Automatic coupling and finite element discretization of the Navier–Stokes and heat equations, in preparation. [3] M. Braack, A. Ern — A posteriori control of modeling errors and discretization errors, Multiscale Model. Simul. 1 (2003), 221–238.

Thank you for your attention.

[5] T. Chac´ on Rebollo, S. Del Pino, D. Yakoubi — An iterative procedure to solve a coupled two-fluids turbulence model, Math. Model. Numer. Anal. [6] W. D¨ orfler — A convergent adaptive algorithm for Poisson’s equation, SIAM J. Numer. Anal. 33 (1996), 1106–1124. [7] J. Pousin, J. Rappaz — Consistency, stability, a priori and a posteriori errors for Petrov–Galerkin methods applied to nonlinear problems, Numer. Math. 69 (1994), 213–231. [8] R. Verf¨ urth — A Review of A Posteriori Error Estimation and Adaptive Mesh-Refinement Techniques, Wiley & Teubner (1996).

[4] F. Brezzi, J. Rappaz, P.-A. Raviart — Finite dimensional approximation of nonlinear problems, Part I : Branches of nonsingular solutions, Numer. Math. 36 (1980), 1–25.