A glimpse of boundary conditions

independent function w⁽⁰⁾.

In summary, the following observations can be stated:

• The irregular expansion (3.44) contributes together with the regular expan-sion to a complete comprehenexpan-sion of the dynamics generated by the lattice-Boltzmann equation. Thanks to the linearity of the system it is easily possible to combine both expansions to understand the evolution associated to arbi-trary initial values.

• Expansion (3.44) gives evidence for the existence of initial layers. This is a particular class of solutions decaying towards zero over a time interval of magnitude O(ǫ²).

Missing higher order terms in the regular expansion of the initial condition can be compensated by initial layers up to a maximal order depending on the regularity of the initial data.

• Finally, the parallel to section 2.3 should be noticed.

With regard to the first point, it remains to check whether any ansatz using another scaled time t/ǫ^α, α6= 2 can only approximate constant or trivial solutions of the lattice-Boltzmann equation (non-existence of further time scales).

3.5 A glimpse of boundary conditions

The conversion of macroscopic boundary conditions¹⁷ into kinetic ones for the lattice-Boltzmann equation is a non-trivial task and still subject to current research [67, 6]. This difficulty is not of purely discrete nature faced in the context of lattice-Boltzmann schemes. In this section we want to substantiate that the problem is of intrinsic kinetic origin and may arise already at the level of kinetic differential equations. The reason for this lies in the fact that kinetic equations involve other quantities as primary variables than the macroscopic ones. This impedes a straight imposition of macroscopic boundary conditions since those quantities are not di-rectly accessible which shall be kept on a certain value. Moreover, kinetic boundary conditions must mediate between those for the target equation and the requirements of the kinetic equation, that might be quite different from a mathematical point of view.

The force-driven Poiseuille flow through a canal provides a simple but illustrative example. Along the walls no-slip conditions apply while inflow and outflow are modeled by periodic boundary conditions. In section 1.2 it was pointed out, that the flow profile of such a quasi 1D flow is a solution of the diffusion equation.

Concretely, the Poiseuille flow is assumed to occur iny-direction. A symmetric in-terval [−w, w],w >0, on thex-axis may correspond to the cross-section of the canal.

The stationary flow profile [−w, w] ∋ x 7→ v(x) is then solution of a 1D Poisson

17Referring to the target equation.

equation (stationary diffusion equation), where homogeneous Dirichlet boundary conditions reflect the vanishing flow velocity at the walls.

−νv^′′(x) =g v(−w) = 0 =v(w)

⇒ v(x) = g

2ν(w²−x²). (3.46) On the right hand side, a constant gappears representing the strength of the force field which is directed parallel to the y-axis and induces the Poiseuille flow. The solution of (3.46) represents a parabola opened from below with roots at −w and w. Its maximum is located in the center of the canal (x= 0) where the flow reaches its highest speed depending on the ratio of the force g and the viscosity ν.

Let us investigate how (3.46) is approximated in the framework of the kinetic ap-proach using the D1P3 lattice-Boltzmann equation. Thanks to the stationarity of the target problem (3.46) we consider the time-independent lattice-Boltzmann equation¹⁸

1

ǫS∂_xf= _ǫ¹2J₀f+gw. (3.47) In order to avoid a reactive term appearing in (3.2), the reactivitycis set to 0, which results into the collision operator J₀. Equation (3.47) must be complemented by boundary conditions at w and −w.

It is well known that the bounce-back rule realizes a kinetic approximation of the macroscopic no-slip condition for the Stokes and Navier-Stokes equation particularly with regard to the D2P9 model. In section 1.2 it has been shown (see page 49) that the reduction from the D2P9 to the D1P3 model yields a bounce-back-like condition with sign-flipping. So we arrive at the following differential algebraic boundary value problem, which is the componentwise version of (3.47) multiplied by τ ǫ².

τ ǫ∂_xf₁(x) = f₁(x)− _2θ¹u(x) − _2θ¹ǫ²τ g f₁(w) = −f₂(w) 0 = f₀(x)−^θ−_θ¹u(x)−^θ−_θ¹ǫ²τ g

−τ ǫ∂_xf₂(x) = f₂(x)− _2θ¹u(x) − _2θ¹ǫ²τ g f₂(−w) = −f₁(−w)





(3.48) Let us derive a solution of this boundary value problem: Motivated by the solution of the target problem, we assume that the populations and their mass momentuare described by quadratic polynomials. To reduce the number of unknown variables, we exploit the symmetry of the problem.

