SLIDE 1
xxxx
Regularized 13 Moment Equations for Maxwell and Hard Sphere molecules
Henning Struchtrup
Anirudh Rana Alireza Mohammadzadeh
University of Victoria, Canada
Manuel Torrilhon
RWTH Aachen, Germany
SLIDE 2 The case for moment equations
- Microscopic theory: Boltzmann Equation
— microscopic variable: distribution function f (xi, t, ci) (7 independent variables!) — Direct Numerical Solutions are accurate, but numerically expensive — Direct Simulation Monte Carlo is powerful (molecules, reactions), but expensive
- Macroscopic transport equations
— Approximation to Boltzmann ⇒ limited range of validity Knn¿ 1 — collective behavior described by finite number N of macroscopic variables, e.g. ρ (xi, t) — density, vi (xi, t) — velocity, T (xi, t) — temperature σij (xi, t) — stress, qi (xi, t) — heat flux , ... — fast deterministic solutions — explicit equations, analytic solutions give deeper insight into processes
SLIDE 3 Are there meaningful macroscopic continuum approximations?
Burnett/super-Burnett equations
(CE expansion) [Burnett, 1935]
- unstable [Bobylev 1982], stable alternatives available [Bobylev, Söderholm, Zhong]
- boundary conditions not clear
- incomplete Knudsen layers, spurious oscillations
Grad moment equations [Grad 1949]
- unclear how many and which moments must be used
- boundary conditions not clear (for non-linear eqs)
- no Knudsen layers with 13 moments
Regularized 13 moment equations [Struchtrup & Torrilhon, since 2003]
- derived from Boltzmann equation (order of magnitude method)
- third order in Knudsen number (CE expanison =
⇒ super-Burnett)
- complete theory of boundary conditions
- linear equations: stable, H-theorem (incl. boundary conditions)
- phase speeds and damping of ultrasound waves agree to experiments
- smooth shock structures for all Ma, agree to DSMC for Ma<3
- reproduce all rarefaction effects for Kn . 0.5
- accessible for fast analytical and numerical solutions
- same methods =
⇒ R26 eqs [Gu&Emerson 2007]
SLIDE 4 Order of magnitude method [HS 2004] Kn =
mean free path macroscopic lengthscale of process
Step by step derivation of equations from Boltzmann eq.
¡ Kn0¢ : Euler
¡ Kn1¢ : Navier-Stokes-Fourier
¡ Kn2¢ : Grad 13
¡ Kn3¢ : regularized 13 moment equations (R13)
- stable equations at all orders [HS & MT 2003]
- accessible for arbitrary interaction potentials [HS 2005, HS&MT 2013 ]
= ⇒ today
SLIDE 5
Motivation: Dependence of moments on molecular interaction model
Couette flow at Kn=0.25, DSMC with VHS, VSS, Maxwell Hydrodynamic quantities not much affected:
velocity temperature shear stress normal heat flux
SLIDE 6
Motivation: Dependence of moments on molecular interaction model
Couette flow at Kn=0.25, DSMC with VHS, VSS, Maxwell Hydrodynamic quantities not much affected:
velocity temperature shear stress normal heat flux
Rarefaction quantities (Burnett etc.): visible dependence on molecular model
parallel heat flux pressure normal stress
SLIDE 7
Moment method in kinetic theory
kinetic equation: Knudsen number ε ∂f ∂t + ck ∂f ∂xk = 1 εS (f) moment method: replace kinetic equation with equations for moments uA = Z ψA (ci) fdc A = 1, . . . , N moment equations: multiply kinetic equation with ψA (ci) and integrate ∂uA ∂t + ∂FAk ∂xk = 1 εPA
fluxes: FAk = R ψA (ci) ckfdc , productions: PA= R ψA (ci) S (f) dc
Question 1: which moments: ψA (ci) = ?? Question 2: how many moments: N = ?? Question 3: how to close the equations: FAk (uB) = ??, PA (uB) = ?? Q3: Grad closure; Q1, Q2, Q3: Order of Magnitude method
SLIDE 8 Order of magnitude method
[HS 2004]
Step 1:
- set up moment system for arbitrary number of moments N
- close with Grad method
Step 2:
- Chapman-Enskog expansion to find leading ε−order of moments
- linear combination of moments such that number of moments at given
ε−order is minimal
- repeat for next order of magnitude
Step 3:
- use ε−orders to rescale equations for new moments
- use scaling for model reduction to a given order of accuracy
SLIDE 9 Step 1: Grad closure for linear moment equations
hydrodynamic moments ρ = Z fdc , ρvi = Z cifdc , ρu = 3 2ρθ = 1 2 Z C2fdc pij = pδij + σij = Z CiCjfdc , qi = 1 2 Z C2Cifdc general moments: a = 0, 1, . . . , Nn , n = 0, 1, 2, . . .
h· · · i indicates tracefree tensors
ua
i1···in =
Z C2aChi1Ci2 · · · Cinifdc with u0 = ρ , u0
i = 0 , u1 = 3ρθ = 3p , u0 ij = σij , u1 i = 2qi
Grad-type distribution: f|G = ρ √ 2πθ
3 exp
∙ −C2 2θ ¸ Ã 1 + X
n=0 Nn
X
a=0
λa
j1···jnC2aChj1Ci2 · · · Cjni
! coefficients λb
i1···in from inversion of
ua
i1···in =
Z C2aChi1Ci2 · · · Cinif|Gdc = ρn!
Nm
X
b=0
λb
i1···inθa+b+n(2 (a + b + n) + 1)!!
(2n + 1)!!
SLIDE 10 Step 1: Grad closure for linear moment equations (dimless)
Conservation laws ∂ρ ∂t + ∂vk ∂xk = 0 , ∂vi ∂t + ∂ρ ∂xi + ∂θ ∂xi + ∂σik ∂xk = Gi , 3 2 ∂θ ∂t + ∂vk ∂xk + ∂qk ∂xk = 0 equations for higher moments (renumbered), a = 1, . . . , N ∂ ˜ wa ∂t +
N
X
b=1
˜ R(1)
ab
∂˜ ub
k
∂xk − (2a + 3)!!2 (a + 1) 3 ∂qk ∂xk = −1 ε
N
X
b=1
˜ C(0)
ab ˜
wb ∂˜ ua
i
∂t +
N
X
b=1
˜ R(2)
ab
∂˜ ub
ik
∂xk + 1 3 ∂ ˜ wa ∂xi − (2a + 3)!! 3 ∂σik ∂xk + (2a + 3)!!a 3 ∂θ ∂xi = −1 ε
N
X
b=1
˜ C(1)
ab ˜
ub
i
∂˜ ua
ij
∂t + ∂˜ ua
ijk
∂xk + 2 5 ∂˜ ua
hi
∂xji + 2 (2a + 3)!! 15 ∂vhi ∂xji = −1 ε
N
X
b=1
˜ C(2)
ab ˜
ub
ij
∂˜ ua
ijk
∂t + ∂˜ ua
ijkl
∂xl + 3 7
N
X
b=1
˜ R(2)
ab
∂˜ ub
hij
∂xki = −1 ε
N
X
b=1
˜ C(3)
ab ˜
ub
ijk
etc. ˜ wa = ˜ ua − ˜ ua
E: non-equilibrium part of scalar moments
˜ R(n)
ab : coefficients from Grad closure (closure coefficients only for b = N)
˜ C(n)
ab : production matrices: collision term + Grad closure [Gupta & Torrilhon, RGD28]
ε : Knudsen number
SLIDE 11 Step 2: O of M method: CE expansion
Chapman-Enskog expansion of moments: for now, interested in zeroth and first order φ = φ0 + εφ1 + ε2φ2 + · · ·
- nly vector and 2-tensor moments have first order contributions:
˜ ub
i,1 = −κb
∂θ ∂xi , ˜ ub
ij,1 = −μb
∂vhi ∂xji
κa =
N
X
b=1
h ˜ C(1)i−1
ab (2b + 3)!!b
3 , μa =
N
X
a=1
h ˜ C(2)i−1
ab (2b + 3)!! 2
15
Navier-Stokes-Fourier qNSF
i
= −ˆ κ ∂θ ∂xi , σNSF
ij
= −2ˆ μ∂vhi ∂xji heat conductivity: ˆ κ = εκ1
2 , viscosity ˆ
μ = ε (dimless), Prandtl: Pr = 5
2 ˆ μ ˆ κ = 5 κ1
= ⇒ ρ, vi, θ are O ¡ ε0¢ = ⇒ ˜ ub
i, ˜
ub
ij are O
¡ ε1¢ = ⇒ all others are at least O ¡ ε2¢
SLIDE 12
Step 2: O of M method: Reconstructing variables
to first order: linear dependence (a = 2, . . . , N) ˜ ua
i,1 = −κa
∂θ ∂xi = 2κa κ1 qi,1 ˜ ua
ij,1 = −μa
∂vhi ∂xji = μa μ1 σij,1 define second order variables (a = 2, . . . , N) ˜ wa
i = ˜
ua
i − 2κa
κ1 qi ˜ wa
ij = ˜
ua
ij − μa
μ1 σij = ⇒ ρ, vi, θ are O ¡ ε0¢ = ⇒ qi, σij are O ¡ ε1¢ = ⇒ all others are at least O ¡ ε2¢ moments describe collective behavior = ⇒ new moments ˜ wa
i , ˜
wa
ij depend on interaction C(n) ab
SLIDE 13 Step 2: O of M method: Next CE expansion
higher moments to leading order O ¡ ε2¢ ˜ wb = −ζbε∂qk ∂xk ˜ wb
i = −ϑbε∂σik
∂xk − ηb ∙ qi + 5ε 2 Pr ∂θ ∂xi ¸ ˜ wb
ij = −ϕbε ∂qhi
∂xji − φb ∙ σij + 2ε∂vhi ∂xji ¸ ˜ ub
ijk = −ξbε∂σhij
∂xki
coefficients depend on interaction, closure
ζc =
N
X
a=1
h ˜ C(0)i−1
ca
" 2
N
X
b=1
˜ R(1)
ab
κb κ1 − (2a + 3)!!2 (a + 1) 3 # ϑb =
N
X
a=2
h ˜ D(1)i−1
ba
" N X
c=1
˜ R(2)
ac
μc μ1 − (2a + 3)!! 3 − κa κ1 ∙μ2 μ1 − 5 ¸# , ηb = 2 κ1
N
X
a=2
h ˜ D(1)i−1
ba
∙a (2a + 3)!! 3 − 5κa κ1 ¸ , ˜ D(1)
ab =
∙ ˜ C(1)
ab − κa
κ1 ˜ C(1)
1b
¸ ϕb = 4 5
N
X
a=2
h ˜ D(2)i−1
ba
∙κa κ1 − μa μ1 ¸ , φb = 2 μ1
N
X
a=2
h ˜ D(2)i−1
ba
∙(2a + 3)!! 15 − μa μ1 ¸ , ˜ D(2)
ab =
∙ ˜ C(2)
ab − μa
μ1 ˜ C(2)
1b
¸ ξb = 3 7
N
X
a=1
h ˜ C(3)i−1
ba N
X
c=1
˜ R(2)
ac
μc μ1
SLIDE 14 Step 3: O of M method: Model reduction by order
new moment equations for qi, σij, ˜ wa, ˜ wa
i , ˜
wa
ij, ˜
ua
ijk (after eliminating time derivatives of qi and σij)
variables scaled by Knudsen order σij = ε¯ σij , qi = ε¯ qi , ˜ wa = ε2 ¯ wa , ˜ wa
i = ε2 ¯
wa
i ,
˜ wa
ij = ε2 ¯
wa
ij , ˜
ua
ijk = ε2¯
ua
ijk ,
˜ ua
ijkl = ε3¯
ua
ijkl .
Knudsen scaled equations (from now on: N = 3)
ε ∙∂¯ qi ∂t + a(1,1)∂¯ σik ∂xk ¸ +ε2 ∙1 2 ∂ ¯ w2
ik
∂xk + 1 6 ∂ ¯ w1 ∂xi ¸ = −a(1,0) ∙ ¯ qi + 5 2 Pr ∂θ ∂xi ¸ −ε1 2
3
X
b=2
˜ C(1)
1b
∙ ¯ wb
i + ϑb
∂¯ σik ∂xk + ηb ∙ ¯ qi + 5 2 Pr ∂θ ∂xi ¸¸ ε ∙∂¯ σij ∂t + a(2,1) ∂¯ qhi ∂xji ¸ + ε2∂¯ u1
ijk
∂xk = −a(2,0) ∙ ¯ σij + 2 ∂vhi ∂xji ¸ − ε
3
X
b=2
˜ C(2)
1b
∙ ¯ wb
ij + ϕb
∂¯ qhi ∂xji + φb ∙ ¯ σij + 2 ∂vhi ∂xji ¸¸ ε2 " ∂ ¯ wa ∂t +
3
X
b=2
˜ R(1)
ab
∂ ¯ wb
k
∂xk # = −ε
3
X
b=1
˜ C(0)
ab
∙ ¯ wb + ζb ∂¯ qk ∂xk ¸ (a = 1, 2, 3) ε2 " ∂ ¯ wa
i
∂t + 1 3 ∂ ¯ wa ∂xi − κa κ1 1 3 ∂ ¯ w1 ∂xi +
3
X
b=2
˜ R(2)
ab
∂ ¯ wb
ik
∂xk − κa κ1 ∂ ¯ w2
ik
∂xk # = −ε
3
X
b=2
˜ D(1)
ab
∙ ¯ wb
i + ϑ(1) b
∂¯ σik ∂xk + ηb ∙ ¯ qi + 5 2 Pr ∂θ ∂xi ¸¸ (a = 2, 3) ε2 " ∂ ¯ wa
ij
∂t + ∂¯ ua
ijk
∂xk − μa μ1 ∂¯ u1
ijk
∂xk + 2 5 ∂ ¯ wa
hi
∂xji # = −ε
3
X
b=2
˜ D(2)
ab
∙ ¯ wb
ij + ϕb
∂¯ qhi ∂xji + φb ∙ ¯ σij + 2 ∂vhi ∂xji ¸¸ (a = 2, 3) ε2 " ∂¯ ua
ijk
∂t + 3 7
3
X
b=2
˜ R(2)
ab
∂ ¯ wb
hij
∂xki # + ε3∂¯ ua
ijkl
∂xl = −ε
N
X
b=1
˜ C(3)
ab
∙ ¯ ub
ijk + ξb
∂¯ σhij ∂xki ¸ (a = 1, 2, 3)
SLIDE 15 Step 3: O of M method: Model reduction by order
consider contribution to conservation laws of order O ¡ ελ¢ first order (λ = 1): leading order terms for qi, σij (underlined)
ε ∙∂¯ qi ∂t + a(1,1)∂¯ σik ∂xk ¸ + ε2 ∙1 2 ∂ ¯ w2
ik
∂xk + 1 6 ∂ ¯ w1 ∂xi ¸ = −a(1,0) ∙ ¯ qi + 5 2 Pr ∂θ ∂xi ¸ − ε1 2
3
X
b=2
˜ C(1)
1b
∙ ¯ wb
i + ϑb
∂¯ σik ∂xk + ηb ∙ ¯ qi + 5 2 Pr ∂θ ∂xi ¸¸ ε ∙∂¯ σij ∂t + a(2,1) ∂¯ qhi ∂xji ¸ + ε2∂¯ u1
ijk
∂xk = −a(2,0) ∙ ¯ σij + 2 ∂vhi ∂xji ¸ − ε
3
X
b=2
˜ C(2)
1b
∙ ¯ wb
ij + ϕb
∂¯ qhi ∂xji + φb ∙ ¯ σij + 2 ∂vhi ∂xji ¸¸
= ⇒ Navier-Stokes-Fourier (again): qi = ε¯ qi = − 5 2 Prε ∂θ ∂xi σij = ε¯ σij = −2ε∂vhi ∂xji
SLIDE 16 Step 3: O of M method: Model reduction by order
consider contribution to conservation laws of order ελ second order (λ = 2): also next order terms for qi, σij (double underlined)
ε ∙∂¯ qi ∂t + a(1,1)∂¯ σik ∂xk ¸ + ε2 ∙1 2 ∂ ¯ w2
ik
∂xk + 1 6 ∂ ¯ w1 ∂xi ¸ = −a(1,0) ∙ ¯ qi + 5 2 Pr ∂θ ∂xi ¸ − ε1 2
3
X
b=2
˜ C(1)
1b
∙ ¯ wb
i + ϑb
∂¯ σik ∂xk + ηb ∙ ¯ qi + 5 2 Pr ∂θ ∂xi ¸¸ ε ∙∂¯ σij ∂t + a(2,1) ∂¯ qhi ∂xji ¸ + ε2∂¯ u1
ijk
∂xk = −a(2,0) ∙ ¯ σij + 2 ∂vhi ∂xji ¸ − ε
3
X
b=2
˜ C(2)
1b
∙ ¯ wb
ij + ϕb
∂¯ qhi ∂xji + φb ∙ ¯ σij + 2 ∂vhi ∂xji ¸¸
qi, σij balances not enough, needs also leading contributions for ¯ wb
i, ¯
wb
ij
¯ wb
i = −ϑb
∂¯ σik ∂xk − ηb ∙ ¯ qi + 5 2 Pr ∂θ ∂xi ¸ ¯ wb
ij = −ϕb
∂¯ qhi ∂xji − φb ∙ ¯ σij + 2 ∂vhi ∂xji ¸ note: equations and coefficients a(α,β) constructed such that ¯ wb
i, ¯
wb
ij cancel to leading order
SLIDE 17 Step 3: O of M method: Model reduction by order
consider contribution to conservation laws of order ελ second order (λ = 2) generalized Grad 13 moment equations (linear): ∂qi ∂t + a(1,1)∂σik ∂xk = −a(1,0) ε ∙ qi + 5 2 Prε ∂θ ∂xi ¸ ∂σij ∂t + a(2,1) ∂qhi ∂xji = −a(2,0) ε ∙ σij + 2ε∂vhi ∂xji ¸
coefficients:
a(1,0) = 5 κ1 − 1 2
3
X
b=2
˜ C(1)
1b ηb
a(1,1) = μ2 2μ1 − 5 2 − 1 2
3
X
b=2
˜ C(1)
1b ϑb ,
a(1,2) = φ2 2 − 1 2
N
X
a,b=2
˜ C(1)
1b
h ˜ D(1)
ba
i−1 " 3 X
c=2
˜ R(2)
ac φc − κa
κ1 φ2 − a(2,0)ϑa − a(1,1)ηa # a(1,3) = ϕ2 2 − 1 2
N
X
a,b=2
˜ C(1)
1b
h ˜ D(1)
ba
i−1 " 3 X
c=2
˜ R(2)
ac ϕc − κa
κ1 ϕ2 − ϑaa(2,1) # a(1,4) = ζ1 6 − 1 2
N
X
a,b=2
˜ C(1)
1b
h ˜ D(1)
ba
i−1 ∙ζa 3 − κa κ1 ζ1 3 − 5 3 Prηa ¸
SLIDE 18 Step 3: O of M method: Model reduction by order
consider contribution to conservation laws of order ελ third order (λ = 3): full eqs for qi, σij, next order terms for ˜ w2
i , ˜
w3
i , ˜
w2
ij, ˜
w3
ij ∂qi ∂t +a(1,1)∂σik ∂xk +1 2 ∂ ˜ w2
ik
∂xk −ζ1 6 ε ∂ ∂xi ∂qk ∂xk = −1 εa(1,0) ∙ qi + 5 2 Prε ∂θ ∂xi ¸ −1 ε 1 2
N
X
b=2
˜ C(1)
1b
∙ ˜ wb
i + ϑbε∂σik
∂xk + ηb ∙ qi + 5 2 Prε ∂θ ∂xi ¸¸ ∂σij ∂t + a(2,1) ∂qhi ∂xji − ξ1ε ∂ ∂xk ∂σhij ∂xki = −1 εa(2,0) ∙ σij + 2ε ∂vhi ∂xji ¸ − 1 ε
3
X
b=2
˜ C(2)
1b
∙ ˜ wb
ij + ϕbε ∂qhi
∂xji + φb 1 ε ∙ σij + 2ε ∂vhi ∂xji ¸¸ ∂ ˜ wa
i
∂t − ∙ζa 3 − ζ1 3 κa κ1 ¸ ε ∂ ∂xi ∂qk ∂xk +
3
X
b=2
˜ R(2)
ab
∂ ˜ wb
ik
∂xk −κa κ1 ∂ ˜ w2
ik
∂xk = −1 ε
3
X
b=2
˜ D(1)
ab
∙ ˜ wb
i + ϑbε∂σik
∂xk + ηb ∙ qi + 5 2 Prε ∂θ ∂xi ¸¸ (a = 2, 3) ∂ ˜ wa
ij
∂t − ∙ ξa − ξ1 μa μ1 ¸ ε ∂ ∂xk ∂σhij ∂xki + 2 5 ∂ ˜ wa
hi
∂xji = −1 ε
3
X
b=2
˜ D(2)
ab
∙ ˜ wb
ij + ϕbε ∂qhi
∂xji + φb ∙ σij + 2ε ∂vhi ∂xji ¸¸ (a = 2, 3)
+ conservation laws closed set for 29 variables: ρ, vi, θ, σij, qi, ˜ w2
i , ˜
w3
i , ˜
w2
ij, ˜
w3
ij
for Maxwell molecules: ˜ C(n)
1b are lower triangular, no coupling to ˜
wa
i , ˜
wa
ij on the right
= ⇒ 13 variables: ρ, vi, θ, σij, qi for non-Maxwell molecules: ˜ C(n)
1b are fully occupied
= ⇒ use Knudsen order arguments for further redcution to 13 variables: ρ, vi, θ, σij, qi
SLIDE 19 Step 3: O of M method: Model reduction by order
further reduction at third order: leading order CE expansion of ˜ w2
i , ˜
w3
i , ˜
w2
ij, ˜
w3
ij
some detailed steps . . . Example: ∂ ∂t ∙ qi + 5 2 Prε ∂θ ∂xi ¸ required to second order O ¡ ε2¢ time derivative of Burnett heat flux (Burnett extracted from Moment system = ⇒ later) ∂ ∂t ∙ qi + 5 2 Prε ∂θ ∂xi ¸ = ε2 a(1,0) ∙ 2a(1,1) ∂2 ∂xk∂xhi ∂vki ∂t − 5 3 Pr ∂2 ∂xk∂xi ∂vk ∂t ¸ + O ¡ ε3¢ use conservation laws (Euler suffices) (with q(1)
k
= − 5
2 Pr ∂θ ∂xk)
∂ ∂t ∙ qi + 5 2 Prε ∂θ ∂xi ¸ = ε2 a(1,0) ∙ −2a(1,1) + 5 2 Pr ¸ ∂ ∂xk ⎛ ⎝ ∂2ρ ∂xhi∂xki − 2 Pr 5 ∂q(1)
hi
∂xki ⎞ ⎠ + O ¡ ε3¢ space derivative of Burnett stress (with q(1)
k
= − 5
2 Pr ∂θ ∂xk)
∂ ∂xk ∙ σik + 2ε∂vki ∂xhi ¸ = ε2 a(2,0) ⎡ ⎣−2 ∂ ∂xk ∂2ρ ∂xhi∂xki + µ4 Pr 5 − a(2,1) ¶ ∂ ∂xk ∂q(1)
hi
∂xki ⎤ ⎦ + O ¡ ε3¢ Combine: Burnett used “forwards, backwards” ∂ ∂t ∙ qi + 5 2 Prε ∂θ ∂xi ¸ = £ 2a(1,1) −
5 2 Pr
¤ 2a(1,0) ∙ a(2,1)ε ∂ ∂xk ∂qhi ∂xki + a(2,0) ∂ ∂xk ∙ σik + 2ε∂vki ∂xhi ¸¸ + O ¡ ε3¢
SLIDE 20
Step 3: O of M method: Model reduction by order
final equations at third order: Regularized 13 moment equations (linear)
∂ρ ∂t + ∂vk ∂xk = 0 ∂vi ∂t + ∂ρ ∂xi + ∂θ ∂xi + ∂σik ∂xk = Gi 3 2 ∂θ ∂t + ∂vk ∂xk + ∂qk ∂xk = 0 ∂qi ∂t +a(1,1)∂σik ∂xk −a(1,2) ∂ ∂xk ∙ σik + 2ε ∂vhi ∂xki ¸ −a(1,3)ε ∂ ∂xk ∂qhi ∂xki −a(1,4)ε ∂ ∂xi ∂qk ∂xk = −a(1,0) ε ∙ qi + 5 2 Prε ∂θ ∂xi ¸ ∂σij ∂t +a(2,1) ∂qhi ∂xji −a(2,2) ∂ ∂xhi ∙ qji + 5 2 Prε ∂θ ∂xji ¸ −a(2,3)ε ∂ ∂xk ∂σhij ∂xki −a(2,4)ε ∂2σij ∂xk∂xk = −a(2,0) ε ∙ σij + 2ε ∂vhi ∂xji ¸ coefficients a(α,β) determined through Grad closure, collision term hyperbolic + diffusion
SLIDE 21 Regularized 13 moment equations coefficients a(α,β) determined through Grad closure, collision term
κa =
N
X
b=1
h ˜ C(1)i−1
ab (2b + 3)!!b
3 , μa =
N
X
a=1
h ˜ C(2)i−1
ab
2 (2b + 3)!! 15 ζc =
N
X
a=1
h ˜ C(0)i−1
ca
" 2
N
X
b=1
˜ R(1)
ab
κb κ1 − (2a + 3)!!2 (a + 1) 3 # , ϑb =
N
X
a=2
h ˜ D(1)i−1
ba
" N X
c=1
˜ R(2)
ac
μc μ1 − (2a + 3)!! 3 − κa κ1 ∙μ2 μ1 − 5 ¸# ηb = 2 κ1
N
X
a=2
h ˜ D(1)i−1
ba
∙a (2a + 3)!! 3 − 5κa κ1 ¸ , ˜ D(1)
ab = ˜
C(1)
ab − κa
κ1 ˜ C(1)
1b
, ϕb = 4 5
N
X
a=2
h ˜ D(2)i−1
ba
∙κa κ1 − μa μ1 ¸ φb = 2 μ1
N
X
a=2
h ˜ D(2)i−1
ba
∙(2a + 3)!! 15 − μa μ1 ¸ , ˜ D(2)
ab = ˜
C(2)
ab − μa
μ1 ˜ C(2)
1b
, ξb = 3 7
N
X
a=1
h ˜ C(3)i−1
ba N
X
c=1
˜ R(2)
ac
μc μ1 a(1,0) = 5 κ1 − 1 2
3
X
b=2
˜ C(1)
1b ηb , a(1,1) = μ2
2μ1 − 5 2 − 1 2
3
X
b=2
˜ C(1)
1b ϑb , a(1,4) = ζ1
6 − 1 2
N
X
a,b=2
˜ C(1)
1b
h ˜ D(1)
ba
i−1 ∙ζa 3 − κa κ1 ζ1 3 ¸ a(1,2) = φ2 2 − 1 2
N
X
a,b=2
˜ C(1)
1b
h ˜ D(1)
ba
i−1 ( 3 X
c=2
˜ R(2)
ac φc − κa
κ1 φ2 − a(2,0) ∙ ϑa − ηa ∙a(1,1) a(1,0) − 5 4 Pr a(1,0) ¸¸) a(1,3) = ϕ2 2 − 1 2
N
X
a,b=2
˜ C(1)
1b
h ˜ D(1)
ba
i−1 ( 3 X
c=2
˜ R(2)
ac ϕc − κa
κ1 ϕ2 − a(2,1) ∙ ϑa − ηa ∙a(1,1) a(1,0) − 5 4 Pr a(1,0) ¸¸) a(2,0) = 2 μ1 −
3
X
b=2
˜ C(2)
1b φb , a(2,1) = 4
5 −
2
X
b=2
˜ C(2)
1b ϕb ,
a(2,2) =
3
X
a,b=2
˜ C(2)
1b
h ˜ D(2)i−1
ba
½ −2 5ηa + a(1,0) ∙ ϕa + φa ∙ 2 Pr a(2,0) − a(2,1) a(2,0) ¸¸¾ a(2,3) = ξ1 −
3
X
a,b=2
˜ C(2)
1b
h ˜ D(2)
ba
i−1 ½ ξa + ϑa − ξ1 μa μ1 − 5 2a(1,1) ∙ ϕa + φa ∙ 2 Pr a(2,0) − a(2,1) a(2,0) ¸¸¾ a(2,4) = 1 3
3
X
a,b=2
˜ C(2)
1b
h ˜ D(2)
ba
i−1 ½ ϑa − 5 2a(1,1) ∙ ϕa + φa ∙ 2 Pr a(2,0) − a(2,1) a(2,0) ¸¸¾
SLIDE 22 Regularized 13 moment equations
coefficients a(α,β) determined through Grad closure, collision term
∂qi ∂t + a(1,1)∂σik ∂xk − a(1,2) ∂ ∂xk ∙ σik + 2ε ∂vhi ∂xki ¸ − a(1,3)ε ∂ ∂xk ∂qhi ∂xki − a(1,4)ε ∂ ∂xi ∂qk ∂xk = −1 εa(1,0) ∙ qi + 5 2 Prε ∂θ ∂xi ¸ ∂σij ∂t + a(2,1) ∂qhi ∂xji − a(2,2) ∂ ∂xhi ∙ qji + 5 2 Prε ∂θ ∂xji ¸ − a(2,3)ε ∂ ∂xk ∂σhij ∂xki − a(2,4)ε ∂2σij ∂xk∂xk = −1 εa(2,0) ∙ σij + 2ε ∂vhi ∂xji ¸ Maxwell hard spheres BGK Pr
2 3 = 0.6667
0.660851 1 a(1,0)
2 3 = 0.6667
0.650061 1 a(1,1) 1 0.786941 1 a(1,2) 0.186614 a(1,3)
12 5 = 2.4
2.10417
14 5 = 2.8
a(1,4) 2 1.47742
4 3 = 1.3333
a(2,0) 1 0.98632 1 a(2,1)
4 5 = 0.8
0.631247
4 5 = 0.8
a(2,2) 0.0925568 a(2,3) 2 2.15033 3 a(2,4) −0.102588
SLIDE 23 R 13 distribution function
f = ρ m exp ³ −
ˆ C2 2
´ √ 2πθ
3
∙ 1 + Aμ ˆ Chi ˆ Cjiσij + Aκ ˆ Ciqi − Aζ µ ε∂qk ∂xk ¶ − Aϕ ˆ Chi ˆ Cjiε ∂qhi ∂xji − Aξ ˆ Chi ˆ Cj ˆ Cki µ ε∂σhij ∂xki ¶ −Aϑ ˆ Ci µ ε∂σik ∂xk ¶ − Aη ˆ Ci µ qi + εκ1 2 ∂θ ∂xi ¶ − Aφ ˆ Chi ˆ Cji µ σij + εμ1 ∂vhi ∂xji ¶¸
Maxwell molecules: Sonine polynomials Aμ = 1 2 , Aκ = −1 + ˆ C2 5 , Aζ = 3 2 − ˆ C2 + ˆ C4 10 Aϕ = 6 35 ˆ C2 − 6 5 , Aξ = 1 3 , Aϑ = Aη = Aφ = 0 hard spheres (3 × 3 matrices, round-off): Aμ = 0.644446 − 0.0259338C2 + 0.00058875C4 Aκ = −1.21624 + 0.296686C2 − 0.00909016C4 + 0.000161778C6 Aζ = 1.7532 − 1.35824C2 + 0.200694C4 − 0.00770945C6 + 0.000127611C8 Aϕ = −1.59174 + 0.292273C2 − 0.007209C4 Aξ = 0.594947 − 0.0412998C2 + 0.0009871C4 Aϑ = −0.0419332 + 0.0222017C2 − 0.00274908C4 + 0.0000861661C6 Aη = 0.154386 − 0.0711059C2 + 0.0070829C4 − 0.000148437C6 Aφ = −0.138697 + 0.0265345C2 − 0.000746734C4
SLIDE 24 Burnett and super-Burnett equations
CE expansion of generalized R13 (linear)
qi = − 5 2 Prε ∂θ ∂xi +ε2 ∙θB
4
2 ∂2vi ∂xk∂xk + 2 3 µθB
4
4 − θB
2
¶ ∂2vk ∂xk∂xi ¸ −ε3 ∙ θsB
1
∂3ρ ∂xi∂xk∂xk + θsB
2
∂3θ ∂xk∂xk∂xi ¸ σij = −2ε ∂vhi ∂xji −ε2 ∙ $B
2
∂2ρ ∂xhi∂xji + ¡ $B
2 − $B 3
¢ ∂2θ ∂xhi∂xji ¸ +ε3 ∙ $sB
1
∂2 ∂xhi∂xji ∂vk ∂xk − $sB
2
∂2 ∂xk∂xk ∂vhi ∂xji ¸ Maxwell hard spheres BGK Pr
2 3 = 0.6667
0.660851 1 θB
2 45 8 = 5.625
5.81945
5 2 = 2.5
θB
4
3 2.42113 2 $B
2
2 2.02774 2 $B
3
3 2.42113 2 θsB
1 5 8 = 0.625
2.23673 −1 θsB
2 157 16 = 9.813
5.81182
25 6 = 4.1667
$sB
1 5 3 = 1.6667
0.493778 −2
3 = −0.6667
$sB
2 4 3 = 1.3333
1.16695 2 Burnett coefficients agree with literature super-Burnett coefficients for non-Maxwellian molecules: first time
SLIDE 25 Example: Dispersion and Damping of Ultrasound
left: inverse phase velocity q
5 3RT/vph
dashed: Maxwell molecules right: reduced damping α/ω continuous: hard spheres dots: measurement by Meyer&Sessler Navier-Stokes Fourier
1 2 3 4 5 6 0.4 0.5 0.6 0.7 0.8 0.9 1.0 1w 5 3 vph
NSF
1 2 3 4 5 6 0.00 0.05 0.10 0.15 0.20 0.25 0.30 1w
NSF
viscosity: μMM = μHS Prandtl number: PrMM = 0.667, PrHS = 0.661 = ⇒ close agreement between MM, HS; insufficient for higher ω
SLIDE 26 Example: Dispersion and Damping of Ultrasound
left: inverse phase velocity q
5 3RT/vph
dashed: Maxwell molecules right: reduced damping α/ω continuous: hard spheres Burnett / super-Burnett: super-Burnett should be better than Burnett!!
1 2 3 4 5 6 0.4 0.5 0.6 0.7 0.8 0.9 1.0 1w 5 3 vph
Burnett
1 2 3 4 5 6 0.00 0.05 0.10 0.15 0.20 0.25 0.30 1w
Burnett
1 2 3 4 5 6 0.4 0.5 0.6 0.7 0.8 0.9 1.0 1w 5 3 vph
super -Burnett
1 2 3 4 5 6 0.00 0.05 0.10 0.15 0.20 0.25 0.30 1w
super -Burnett
mode with lowest damping only; unstable modes not shown . . .
SLIDE 27 Example: Dispersion and Damping of Ultrasound
left: inverse phase velocity q
5 3RT/vph
dashed: Maxwell molecules right: reduced damping α/ω continuous: hard spheres Grad13 / R13
1 2 3 4 5 6 0.4 0.5 0.6 0.7 0.8 0.9 1.0 1w 5 3 vph
G13
1 2 3 4 5 6 0.00 0.05 0.10 0.15 0.20 0.25 0.30 1w
G13
1 2 3 4 5 6 0.4 0.5 0.6 0.7 0.8 0.9 1.0 1w 5 3 vph
R13
1 2 3 4 5 6 0.00 0.05 0.10 0.15 0.20 0.25 0.30 1w
R13
mode with lowest damping shown; all modes stable . . .
SLIDE 28
Generalized 13 moment eqs: nonlinear [HS 2004]
2nd order for arbitrary interaction potentials
derived by O of M, matched to Burnett by CE expansion Dqi Dt + qk ∂vi ∂xk + 5 3qi ∂vk ∂xk + a(1,1)θ∂σik ∂xk − a(1,5)θσik ∂ ln p ∂xk − a(1,6)σikqk μ = −a(1,0)p μ ∙ qi + 5 2 Prμ ∂θ ∂xi ¸ Dσij Dt +2σkhi ∂vji ∂xk +σij ∂vk ∂xk +a(2,1) ∂qhi ∂xji +a(2,5)qhi ∂ ln p ∂xji +a(2,6)σkhiσjik μ +a(2,7)qhiqji θμ = −a(2,0)p μ ∙ σij + 2μ∂vhi ∂xji ¸
coefficients a(1,1) = 5 4 1 Pr θ4 θ2 , a(1,5) = −5 4 1 Pr θ3 θ2 , a(1,6) = 3 2 θ5 θ2 − 1 − θ4 2θ2 ω a(2,1) = 4 5 Pr $3 $2 , a(2,5) = 4 5 Pr $4 $2 , a(2,6) = 2 − $6 2$2 , a(2,7) = 8 Pr2 25 ∙$3 $2 ω − $5 $2 ¸ Burnett coefficients [Reinecke&Kremer 1990’s] γ ω $2 $3 $4 $5 $6 θ2 θ3 θ4 θ5 MM 5 1 2 3 3 8 5.625 −3 3 9.75 HS ∞ 0.5 2.028 2.418 0.681 0.219 7.424 5.822 −3.09 2.418 8.286
SLIDE 29 R13 nonlinear [MT&HS in prep]
non-linear 2nd order with linear 3rd order contributions
Dρ Dt + ρ∂vk ∂xk = 0 ρDvi Dt + ∂ρθ ∂xi + ∂σik ∂xk = ρGi 3 2ρDθ Dt + ∂qk ∂xk + p∂vk ∂xk + σij ∂vi ∂xj = 0 Dqi Dt + qk ∂vi ∂xk + 5 3qi ∂vk ∂xk + a(1,1)θ∂σik ∂xk − a(1,2) ∂ ∂xk µ θ ∙ σik + 2μ ∂vhi ∂xki ¸¶ − a(1,3) ∂ ∂xk µμ ρ ∂qhi ∂xki ¶ − a(1,4) ∂ ∂xi µμ ρ ∂qk ∂xk ¶ − a(1,5)θσik ∂ ln p ∂xk − a(1,6)σikqk μ = −a(1,0) p μ ∙ qi + 5 2 Prμ ∂θ ∂xi ¸ Dσij Dt + 2σkhi ∂vji ∂xk + σij ∂vk ∂xk + a(2,1) ∂qhi ∂xji − a(2,2) ∂ ∂xhi ∙ qji + 5 2 Prμ ∂θ ∂xji ¸ − a(2,3) ∂ ∂xk µμ ρ ∂σhij ∂xki ¶ − a(2,4) ∂ ∂xk µμ ρ ∂σij ∂xk ¶ + a(2,5)qhi ∂ ln p ∂xji + a(2,6)σkhiσjik μ + a(2,7)qhiqji θμ = −a(2,0) p μ ∙ σij + 2μ ∂vhi ∂xji ¸
a(1,0) a(1,1) a(1,2) a(1,3) a(1,4) a(1,5) a(1,6) Pr MM 0.66667 1.0 0.0 2.4 2.0 1.0 1.33333 0.66667 HS 0.65006 0.78694 0.17693 2.09248 1.50489 1.00504 0.98678 0.66085 a(2,0) a(2,1) a(2,2) a(2,3) a(2,4) a(2,5) a(2,6) a(2,7) MM 1.0 0.8 0.0 2.0 0.0 0.0 0.0 0.0 HS 0.98632 0.63125 0.09466 2.19368 0.11447 0.231279 0.35548 0.10270
SLIDE 30 Boundary conditions
from Maxwell BC and Grad closure Z
R3
Ψ (C) fR13 (C) dC = 2χ 2 − χ
∞
Z Z
R2
h Ψ (CW − V) fwall (CW) − Ψ (C) f(even)
R13
(C) i dCtdCn BC’s for moments rπ 2θ2 − χ χ σtn = b(1,0)PVt + b(1,1)qt − μ ρ µ b(1,2)∂σhtn ∂xni + b(1,3)∂σtk ∂xk ¶ + b(1,4)∆qt rπ 2θ2 − χ χ qn = b(2,0)P∆θ + b(2,1)θσnn − μ ρ µ b(2,2) ∂qhn ∂xni + b(2,3) ∂qk ∂xk ¶ − b(2,4)θ∆σnn − b(2,5)PV 2
t
rπ 2θ2 − χ χ µμ ρ ∂qht ∂xni + b(3,7)θ∆σtn ¶ = b(3,0)θPVt−b(3,1)θqt+μ ρθ µ b(3,2)∂σhtn ∂xni − b(3,3)∂σtk ∂xk ¶ +b(3,4)θ∆qt−b(3,5)P∆θVt−b(3,6)PV 3
t
rπ 2θ2 − χ χ μ ρ ∂σhnn ∂xni = b(4,0)P∆θ − b(4,1)θσnn + μ ρ µ b(4,2) ∂qhn ∂xni − b(4,3) ∂qk ∂xk ¶ + b(4,4)θ∆σnn − b(4,5)PV 2
t
rπ 2θ2 − χ χ μ ρ ∂σhns ∂xsi = −b(5,1)θσss + b(5,2)μ ρ ∂qhs ∂xsi + b(5,3)θ∆σss + b(5,4)PV 2
t
where: Vt = vt − v(W)
t
, ∆θ = θ − θW ∆σij = σij − σ(NSF)
ij
= σij + 2ε ∂vhi ∂xji ∆qi = qi − q(NSF)
i
= qi + ε 5 2 Pr ∂θ ∂xi P = ρW p θWθ = b(0,0)ρθ + b(0,1)σnn + μ ρ µ b(0,2) ∂qhn ∂xni + b(0,3) ∂qk ∂xk ¶ + b(0,4)θ∆σnn
SLIDE 31
Boundary conditions
more coefficients
MM HS b(0,0) 1.0 1.0 b(0,1) 0.5 0.51710 b(0,2) 0.17143 0.18414 b(0,3) 0.1 0.09829 b(0,4) 0.0 0.01533 MM HS b(1,0) 1.0 1.0 b(1,1) 0.2 0.20521 b(1,2) 1.0 1.03061 b(1,3) 0.0 0.00068 b(1,4) 0.0 0.00352 MM HS b(2,0) 2.0 2.0 b(2,1) 0.5 0.41802 b(2,2) 0.85714 0.87743 b(2,3) 0.8 0.73977 b(2,4) 0.0 0.07219 b(2,5) 0.5 0.5 MM HS b(3,0) 0.20833 0.12539 b(3,1) 0.45833 0.47552 b(3,2) 0.20833 0.09633 b(3,3) 0.0 0.00299 b(3,4) 0.0 0.01784 b(3,5) 0.20833 0.12539 b(3,6) 0.20833 0.21976 b(3,7) 0.0 0.08057 MM HS b(4,0) 0.2 0.20770 b(4,1) 0.7 0.71067 b(4,2) 0.17143 0.17326 b(4,3) 0.08 0.07683 b(4,4) 0.0 0.01406 b(4,5) 0.3 0.31155 MM HS b(5,1) 0.5 0.50272 b(5,2) 0.17143 0.17626 b(5,3) 0.0 0.01437 b(5,4) 0.25 0.259628
SLIDE 32 Example: Couette flow ε = 0.1, χ = 1
preliminary results
solution of semi-linear R13, similar to [PT,MT&HS 2009] shear stress σxy = −C0 Kn velocity vx (y) = −a(2,0) a(2,0) 1 KnC1 y − a(2,1) − a(2,2) 2a(2,0) qx (y) heat flux parallel to wall qx(y) = −1 + a(1,6) a(1,0) Kn2 C3
0 y + C2 sinh
⎡ ⎣ s 2a(1,0)a(2,0) a(1,3)a(2,0) − a(1,2)a(2,1) + a(1,2)a(2,2) y Kn ⎤ ⎦ shear stress velocity heat flux parallel to wall
constants of integration Cα from BC
SLIDE 33 Example: Couette flow ε = 0.1, χ = 1
preliminary results
solution of semi-linear R13, similar to [PT,MT&HS 2009] temperature θ (y) = C4 − Pr 5 C2
0 y2 − 2 Pr a(1,1) − a(1,2)
5a(1,0) σyy pressure p = P0 − σyy (y) normal stress
σyy (y) = −2 + 2a(2,1) + a(2,6) 3a(2,0) Kn2 C2
0−2a(2,7)
3a(2,0) Kn2 C4
0 y2+C3 cosh
Ãs 15a(1,0)a(2,0) 9a(2,3) + 15a(2,4) + 10a(2,2)(a(1,2) − a(1,1)) y Kn ! temperature pressure normal stress constants of integration Cα from BC
SLIDE 34
Shocks: Comparison with DSMC results [MT & HS 2004]
Failure of NSF, Burnett, super-Burnett, Grad13
SLIDE 35
Shocks: Comparison with DSMC results [MT & HS 2004]
Success of R13
SLIDE 36 Lid-driven cavity flow (MM) [A.Rana, MT, HS 2013]
velocity streamlines and stress contours Kn = 0.08, vlid = 50m
s
SLIDE 37
Lid-driven cavity flow (MM) [A.Rana, MT, HS 2013]
velocity on horizontal center line temperature on vertical center line
SLIDE 38 Lid-driven cavity flow (MM) [A.Rana, MT, HS 2013]
temperature contours and heat flux streamlines Kn = 0.08, vlid = 50m
s
DSMC, R13: heat flux from hot to cold!
SLIDE 39
2D Bottom heated plate (MM) [AR, AM, HS submitted]
comparison of NSF, DSMC, R13 Kn=0.05: heat flux and shear stress Kn=0.05: temperature and streamlines
SLIDE 40
2D Bottom heated plate (MM) [AR, AM, HS submitted]
comparison of NSF, DSMC, R13: temperature and streamlines Kn=0.13 Kn=0.3
SLIDE 41
2D Bottom heated plate (MM) [AR, AM, HS submitted]
Normal and total heat flux at bottom plate Normal heat flux at center line
SLIDE 42 Summary and outlook
- regularized 13 moment eqs for non-Maxwellian molecules
- order of magnitude reduction of large Grad-type moment system
- proper moments to describe collective behavior (depend on C(n)
ab )
- coefficients for hard sphere molecules and Maxwell molecules
- first time super-Burnett coefficients for hard sphere molecules
- other interaction potentials: compute production matrices C(n)
ab
- extension to non-linear equations
- boundary conditions for moments from Grad distribution
- visible differences between HS/MM in rarefaction effects
Outlook
- transpiration flow, heat transfer, . . .
- 2-D / 3-D flows
- other interaction potentials
- polyatomic, mixtures, . . .
