Exact linear hydrodynamics from the Boltzmann equation

0.15 MB·English

by  I. V. Karlin

Exactlinearhydrodynamics from the Boltzmannequation
I.V. Karlin,1,∗ M. Colangeli,2 and M. Kro¨ger2,†
1AerothermochemistryandCombustionSystemsLab, ETHZu¨rich, CH-8092Zu¨rich, Switzerlandand
School ofEngineeringSciences, UniversityofSouthampton, SO171BJ,Southampton, UnitedKingdom
2PolymerPhysics,DepartmentofMaterialsandMaterialsResearchCenter,ETHZu¨rich,CH-8093Zu¨rich,Switzerland
(Dated: 2008-01-1809:47:09)
Exact(toallordersinKnudsennumber)equationsoflinearhydrodynamicsarederivedfromtheBoltzmann
8 kineticequationwiththeBhatnagar-Gross-Krookcollisionintegral.Theexacthydrodynamicequationsarecast
0 inaformwhichallowsustoimmediatelyprovetheirhyperbolicity,stability,andexistenceofanH-theorem.
0
2
PACSnumbers: 51.10.+y(Kinetictheory)05.20.Dd(Kinetictheory)

n Keywords:Kinetictheory;heattransfer;hydrodynamics;hyperbolicequations;Boltzmannequation
a
J
Hydrodynamicsassumesthatastateofafluidissolelyde-
8 0
1 scribed by five fields: density, momentum and temperature.
DerivationoftheNavier-Stokes-Fourier(NSF)hydrodynamic −0.2
Re(ω )
] equations from the Boltzmann kinetic equation as the first- ac
h Im(ω ) / 10
c order approximation in the Knudsen number (ratio between −0.4 Re(ωac) hydrodynamic modes
e mean free path and a flow scale) by Enskog and Chapman −0.6 Im(ωaacc) / 10 (real and imaginary parts)
m is a textbook example of a success of statistical physics [1]. Re(ω )
diff
- Recentrenewedinteresttotheproblemsbeyondthestandard −0.8 Im(ωdiff)
t hydrodynamicsisdue,inparticular,toflowsimulationandex- Re(ω )
a shear
t perimentsatamicro-andnano-scale[2–4]. However,almost −1 Im(ωshear)
.s acenturyofefforttoextendthehydrodynamicdescriptionbe- 0 0.5 1 k 1.5 2 2.5
t
a yondtheNSFapproximationfailedeveninthecaseofsmall
m deviationsaroundthe equilibrium. Inordertoappreciatethe
FIG. 1: Exact hydrodynamic modes ω of the Boltzmann-BGK kinetic
problem, let us remind that, in the NSF approximations, the
- equation asafunction ofwavenumberk (twocomplex-conjugated acous-
d decay rate of the hydrodynamic modes is quadratic in the
ticmodesωac,twicedegeneratedshearmodeωshearandthermaldiffusion
on wavevector,Re(ω) ∼ −k2,andisunbounded. Ontheother modeωdiff).Thenon-positivedecayratesRe(ω)attainthelimitofcollision
hand, Boltzmann’s collision term features equilibration with frequency(−1)ask→∞.
c
finite characteristicrates. This“finite collisionfrequency”is
[
obviously incompatible with the arbitrary decay rates in the

1 NSFapproximation:intuitively,hydrodynamicmodesatlarge
v popular in applications [9], and features a single relaxation
k cannot relax faster than the collision frequency. Now, the
2 rate. Theresultforthehydrodynamicmodesisdemonstrated
classicalmethodofEnskogandChapmanextendsthehydro-
3 inFig.1. ItisclearfromFig.1thattherelaxationofnoneof
dynamicsbeyondtheNSFinsuchawaythatthedecayrateof
9
thehydrodynamicmodesisfasterthanω = 1whichisthe
2 thenextorderapproximations(Burnettandsuper-Burnett)are −
collision frequencyin the units adopted in this paper. Thus,
. polynomialsofhigherorderink. Insuchanextension,relax-
1 theresultfortheexacthydrodynamicsindeedcorrespondsto
ation rate may becomecompletelyunphysical(amplification
0 theaboveintuitivepicture.Below,weapplythemethodofin-
8 instead of attenuation), as first shown by Bobylev [5] for a
variantmanifold [10] to derive the hydrodynamicequations.
0 particularcaseofMaxwellmolecules.Thisindicatesinability
Thenon-perturbativederivationismadepossiblewithanop-
: oftheChapman-Enskogmethodtotackletheaboveproblem,
v timalcombinationof analyticalandnumericalapproachesto
andnon-perturbativeapproachesaresought. Theproblemof
i
X exact hydrodynamics has been studied in depth recently for solvetheinvarianceequation.
PointofdepartureisthelinearizedBoltzmann-BGKequa-
r toy(finite-dimensional)models-momentsystemsofGrad-in
a tionforthedeviation f =f fGMofthedistributionfunc-
[6–8],andmanyremarkableresultswereobtained. Inpartic- tionf fromaglobalM△axwelli−anfGM(c2) = π−3/2e−c2. In
ular,in[6,7]itwasshownthattheexacthydrodynamicequa-
thereciprocalspace,itreads,
tions are hyperbolicand stable for all wave numbers. How-
ever,for“true”kineticequationssuchquestionsremainopen. ∂ f = ik c f δf; δf =f fLM, (1)
t
△ − · △ − −
InthisLetterwederiveexacthydrodynamicequationsfrom
withthewavevectork = ke defininge ,k k,peculiar
thelinearizedBoltzmannequationwiththeBhatnagar-Gross- k k ≡ | |
velocity c and time t. All quantities are considered dimen-
Krook (BGK) collision term. This kinetic equation remains
sionless, i.e.,reducedwiththeunitsoftherelaxationtimeτ,
the thermalvelocity 2kBT/mand mass m of the particle.
In(1),thelinearizedlpocalMaxwellianisfLM =fGM(1+ϕ0)
∗†UElRecLt:rownwicw.adcdormespsl:[email protected] ewrhagerees1ar+edϕe0fi=nehd1fiofr+ar2bcitr·ahrcyifω+vi32a(cω2−f 23)hc2ω−f32di3fc.,Aanvd-
h i ≡ R
2
weintroducepertinentquantitieswhichcharacterizedeviation σ1k σ2k σ3k σ4
fromtheglobalequilibrium:n 1 1(densityperturba-
tion),u ≡ hcif (velocitypertur≡bahtioinf)−andT ≡ 32hc2− 23if hλ−kkδ2XB1i hλikkδAX2i hλ−kkδ2XC3i hckick⊥DδY4i
(temperatureperturbation). Sincethescalarproductbetween
kandcappearsin(1),thedistributionfunctionofferssymme- B0 =−31 A0 =−23 C0 =0 D0 =−21
trywithrespecttothee -axis,whichisnotuniaxialincaseu real,⊕ imag,⊕ real,⊕ imag,⊖
k
is not collinear with ek. We denote the two components of q1k q2k q3k q4
the mean velocity as uk = u·ek and u⊥ = u·e⊥, where hγkδX1i hγkδX2i hγkckδX3i h(c2− 25)c⊥δY4i
theunitvectore⊥belongstotheintersectionoftheplaneper- ikX −k2Z ikY −k2U
pue=ndiuckuelar+tou⊥kean.dEqthueatpiolannseofspcahnannegdefboyrmkoamndenuts, sωo tharaet X0 =0 Z0 =−16 Y0 =−45 U0 = 21
obtainedkbyinteg⊥rationoftheweighted(1)overd3cash i imag,⊖ real,⊖ imag,⊖ real,⊕
∂ ω = ik cω ω . (2) TABLEI:Symmetryadaptedcomponentsof(nonequilibrium)stress
th i△f − ·h i△f −h iδf tensor σ andheat fluxq, introducedin(7a)and(7b), respectively.
In order to calculate such averages, we can switch to spher- Row2:Microscopicexpressionofthesecomponents(averagingwith
ical coordinates. For each (at present arbitrary) wave vec- theglobalMaxwellian).Short-handnotationused:λk =c2−c2 and
k 3
tor, we choose the coordinate system in such a way that its γk = (c2− 5)c . Row3: Expressionofthecomponents interms
2 k
z-direction aligns with e . We can then express c in terms of(asweshow,real-valued)functionsA–Z (seetext). Row4: Val-
k
of its norm c, a vertical variable z and plane vector e (az- uesoffunctionsA–Zatk=0.Thesevaluesrecoverhydrodynamic
φ
imuthal angle e e = cosφ) for the present purpose as equationsuptoBurnettapproximation. Row5: Paritywithrespect
c/c = √1 z2φe· +⊥ze . We could have equally chosena toz–symmetric(⊕)orantisymmetric(⊖)–ofthepartofthecor-
− φ k responding δX entering the averaging in row 2, and whether this
fixedcoordinatesystemintheplaneorthogonaltok,andtwo
partisimaginaryorreal-valued(seeFig.2). Row3isanimmediate
fieldsinsteadofu⊥plusanangle,viz.u⊥+iu⊥ =eiφu⊥.Due
x y consequenceofrow5.
to isotropy, u⊥ alone fully represents the twice degenerated
(shear)dynamics. Inordertosimplifynotationandcompute
thedynamicsofallfivefieldsweintroduceafour-dimensional whichisanonlinearintegralequationfortheunknownfields
v(sxiemckt,poxlre4)foowfrmihth,yϕdxr0ko=d≡yXn(an0m,uixck.,fiTTeh)ledasvn,edcxtxo4r≡X=0u(c⊥xa.n1T,ixmh2em,nxeϕ3d0,iaxtta4ek)lyesb=ae δoHXfe(r,2e,)b,eωfcoairussaseismu∂ititlxaabrhlwyasicthhtooXsbe0enarvenepdcltadocireffωderbsfyuflrfitohlmleinrXgigh0hωtmih△aainfndly=sbidxee-.
readoff,wehaveX0(c·,z)=(1,2c ,c2 3,2c ),wherewe causeofconventionsforprefactorsinthetemperatureandve-
introduced,forlateruse,theabbrevikation−s 2 φ locitydefinitions,ω =(1,c ,2(c2 3),c ). Withinthesame
k 3 −2 φ
c eigen-closure,Eq.(2)islinearinxandhencewrittenas
c c e , c c e , c φ , (3)
k ≡ · k φ ≡ · ⊥ ⊥ ≡ e⊥·eφ ∂tx=M x. (6)
·
suchthatik c = ikc , c = c e +c e with c = cz and
· k ⊥ φ k k k The matrix M solely depends on the non-hydrodynamic
c⊥ = c√1−z2, contrasted by cφ (and eφ), do not depend fields, the heat flux q c(c2 52) f and the stress tensor
ontheazimuthalangle. Similarly,weintroduceyetunknown ≡ h − i
σ cc , where s denotesthe symmetric tracelesspart
tfiheelddsisδtXrib(uct,ikon)wfuhniccthiocnh,aδraϕct=erizδeft/hfeGnMonienquteirlimbrsiuomf tphaerthoyf- of≡a tehnsorisf[6, 11], s = 12(s+sT)− 13tr(s)I. Using(4),
thestresstensorandheatfluxuniquelydecomposeasfollows
drodynamicfieldsxthemselves,
δϕ=δX x=δX1n+δX2uk+δX3T +δX4u⊥, (4) 3
· σ = σk e e +σ⊥2e e , (7a)
Twhhiesre“eδiXge4n”f-accltoosruizrees(a4s) δwXhi4c(hc,fzo,rφm)al=ly a2nδdY4v(ecr,yz)geenφe·rael⊥ly. q = qke2 +kqk⊥e , k ⊥ (7b)
k ⊥
addressesthefact,thatwewishtonotincludeotherthanhy-
drodynamicvariablesimpliesa closurebetweenmomentsof withthemomentsσk =(σ1k,σ2k,σ3k)·xkandσ⊥ =σ4x4,and
similarlyforq(seeRow2ofTab.I).Theprefactorsarisefrom
thedistributionfunction,tobeworkedoutindetailbelow. It theidentities e e : e e = 2 and e e : e e = 1. The
assumes the existence of an invariant manifold, and the hy- k k k k 3 k k k ⊥ 2
appearanceof δY4 rather than δX4 in the expression for the
drodynamic fields as slow variables which leave the higher
orthogonalmoment(inTab.I)reflectsthefactthatwehaveal-
moments “slaved”. In order to ensure these contributionsto readyintegratedouttheangularvariable, 2πe e e dφ=
notinterferewiththelocalMaxwellian,onehasthefreedom R0 φ φ· ⊥
torequire X0 δf = 0withoutproducinganylimitation,i.e., πe⊥. We note in passing that, while the stress tensorhas, in
h i general, threedifferenteigenvalues,in the presentsymmetry
bykeepingxtobedefinedthroughthelocalMaxwellianpart
adaptedcoordinatesystemitexhibitsavanishingfirstnormal
of the distribution function. Using the above form for δf in
(1), and using the canonicalabbreviation X X0 +δX, stressdifference. Sincetheintegralkernelsofallmomentsin
△ ≡ (7) do not depend on the azimuthalangle, these are actually
yields
two-dimensionalintegrals over c [0, ] and z [ 1,1],
X ∂ x= ikc X x δX x, (5) weightedby2πc2fGM(c2)δX . ∈ ∞ ∈ −
△ · t − k△ · − · µ
3
Stresstensorandheatfluxcanyetbewritteninanalterna- [6,7]andbelow),oncethefunctionsA–Z areexplicitlyeval-
tive form which is defined by Row 3 of Tab. I. As we will uated. For that, we still requireδX. Combining(5)and (6),
seelateron,duetobasicsymmetryconsiderations,thehereby andrequiringthatthe resultholdsforanyx(invariancecon-
introducedfunctions A–Z are real-valued. We postpone the dition),weobtainaclosed,singularintegralequation(invari-
related proof, and proceed by using these functions A–Z to anceequation)forcomplex-valuedδX,
splitMintopartsasM=Re(M) iIm(M)with
− δX=X0 (M+[ikc +1]I)−1 X0. (9)
· k −
0 0 0 0
  Notice that δX vanishes for k = 0, and that (9) is supple-