So we require that u and the rest population f₀ are invariant under reflection at the origin: u(x) = u(−x) and f₀(x) = f₀(−x). This means that the zeros of the associated quadratic polynomials are symmetrically situated with respect to 0. In particular linear terms are excluded. In opposition, left and right moving popu-lations should obey the relation f₁(−x) = f₂(x). All assumptions are taken into account by the following ansatz

u(x) = a(n−x)(n+x) = −ax²+an² f0(x) = b(p−x)(p+x) = −bx²+bp²

f₁(x) = d(q−x)(r+x) = −dx²+d(q−r)x+dqr f₂(x) = d(r−x)(q+x) = −dx²+d(r−q)x+dqr

18In this section we do not stress theǫ-dependence by index notation. So we writef,ˆf instead offǫ,ˆfǫ.

3.5. A glimpse of boundary conditions 151

where a, b, d ≥ 0 and n, p, q, r ≥ 0 are respectively the amplitudes and roots of the unknown polynomials. It turns out that these seven coefficients are uniquely determined by (3.48).

• First, the boundary conditions for f₁ and f₂ yield: rq=w².

• Performing a comparison of coefficients in front of the x²-term and the con-stant term occurring in the algebraic equation forf⁽⁰⁾ ends up with the rela-tions:

b= ^θ−_θ¹a and a(z²−p²) =ǫ²τ g. (3.49)

• Similarly, we extract from the equation forf₁

d= _2θa¹ r−q = 2τ ǫ a(p²−w²) = 2aτ²ǫ²−gτ²ǫ² (3.50)

• The relation u=f₁+f₀+f₂ imposes the equation

(θ−1)n²+w² =θp². (3.51) Adding the right equations of (3.49) and (3.50), and solving for p gives p =

√w²+ 2τ ǫ².

Combining this with (3.51) producesn=q

w²+ 2^θ−_θ¹τ²ǫ². Utilizing a= ¹₂^θ_τg leads to the final result foru:

u(x) = _2ν^g

w²−x²+ 2^θ−_θ¹ǫ²τ²

=v(x) + ^g_ν^θ−_θ¹τ²ǫ² (3.52) Comparing the mass moment u with the solution v of the target problem (3.46) confirms the statement of corollary 3.4 i.e.:

hf,1i =u−−→^ǫ↓0 v more precisely v=u+ O(ǫ²).

Formula (3.52) corroborates that homogeneous Dirichlet boundary conditions are not exactly satisfied. This is symptomatic for other macroscopic boundary condi-tions, whose translation into kinetic boundary conditions is possibly less evident.

Returning to (3.52), the error term vanishes remarkably if θ attains its limit value 1. Hence, for the example of the Poiseuille flow, the solution of the target equation is produced exactly if the D1P2 model is used as a special case of the D1P3 model.

In general, kinetic boundary conditions are obtained as conditions on the the mo-ments of the population functionf.

bounce-back with sign-flipping f₁ =−f₂ hs²,fi= 0 → Dirichlet b.c. (3.53) bounce-back without sign-flipping f₁ = f₂ hs,fi= 0 → Neumann b.c. (3.54) In order to see what kind of macroscopic boundary condition ensues from (3.53), we insert the comparison function. More precisely the boundary conditions for the

mass moments u⁽⁰⁾, u⁽²⁾, ...of the asymptotic order functions are chosen such that the comparison function ˆf^[2n+3] satisfies the boundary condition with a residual as small as possible. Explicitly let us consider ˆf^[3] which is constructed alone by the solutionv=u⁽⁰⁾of the target equation (3.32). Using the expressions forf⁽⁰⁾,f⁽¹⁾,f⁽²⁾ and f⁽³⁾ in (3.38) we obtain

hˆf^[3],s²i= ¹_θu⁽⁰⁾ + 0·ǫ + ^θ_θ⁻2¹τ²∂²_xu⁽⁰⁾ǫ² + 0·ǫ³

Obviously, the right hand side is minimized if we require u⁽⁰⁾ to vanish at the boundaries, which yields Dirichlet conditions foru⁽⁰⁾. Hence (3.53) may be used to approximate the target problem (3.32) endowed with homogeneous Dirichlet instead of periodic boundary conditions, i.e. v(t,0) = v(t, L) = 0. ˆf^[3] leaves a residual of magnitude O(ǫ²) since it is not possible to additionally require ∂_x²u⁽⁰⁾ = 0 which would result in an ill-posed problem foru⁽⁰⁾. To obtain a higher order of consistency, we need another quantity that can compensate the ∂_x²u⁽⁰⁾-term. This is achieved by taking ˆf^[5] where u⁽²⁾ occurs in the associated asymptotic order function f⁽²⁾:

f⁽²⁾ =u⁽²⁾w+τ²∂_x²u⁽²⁾(s²w−¹_θw)

⇒ hˆf^[5],1i= ¹_θu⁽⁰⁾ + 0·ǫ + ¹_θu⁽²⁾+^θ_θ⁻₂¹τ²∂_x²u⁽⁰⁾

ǫ² + 0·ǫ³ So we may demand inhomogeneous Dirichlet conditions for u⁽²⁾ by setting

u⁽²⁾ =−^θ−θ¹τ²∂²_xu⁽⁰⁾ at the boundary.

In this way, ˆf^[5] reaches a consistency order of at least 4, because hf⁽³⁾,s²i must vanish due to lemma 3.4 (f⁽³⁾has only terms containingswhich is orthogonal to s²).

Let us finally note that the mass moment (3.52) of the example is recovered as sum of u⁽⁰⁾ =v and u⁽²⁾ =−^θ−_θ¹τ²∂_x²v = ^θ−_θ¹τ²_ν^g. Further terms do not appear in the expansion ofubecause they involve derivatives that vanish for the time-independent quadratic polynomial v.

The second kinetic boundary condition (3.54) is analyzed analogously. Here we obtain

hˆf^[3],si= 0 − ^τ_θ∂xu⁽⁰⁾ǫ + 0·ǫ² − τ^3θ_θ⁻2¹∂_x³u⁽⁰⁾−τ^{2 1}_θ∂t∂xv ǫ³ .

The residual becomes at most of magnitude O(ǫ³) if we set∂_xu⁽⁰⁾ = 0 at the bound-ary. Hence we can employ (3.54) to approximate target problems with Neumann¹⁹ boundary conditions. By using ˆf^[5], a higher order of consistency can be attained which yields an inhomogeneous Neumann condition for u⁽²⁾.

Finally, an objection with regard to condition (3.53) may be thrown in: It is true that bounce-back type conditions are well adapted to the hyperbolic character of lattice-Boltzmann equations, since they couple populations referring to out-going and incoming characteristics at the boundary points. But from the perspective of

19Observe that−ν∂xvis interpreted as the diffusive flux which is inversely proportional to the gradient of the density to level out spatial density fluctuations.

3.5. A glimpse of boundary conditions 153

the target equation, (3.53) is not a natural approach to realize Dirichlet boundary conditions upon the mass momentu. Instead, it would appear much more intuitive to require directly

hf,1i= 0. (3.55)

Since 1 ∈ F is orthogonal to s²w− ¹_θw in opposition to s², we infer hˆf^[3],1i = 0 if u⁽⁰⁾ = 0. Similarly we get hˆf^[5],1i = 0 if u⁽⁰⁾ = 0 and u⁽²⁾ = 0 at the boundary.

This time, homogeneous Dirichlet conditions are also found for the second and even for higher order functions of the mass moment.

Against this backdrop, the exceptional role of the D1P2 (θ= 1) model with respect to the accuracy of (3.52) becomes immediately comprehensible, because condition (3.53) and (3.55) coincide in this case as s² =1.

Remark 3.5. Both condition (3.53) and (3.55) are suitable for the D1P3 lattice-Boltzmann algorithm because there is only one population at each boundary node which is not updated by the standard evolution step. This population can be determined either by the bounce-back prescription or by solving a linear equation such that the mass moment remains zero. In two dimensions there are usually more than one population (even along a rectilinear boundary) where boundary conditions must be applied. In this case condition (3.55) yields only a single equation which is not enough to determine several unknown populations. In this situation bounce-back offers a reasonable possibility to cope with the under-determined problem.

Boundary conditions and energy estimate

Let us finally check whether the energy estimate in section 3.1 also works if we assume bounce-back type boundary conditions like (3.53) or (3.54) instead of pe-riodic boundary conditions. The essential step in the proof of theorem 3.1, where the periodicity of f and ˆf has been exploited, is the elimination of the integral

I :=

Z L 0

D

B f(t, x)−ˆf(t, x)

, S∂x f(t, x)−ˆf(t, x)E dx