Re(M) = k2 320X A0 230Y 00 , (8) mvaennistehdinwgitLhagthraenbgaesimcuclotinpsltirearisnt(mhXat0riixδf) =X00δ,Xor,ewquhailclhy,,
  h i
 0 0 0 D  however, is automatically dealt with if we only evaluate
anisotropic (irreducible) moments with δf, such as those
0 1 0 0
  listed in Tab. I. The implicit equation (9) is identical with
Im(M) = k 12−0k2B 32(1−0k2Z) 12−0k2C 00 . trheseuelti,gwenit-hclMosufrreo(m4)(,8a)n,dcki,sXo0u,ramndaiAn–aZnddperfiancetidcainllyanudsejufustl
 0 0 0 0  before(3)andTab.I,respectively.
Weiterativelycalculate(i)δXdirectlyfrom(9)foreachk
intermsofM,(ii)subsequentlycalculatemomentsfromδX
bynumericalintegration. Importantly,thefix pointoftheit-
eration(i)-(ii)-(i)-..isuniqueforeachk,i.e.,doesnotdepend
ontheinitialvaluesformomentsA–Z,aslongaswechoose
real-valued ones which are consistent with (9), as we prove
in the next paragraph. In addition, two other computational
strategies were implemented: First, we used continuation of
functionsA–Z fromtheirvaluesatk =0tosolve(9)withan
incrementalincreaseofk,wherethesolutionatkwasusedas
theinitialguessfork+dk. Second,weusedalsoacontinu-
ation“backwards”inwhichthesolutionatsomek (obtained
byconvergentiterationswitharandominitialcondition)was
usedastheinitialguessfora solutionatk dk. Both these
−
strategiesreturnedthesamevaluesoffunctionsA–Z ascom-
putedbyiterationsfromarbitraryinitialcondition. Thesolu-
tionδXallowstocalculatethewholedistributionfunctionf
via (4) as illustrated by Fig. 2. For resulting moments for a
widerangeofk-valuesseeFig.3.
Finally, we need to clarify the origin of row 5 in Tab. I
(whichisdirectlyillustratedbyFig.2)anditsimplicationon
thestructureofM(8)whoseentriesare–apriori–complex-
valued functions to be calculated with complex-valued δX.
Wewishtomakeuseofthefactthatallintegralsoverz van-
ish for odd integrands. To this end we introduce abbrevia-
tions ( ) for a real-valued quantity which is even (odd)
⊕ ⊖
with respect to the transformation z z. One notices
X0 = ( , , , ), and we recall th→at A−–Z are integrals
⊕ ⊖ ⊕ ⊕
overeitherevenoroddfunctionsinz,timesa componentof
δX (see Tab. I). Let us prove the consistency of the spec-
ified symmetry of M and the invariance condition: Start by
assumingA–Z to be real-valuedfunctions. ThenM =
µν
⊕
if µ + ν is even, and M = i otherwise. This implies
FIG.2: (Coloronline)Sampledistributionfunctionf(c,k)atk=1,fully µν ⊕
characterizedbythefourquantitiesδX1,2,3(c,z)andδY4(c,z).Shownhere δX1 = ⊕ + i⊖, δX2 = ⊖ + i⊕, δX3 = ⊕ + i⊖, and
areboththeirreal(left)andimaginaryparts(rightcolumn). Inordertoim- δX4 = +i , i.e., different symmetry properties for real
⊕ ⊖
provecontrast,weactuallyplotln|1+fGMδXµ|multipliedbythesignof andimaginaryparts. Withthese“symmetry”expressionsfor
δXµ.Samecolorcodeforallplots,rangingfrom−0.2(red)to+0.2(blue). X0,δX,andMathand,wecaninsertintotherighthandside
of the equation, δX = (X0 +δX) (M+i I), which is
· ⊖
identicalwiththeinvarianceequation(9). Thereareonlytwo
Note that the checkerboard structure of the matrix M (8) cases to consider, because M has a checkerboard structure,
is particularlyusefulfor studyingpropertiesof the hydrody- i.e.,onlytwotypesofcolumns: Columnsµ = 1andµ = 3:
namic equations(6), such as hyperbolicityand stability (see δXµ = +i becauseM1−3,4 = 0; Columnsµ 2,4 :
⊕ ⊖ ∈ { }
4
δXµ = +i ifMµ,1−3 =0(whichisthecaseforcolumn ∂tx′ =M′ x′,withM′ =T M T−1ismanifestlyhyper-
4)and ⊕+i ⊖ifMµ,4 =0(whichisthecaseforcolumn2). bolicandsta·ble;Im(M′)issym· me·tric,Re(M′)issymmetric
⊖ ⊕
and non-positive semi-definite. The corresponding transfor-
mationmatrixTcanbeeasilyreadofftheresultsobtainedin
0.5 [6,7]forGrad’ssystemssincethestructureofthematrixM
(8) is identical to the one studied in [6, 7]. We have explic-
itlyverifiedthatmatrixT(equations(21)–(23)inRef.[7]and
0
Z
o (13)inRef.[6])withthefunctionsA–Zderivedhereinisreal-
A t valued and thus render the transformedhydrodynamicequa-
nts −0.5 A tionsmanifestlyhyperbolicandstable.Wenotethatthisresult
e B
m – hyperbolicity of exact hydrodynamic equations – strongly
C
o
m D supportsarecentsuggestionbyBobylevtoconsiderahyper-
−1 U
X bolicregularizationof the Burnettapproximation[12]. Sim-
Y ilarly, using the hyperbolicity, an H-theorem is elementary
Z
−1.5 provenasin[6,12]. Finally,usingtheaccuratedataforfunc-
0 0.5 1 1.5 2 2.5 3
k tionsA–Z, wecanwriteanalyticapproximationsforthehy-
drodynamicequations(6)insuchawaythathyperbolicityand
FIG.3:MomentsA–Zvs.wavenumberkobtainedwiththesolutionof(9). stabilityisnotdestroyedinsuchanapproximation(see[7]).
In conclusion, we derived exact hydrodynamic equations
Wehavethusshownthatbothsidesoftheinvarianceequa- fromthelinearizedBoltzmann-BGKequation.Themainnov-
tion(9)haveequalsymmetryproperties,andthatδXwiththe elty is the numericalnon-perturbativeprocedureto solve the
specified symmetriesis consistentwith real-valuedmoments invariance equation. In turn, the highly efficient numerical
A–Z. The proof implies, that any iteratively obtained solu- approach is made possible by choosing a convenient coor-
tion, if it exists, starting with arbitrary real-valued moments dinate system and establishing symmetries of the invariance
A–Z in (9) to evaluate δX must convergeto real-valuedso- equation. Theinvariantmanifoldin the spaceof distribution
lution A–Z. Since the solution is smoothly varying with k, functionsistherebycompletelycharacterized,thatis,notonly
and since A–Z at k = 0 are known and are real-valued, the equationsof hydrodynamicsare obtained but also the corre-
momentsmustbereal-valuedoverthewholek-space. sponding distribution function is made available. The perti-
WiththeresultforthefunctionsA–Z athand,theextended nentdatacanbeused,inparticular,asamuchneededbench-
hydrodynamic equations are closed. Let us briefly discuss markforcomputation-orientedkinetictheoriessuchaslattice
the pertinent properties of this system. First, the general- Boltzmannmodels, as well as forconstructingnovelmodels
ized transport coefficients are given by the nontrivial eigen- usingquadraturesin the velocityspace [13, 14]. Finally, we
values of k−2Re(M): λ2 = A (elongation viscosity), have established a novel non-perturbativecomputationalap-
λ3 = −32Y−(thermal diffusivity),−and λ4 = −D (shear vis- proachtofindinginvariantmanifoldsofkineticequations.
cosity). All these generalized transportcoefficients are non-
negative (see Fig. 3). Second, computing the eigen-values The present approach can be extended to the Boltzmann
of matrix M we obtain the dispersion relation ω(k) of the equationwith othercollision terms. Theabovederivationof
correspondinghydrodynamicmodesalreadypresentedinFig. hydrodynamicsisdoneunderthe standardassumptionof lo-
1. Third, a suitable transform of the hydrodynamic fields, calequilibrium[1], howeverthe assumptionitself is open to
x′ = T x, where T is a real-valued matrix, can be es- furtherstudy[15][WethankH.C.O¨ttingerforthisimportant
·
tablishedsuchthatthetransformedhydrodynamicequations, remark].I.V.K.acknowledgessupportofCCEM-CH.
[1] S. Chapman, T. G. Cowling, The Mathematical Theory of [8] A.N.Gorban,I.V.Karlin,Phys.Rev.Lett.77(1996)282.
NonuniformGases (CambridgeUniv.Press,NewYork,1970). [9] C.Cercignani,TheoryandApplicationoftheBoltzmannEqua-
[2] A.Beskok, G.E.Karniadakis, Microflows: Fundamentalsand tion (ScottishAcademicPress,Edinburgh,1975).
Simulation (Springer,Berlin,2001). [10] A.N.Gorban,I.V.Karlin,Transp.Th.Stat.Phys.23(1994)559.
[3] S. Ansumali, I.V. Karlin, S. Arcidiacono, A. Abbas, N.I. [11] M. Kro¨ger, Models for Polymeric and Anisotropic Liquids
Prasianakis,Phys.Rev.Lett.98(2007)124502. (Springer,Berlin,2005).
[4] D.M.Karabacak,V.Yakhot,K.L.Ekinci,Phys.Rev.Lett.98 [12] A.V.Bobylev,J.Stat.Phys.124(2006)371.
(2007)254505. [13] S.S. Chikatamarla, I.V. Karlin, Phys. Rev. Lett. 97 (2006)
[5] A.V.Bobylev,Sov.Phys.Dokl.27(1982)29. 190601.
[6] M.Colangeli, I.V.Karlin, M.Kro¨ger, Phys.Rev.E76(2007) [14] X.Shan,X.He,Phys.Rev.Lett.80(1998)65.
022201. [15] H.C. O¨ttinger, H. Struchtrup, Multiscale Model. Simul. 6
[7] M.Colangeli, I.V.Karlin, M.Kro¨ger, Phys.Rev.E75(2007) (2007)53.
051204.