= Z _L

0

1 2

d dx

D

B f(t, x)−ˆf(t, x)

, S f(t, x)−ˆf(t, x)E dx.

in (3.11). We consider now the boundary conditions

f₂(t,0) =±f₁(t,0) and f₂(t, L) =±f₁(t, L),

which are supposed to be satisfied likewise by the comparison function ˆf up to a residual of magnitude O(ǫ^α). Then

ˆf₂(t,0) =±ˆf₁(t,0) + O(ǫ^α) and ˆf₂(t, L) =±ˆf₁(t, L) + O(ǫ^α) Evaluating the integral by the fundamental theorem of calculus we obtain

I = 1 2

DB f−ˆf

, S f−ˆfE^(t,L)

(t,0)

which becomes more explicitly

I = 1 2

*



2θ f₁−ˆf₁

θ

θ−1 f₀−ˆf₀ 2θ f2−ˆf2



,



−1· f1−ˆf1 0· f₀−ˆf₀ +1· f₂−ˆf₂



 +

(t,L)

(t,0)

=θh

(f₂−ˆf₂)²−(f₁−ˆf₁)²i(t,L) (t,0)

=θh

(±1)² f₁−ˆf₁+ O(ǫ²)2

− f₁−ˆf₁2i(t,L) (t,0)

= O(ǫ^α).

In contrast to the periodic case the integral usually does not become zero since it is not possible to construct comparison functions, which satisfy the non-periodic boundary conditions exactly. However, the integral vanishes in the limit of ǫ↓ 0.

Therefore it acts like a small additional source term, which shows that the an energy estimate similar to theorem 3.1 is basically possible.

Outlook – What else?

The following aspects of the presented D1P3 lattice-Boltzmann equation might be interesting for further research.

• Lemma 3.1 does not refer to an equilibrium operator which includes advection, e.g. E_ǫf = (1−ǫ²τ c)hf,1iw+ǫaθhf,1iws, because thenBE_ǫ is not symmetric.

Here another approach might be necessary to close this gap.

• The energy estimate should be performed carefully for non-periodic bound-ary conditions particularly with respect to the resulting order of asymptotic similarity. Furthermore norms stronger thanL² should be taken into account to get also asymptotic similarity for the derivatives of f and ˆf. What are the boundary conditions for∂_xf, ∂_x²f, . . .?

• Instead of (3.53) and (3.54) one might consider the boundary conditionhf,s²−

1

θi= 0 which does not appear suspicious as it is analogue to the other moment conditions. However, the asymptotic analysis reveals quickly that it does not result in a macroscopic boundary conditions for u⁽⁰⁾ or ∂_xu⁽⁰⁾ but for

∂_x²u⁽⁰⁾. This leads to an ill-posed problem for u⁽⁰⁾, whose second spatial derivative cannot be prescribed as being a solution of a second order parabolic equation. If the boundary conditionhf,s²−¹_θi= 0 is well-posed for the lattice-Boltzmann equation, how does the solution look like? Certainly, it cannot be approximated by means of a regular expansion alone (arising of boundary layers).

• Another question concerns the natural translation of Neumann boundary con-ditions. So∂_xv= 0 suggests 0 =∂_xu=∂_xhf,1i=h∂_xf,1i. Where is the snag?

• Classical lattice-Boltzmann algorithms work on uniform (mostly quadratic or cubic) grids. Local grid refinement requires the coupling of grids with different

3.5. A glimpse of boundary conditions 155

mesh spacings (coarse grid↔fine grid). The grid spacing of the discrete algo-rithm is related to the sacling parameterǫof the lattice-Boltzmann equation (see end of subsection 1.1.2, page 41). Therefore it might be interesting from the point of grid coupling to study the model problem for a locally constant but discontinuous parameter function which replacesǫ.

• In the course of the computations particularly in section 3.2 it became notice-able, that the full information needed to determine an order function univo-cally was not directly found in the equations pertaining to this order but two orders later. This obliges to continue the expansion, which makes stronger regularity assumptions necessary concerning initial data and sources, as al-ready observed at the end of section 2.3. What can be said, if these regularity assumptions are not available? We will encounter related problems in analyz-ing numerical schemes throughout the next chapters.

Chapter 4

Consistency of a D1P3

Im Dokument Analysis of Lattice-Boltzmann Methods : Asymptotic and Numeric Investigation of a Singularly Perturbed System (Seite 161-169)