Appendix69uuu.com

ThermallatticeBoltzmannequationforlowMachnumberows:DecouplingmodelZhaoliGuo,*ChuguangZheng,andBaochangShiNationalLaboratoryofCoalCombustion,HuazhongUniversityofScienceandTechnology,Wuhan,People'sRepublicofChinaT.
S.
ZhaoDepartmentofMechanicalEngineering,HongKongUniversityofScienceandTechnology,Kowloon,HongKongReceived4November2006;published14March2007AlatticeBoltzmannmodelisproposedforsolvinglowMachnumberthermalowswithviscousdissipationandcompressionworkinthedouble-distribution-functionframework.
Adistributionfunctionrepresentingthetotalenergyisdenedbasedonasinglevelocitydistributionfunction,anditsevolutionequationisderivedfromthecontinuousBoltzmannequation.
AlatticeBoltzmannequationmodelwithclearphysicsandasimplestructureisthenobtainedfromakineticmodelforthedecoupledhydrodynamicandenergyequations.
ThemodelistestedbysimulatingathermalPoiseuilleowandnaturalconvectioninasquarecavity,anditisfoundthatthenumericalresultsagreewellwiththeanalyticalsolutionsand/orthedatareportedinpreviousstudies.
DOI:10.
1103/PhysRevE.
75.
036704PACSnumbers:47.
11.
j,44.
05.
eI.
INTRODUCTIONAlthoughthelatticeBoltzmannequationLBEmethodhasachievedgreatsuccessinsimulatingathermalandiso-thermaluidows,itsapplicationsforthermohydrodynam-icsisstillnotsatisfactory.
ConstructingLBEmodelsforther-malowsremainschallengingintheLBEcommunity,althoughsomeeffortshavebeenmadefromvariousview-points.
Arecentcomprehensivereviewonthistopiccanbefoundelsewhere1.
TheexistingstrategiesforconstructingthermalLBETLBEmodelscanbeclassiedintothreecategories,i.
e,themultispeedapproach,thedouble-distribution-functionDDFapproach,andthehybridapproach.
ThemultispeedapproachisastraightforwardextensionoftheathermalLBE,inwhichonlythevelocitydistributionfunctionVDFisused3–5;theDDFapproachutilizestwodifferentdistribu-tionfunctionsDFs,oneforthevelocityeldandtheotherforthetemperatureorinternalenergyeld;thehybridap-proachissimilartotheDDFapproachexceptthattheenergyequationissolvedbydifferentnumericalmethodse.
g.
,nite-differenceornite-volumemethodsratherthanbysolvingtheLBE1.
BoththemultispeedandDDFapproacheshavesomelimitations2.
Themultispeedmodelsusuallysufferfromseverenumericalinstabilityandthetemperaturevariationsimulatedislimitedtoanarrowrange,andusuallyresultsinaxedPrandtlnumber,althoughsomelaterversionshaveovercomethisproblem6.
ForDDFLBEmodels,althoughtheyexhibitgoodnumericalstabilityandanadjustablePrandtlnumber7–23,mostofthemexceptforthosepro-posedinRefs.
9,13takenoaccountoftheviscousdissi-pationandcompressionwork.
TherstDDFmodelthatincorporatesviscousdissipationandcompressionworkisattributedtoHe,Chen,andDoolenHCD9,whereanadditionalDFisdenedfortheuidtemperatureandisderiveddirectlyfromthemomentoftheVDF.
Thismodelhasattractedmuchattentionsinceitsemer-gence,andhasfoundapplicationsinavarietyofelds15–23.
DespitetheapparentadvantagesoftheHCDTLBEmodel,itiswellrecognizedthatthismodelstillsuffersfromsomedeciencies.
Forinstance,boththeLBEforthetem-peraturedistributionfunctionTDFandthecalculationofthetemperatureincludecomplicatedtermsinvolvingtempo-ralandspatialderivativesofthemacroscopicowvariables,whichmayintroducesomeadditionalerrorsanddoharmtothenumericalstability.
Furthermore,inthederivationoftheequilibriumfortheTDF,anadhotregroupingtechniquewasemployedsothatsomehigh-ordertermscouldbeneglected.
Theregroupingissomewhatarbitrary,anddifferentregroup-ingmethodsmayleadtodifferentequilibria9,13.
Itisnotedthatsomeimprovedversionshavebeenproposedbysomegroups.
Forinstance,asimpliedmodelwasderivedbydroppingthespatialgradientterminthetemperatureLBEforthermaluidswhereviscousheatdissipationandcom-pressionworkareneglected12,whichissimilartootherDDFLBEmodelsthatdonotconsidertheviscousandcom-pressioneffectsontheenergy.
Recently,Shietal.
proposedanotherversionbyregroupingtheTaylorexpansionsofthecontinuousequilibriumfortheTDF13.
Unfortunately,thecomplicatedspatialgradienttermsstillexistinthemodelifviscousdissipationisincluded.
Inthispaper,weaimtoproposeanalternativeTLBEmodel,inwhichviscousdissipationandcompressionworkareconsideredintheDDFframework.
Tothisend,werstintroduceadistributionfunctionthatrepresentsthetotalen-ergyratherthanthetemperatureorinternalenergyasthesecondDFinadditiontotheVDF.
Thenweconstructaki-neticequationforthetotalenergydistributionfunctionTEDFbasedonthekineticequationfortheVDF.
BasedonakineticmodelconstitutedoftwokineticequationsfortheVDFandTEDF,respectively,adiscretevelocitymodelDVMisproposedbychoosinganappropriatediscreteve-locitysetbasedonaHermiteexpansionoftheequilibriumfortheVDFandTEDF.
Furtherdiscretizationsofthetempo-ralandspatialderivativesoftheDVMleadtoourTLBE.
It*Correspondingauthor.
Electronicaddress:zlguo@hust.
edu.
cnPHYSICALREVIEWE75,03670420071539-3755/2007/753/036704152007TheAmericanPhysicalSociety036704-1shouldbenotedthat,unliketheHCDmodel,theuseoftheTEDFenablestheproposedLBEmodeltobesimplerwith-outthecomplicatedspatialgradientterms;andtheexpansionofthecontinuousequilibriumfortheTEDFintoaseriesoftensorHermitepolynomialsinsteadofTaylorseriesallowsthediscreteequilibriumtobedetermineduniquely.
Therestofthepaperisorganizedasfollows.
InSec.
II,theDFforthetotalenergyisintroducedanditskineticequa-tionisconstructedbasedontheBoltzmannequation.
InSec.
III,adiscretevelocitymodelisdevelopedfromthekineticequationsbyexpandingtheenergyDFintoaseriesoftensorHermitepolynomialsandapplyingtheGauss-Hermitquadra-turetothevelocitymomentsoftheDF.
SectionIVpresentsthederivationofthethermalLBEmodelfromtheDVMbyemployingsomestandardnumericalprocedures.
NumericaltestsoftheLBEmodelaremadeinSec.
VbysimulatingthermalPoiseuilleowandnaturalconvectioninasquarecavity,andnallyabriefsummaryisgiveninSec.
VI.
II.
KINETICMODELWITHDIFFERENTMOMENTUMANDENERGYRELAXATIONTIMESA.
KineticmodelsofBhatnagar-Gross-KrooktypeInkinetictheory,amonatomicgasisdescribedbythevelocitydistributionfunctionfx,,tofthemolecules,whichisdenedsuchthatfddxistheprobabilityofndingamoleculemovingwithvelocityatpositionxandtimet.
TheevolutionoftheVDFisgovernedbytheBoltzmannequation24tf+·f+a·f=f,1whereaistheacceleration,andfisacollisionoperatorthatsatisesthefollowingconservationlaws:fd=0,2where=1,,2.
TheoriginalBoltzmanncollisionop-eratorisacomplicatedintegralthatdependsontheinterpar-ticlepotentials.
Inpracticalapplications,fisusuallyap-proximatedbysomesimpliedmodels.
Onewidelyusedapproximationistheso-calledsingle-relaxation-timeorBhatnagar-Gross-KrookBGKmodel25f=1fffeq,3wherefistherelaxationtime,andfeqisthelocalMax-wellianequilibriumdistributionfunctionEDFdenedbyfeq;,u,T=2RTD/2expu22RT,4withDbeingthespatialdimensionandRthegasconstant.
TheuidvariablesintheEDF,i.
e.
,thedensity,velocityu,andtemperatureT,aredenedasthemomentsoff,uDRT2=fdfdu22fd.
5AlthoughtheBGKmodelretainsthemainfeaturesoftheoriginalBoltzmanncollisionoperator,itislimitedtogaseswithaxedPrandtlnumber24,whichhasbeenrecognizedasoneofthemaindefectsofthemodel.
SomeeffortshavebeenmadetoovercomethexedPrandtlnumberproblemoftheBGKmodel.
Forinstance,onecanusetheso-calledel-lipsoidalstatisticalBGKESBGKmodelwheretheMax-wellEDFisreplacedwithananisotropicGaussianEDF26,oruseavelocity-dependentrelaxationtime27.
Alterna-tively,Wooddirectlyintroducedtworelaxationtimesintothenonequilibriumdistributionfunctionafternoticingthatthemomentumandenergyshouldhavedifferenttransporttimescalesduringthecollisionprocess28.
Recently,Heetal.
proposedanotherapproachtoxthePrandtlnumberproblembyintroducinganewvariabletheso-calledinternalenergydistributionfunctionasaninter-nalenergyDF,g=u22f.
6ABGK-typekineticequationforgisthenconstructedbasedontheBoltzmannequation1,whichallowsfortheenergyhavingarelaxationtimescaledifferentfromthatofthemo-mentumtransport:tg+·g=1gggeqfq,7wheregistherelaxationtimefortheenergytransport,andgeqisthecorrespondingenergyEDFdenedbygeq=u222RTD/2expu22RT.
8ThequantityqinEq.
7isgivenbyq=u·tu+·u.
9Inthismodel,thedensityandvelocityarestilldeterminedbythemomentsofthedensityDFfasEq.
5,buttheinternalenergyDRT/2isnowdenedbytheenergyDFg:=gd.
10ThroughaChapman-Enskoganalysis,Heetal.
wereabletoshowthatthemacroscopicequationsderivedfromEqs.
1and7couldhaveaproperPrandtlnumbergiventhatgischosenappropriately9.
Inessence,theapproachusedbytheHCDmodelisiden-ticaltothatproposedbyWoodswheretworelaxationtimesareusedtodistinguishthemomentumandenergytransportduetoparticlecollisions.
ThemaindifferencebetweentheseGUOetal.
PHYSICALREVIEWE75,0367042007036704-2twoapproachesliesintherealizationofthetimescalesepa-ration:InthemethodofWoods,thetworelaxationtimesaredirectlyappliedtotherst-orderapproximationoftheVDFfintheChapman-Enskogexpansion,whiletheHCDmodelseparatestheenergytransportfromthemomentumtransportwithtwotimescalesexplicitly.
ItisdifculttoconstructaLBEmodelbasedonWoods'methoddirectlysincetheChapman-EnskogexpansioncannotbeemployedintheLBEexplicitly.
Onthecontrary,theHCDmodelcanserveasagoodbasefortheLBE9.
However,thetermfqinthekineticequationfortheenergyDFoftheHCDmodelmakesthenalLBEmodelcontainsometermsinvolvingspatialgradientsofboththedensityandvelocity.
ThesetermsnotonlymaketheLBEmodelmorecomputationallyexpensiveandmaydoharmtothenumericalstability,butalsomayleadtosomeunphysicalphenomenainuidsystemscontain-inglargespatialgradients,suchasmultiphaseormulticom-ponentandmicroscaleows.
B.
TotalenergydistributionfunctionanditskineticequationUnliketheHCDmodel,whichusestheinternalenergydistributionfunctiong,hereweintroducethefollowingtotalenergydistributionfunction:h=22f,11fromwhichthetotalenergyEcanbedenedasE=+u22=hd.
12TheevolutionofhcanbeobtainedfromtheBoltzmannequation1asth+·h+a·hf=h,13whereh=2f/2isthecollisionoperatorcharacterizingtheenergychangeduringtheparticlecollisions.
ThekeypointfordevelopingakineticmodelbasedonthetotalenergyDFhistospecifythecollisiontermhinEq.
13withsoundphysics.
Bynotingthatthecontributionofhincludestheinternalenergypartandthemechanicalen-ergypart,werstdecomposehintothesetwoparts:h=i+m,14wherei=u2f/2istheinternalenergypart,andm=hi=22u22fZfisthemechanicalenergypart,withZ=·uu2/2.
AccordingtoWoods'theory28,mshouldhavethesametimescaleasthatoff,andthereforeweapproximateitasm=Zfffeq.
15Fortheinternalenergypart,wecanapproximateitwithanotherBGK-typemodelasintheHCDmodel,i.
e.
,i=g1ggeq.
ButsuchamodelwillintroducetheinternalenergyDFgasanauxiliaryvariable.
Inordertoavoidthisinconvenience,wereplacegwithhZf,andthusobtainthefollowingBGK-typemodel:i=1hhheqZffeq,16whereh=gistherelaxationtimecharacterizingtheinternalenergychangeduringtheparticlecollisions,andheq2feq/2isthecorrespondingEDF.
Assuch,thenalcol-lisionoperatorhisgivenbyh=1hhheq+Zhfffeq,17where1hf=1h1f.
Itisclearthat,ash=f,thesecondtermofhvanishesandthemodelisidenticaltotheoriginalBGKmodel.
Oth-erwise,thesecondtermcanbeviewedasacorrectiontothesingle-relaxation-timemodel.
Aswillbeseenlater,withoutthisterm,themodelgivesincorrectviscousheatdissipationintheenergyequation,althoughthePrandtlnumbercanbetunedtobecorrect.
Insummary,weproposethefollowingtwo-relaxation-timemodelfordescribingathermaluidsystemwithavari-ablePrandtlnumber:tf+·f+a·f=1fffeq,18ath+·h+a·h=1hhheq+Zhfffeq+f·a,18bwherefeq=2RTD/2expu22RT,19aheq=222RTD/2expu22RT.
19bTheuidvariablesaredenedasuE=fdfdhd.
20ThroughtheChapman-Enskogexpansion,wecanobtainthefollowinghydrodynamicequationsattheNavier-StokeslevelseeAppendixAfordetails:t+·u=0,21aTHERMALLATTICEBOLTZMANNEQUATIONFORLOW…PHYSICALREVIEWE75,0367042007036704-3tu+·uu=p+·+a,21btE+·p+Eu=·T+··u+u·a,21cwherep=RTisthepressure,=S2/D·uIS=u+uistheviscousstresstensor,andtheviscosityandthermalconductivityaregivenby=fpand=D+2R2hp=cphp,respectively,wherecp=D+2R/2isthespecicheatcoef-cientatconstantpressure.
Usingthemomentumequation21b,wecandeducethetemperatureequationfromthetotalenergyequation21cascvtT+·uT=·Tp·u+·u,22wherecv=DR/2isthespecicheatcoefcientatconstantvolume.
ThePrandtlnumberofthesystem,Pr=cp/=f/h,canbemadearbitrarybytuningthetworelaxationtimes.
ThisresultisjustthesameasthatoftheWoods'model28.
Wewouldliketopointoutthattheabovekineticmodelcanalsobeextendedtopolyatomicgases.
Insuchacase,theVDFfisalsoafunctionoftherotationaland/orvibrationalenergiesthatcanbeeitherdiscreteorcontinuous27.
Inthecontinuouscase,fcanbeexpressedasf=fx,,,t,whereisavectorcontainingKcomponentscorrespondingtotheinternalfreedoms.
Accordingly,aBGK-typemodelcanbeusedtoapproximatethecollisionoperator29:tf+·f+a·f=1fffeq,23wherefeq=2RTD+K/2expu2+22RT,24andtheuidvariablesarenowdenedasuD+KRT2=fddfddu2+22fdd.
25ByintroducingtworeducedDFs,f=fdandh=2+2/2fd,wecanobtainthefollowingtwokineticequa-tionsfromtheBoltzmannequation1:tf+·f+a·f=1fffeq,26ath+·h+a·h=1hhheq+Zhfffeq+f·a,26bwherefeq=feqd=2RTD/2expu22RT,27aheq=2+22feqd=2+KRT22RTD/2expu22RT.
27bTheuidvariablesaredenednowasuE=fdfdhd.
28ThemacroscopicequationsderivedfromthemodelarethesameasthosegiveninEq.
21,butintheenergyequationthespecicheatscvandcpincludenowthecontributionoftherotationaland/orvibrationalenergies,i.
e.
,cv=D+KR/2andcp=D+K+2R/2.
III.
THEDISCRETEVELOCITYMODELPreviousstudieshaveshownthatwecanderiveaLBEmodelfromagivenkineticmodelfollowingsomestandardprocedures30–32.
Insuchanapproach,adiscretevelocitymodelisrstconstructedbydiscretizingthevelocityspaceofthecontinuouskineticequationintoanitesetofdiscretevelocities,andthentheLBEmodelisobtainedbydiscretiz-ingthetemporalandspatialderivationsoftheDVM,usingsomestandardnumericalschemes.
Thekeypointfordevel-opingLBEmodelsfollowingthisapproachliesintherststep,i.
e.
,thederivationofthediscretevelocitysetsothattheDVMcanmatchtheoriginalkineticmodelwithsufcientaccuracy.
Inthissection,wewillconcentrateonthisstepandpresentaDVMforthermalowsbasedonthekineticmodelproposedinSec.
II.
A.
HermiteexpansionsoftheequilibriumdistributionfunctionsInordertoobtainthecorrecthydrodynamicequations,thevelocityspaceofthekineticmodelmustbediscretizedwithsufcientaccuracy,or,inotherwords,thephysicalsymmetryoftheresultingdiscretevelocitysetshouldbeadequate.
Tothisend,ithasbeensuggestedtoexpandtheEDFfeqaroundthestateatrestundertheconditionoflowMachnumber,eitherbyperformingaTaylorseriesexpansionuptou29,30,31,orbyprojectingfeqontothetensorHermitepolynomialbasisintermsoftheparticlevelocityanduptothesecondorder32.
Forisothermalows,bothexpansionsGUOetal.
PHYSICALREVIEWE75,0367042007036704-4givethesameresults.
Forthermalows,however,thetwomethodswillresultindifferentformulations.
Inthepresentwork,weprefertousetheprojectionmethodbecausetheexpansioncoefcientsobtainedinthiswayarejusttheve-locitymomentsofthedistributionfunction,andthetrunca-tionofhigher-ordertermsdonotdirectlyalterthelower-ordermomentsofthedistributionfunction.
TheHermiteexpansionsoffeqandheqgivenbyEqs.
27aand27bcanbeexpressedasfeq=,TnAnx,tn!
Hn,29aheq=,TnBnx,tn!
Hn,29bwhere,T=12RTD/2exp22RT,and=/RT;Hnarethenth-ordertensorHermitepolyno-mials.
TheexpansioncoefcientsinEq.
29,AnandBn,aregivenbyAn=feqHnd,Bn=heqHnfxid.
30AsseeninAppendixA,thederivationofthehydrody-namicequationsattheNavier-StokeslevelfromthekineticmodelproposedinSec.
IIrequiresthezeroth-throughthird-ordermomentsoffeqandzeroth-throughsecond-ordermo-mentsofheq.
Therefore,inordertoobtainthesameequa-tionsattheNavier-Stokeslevel,itisnecessarytokeepthetermsuptothirdorderintheHermiteexpansionoffeq,andtosecondorderintheexpansionofheq.
Withthesecoef-cients,thetruncatedHermiteexpansionsoffeqandheqcanbewrittenasfeq,3T=,T1+·uRT+12·uRT2u22RT+·u6RT·uRT23u2RT,31heq,2T=,TE+p+E·uRT+p22RTD+p+E2·uRT2u2RT,32ForlowMachows,thethird-orderterminfeq,3canbeneglected,andwecanusethetruncatedexpansionsoffeqandhequptothesecondorder,i.
e.
,feq,2T=,T1+·uRT+12·uRT2u22RT,33aheq,2T=,Tp·uRT+·uRT2u2RT+122RTD+Efeq,2.
33bAccordingly,thetermsassociatedwiththeexternalforceainthekineticmodel18,a·fanda·h,shouldalsobeprojectedontotheHermitebasis.
TheChapman-EnskoganalysisofthekineticmodelseeEqs.
A10andA11inAppendixAindicatesthat,inordertoobtaintheexactNavier-Stokesequations,itisadequatetotruncatetheHer-miteexpansionsofthetwotermsuptothesecondorderandrstorder,respectively.
Withthisinmind,aftersomestan-dardmanipulationsweobtaina·f=,T·aRT+·a·uRT2a·uRT,34aa·h=,TE·aRT.
34bItcanbereadilyveriedthatthethermohydrodynamicequationscorrespondingtothetruncatedEDFs33andtheforcingterms34arejustthesameasthosefortheoriginalunexpandedonesafterneglectingthetermsofOMa3MarepresentstheMachnumber.
Here,wewouldliketopointoutthatifweapplytheHermiteexpansiontotheinternalenergyEDFgeqintheHCDmodel,wecanobtainthefollowingtruncatedexpan-sionuptothesecondorder:geq,2uRT+12·uRT2u22RT+2DRT=feq,22DRT1,35whichissimilartothosegiveninRefs.
9,13whichareobtainedbyregroupingtheTaylorexpansionofgeqheuris-tically.
B.
DiscretizationofthevelocityspaceThediscretevelocitysetcanbeobtainedbychoosingtheabscissaeofasuitableGauss-Hermitequadraturewiththeweightfunction,Tsothattherequiredvelocitymo-mentsofthetruncatedfeqcanbeexactlyevaluated.
How-ever,itisnotedthatthetemperatureappearinginthetrun-catedEDFsisalocallychangedvariable,whichmeansthattheabscissaeoftheGauss-Hermitequadraturearenotxed.
Thediscretevelocitiesobtainedinthisway,sayix,t,willalsodependonthelocaltemperatureandmaychangefrompositiontoposition.
Assuch,theresultantDVMcannotbeconsistentwiththeoriginalkineticmodelwherethecontinu-ousparticlevelocityisindependentoftimeandspace.
Asaresult,wecannotderivethecorrectthermohydrodynamicequationsbecauseoftheincommutabilityofiandthetem-poralandspatialgradients.
Inordertoovercomethisdifculty,wereplacethelocaltemperatureTinthetruncatedEDFswithareferencetem-THERMALLATTICEBOLTZMANNEQUATIONFORLOW…PHYSICALREVIEWE75,0367042007036704-5peratureT0,justlikethestrategyadoptedintheHCDmodel9.
Withthisreplacement,theEDFsfortheVDFandTEDFnowbecomefeq,2T0=,T01+·uRT0+12·uRT02u22RT0,36aheq,2T0=,T0p0·uRT0+·uRT02u2RT0+122RT0D+Efeq,2T0,36bwhere,T0=2RT0D/2exp2/2RT0andp0=RT0.
Accordingly,thelocaltemperatureTappearingintheforcingtermsisalsoreplacedwithT0:a·f=,T0·aRT0+·a·uRT02a·uRT0,37aa·h=,T0E·aRT0.
37bButnoticethatthetemperatureappearinginthetotalenergyEisstillthelocalvalue:E=cvT+u2/2.
Itiseasytoverifythatthezeroth-andrst-ordermomentsoffeq,2T0andthezeroth-ordermomentofheq,2T0arethesameasthoseoftheEDFswiththelocaltemperatureTgivenbyEq.
33,i.
e.
,feq,2T0d=,feq,2T0d=u,38aheq,2T0d=E.
38bThehighermomentsrequiredinthederivationoftheNavier-Stokesequations,however,aredifferentbecauseofthere-placementofTwithT0:feq,2T0d=uu+p0,39afeq,2T0d=p0u+u+u,39bheq,2T0d=p0+Eu,39cheq,2T0d=p0RT0+E+2p0+Euu.
39dBasedonthetwomodiedEDFsgivenbyEq.
36,wecannowdeterminethediscretevelocityseteasilyfromcer-tainGauss-Hermitequadratureswiththeweightfunction,T0.
Thequadratureshouldbeaccurateenoughsothatthevelocitymoments38and39canbeevaluatedexactly.
Forfeq,2T0,sincethethird-ordervelocitymomentneedstobeevaluatedaccuratelyinordertoobtainthehydrodynamicequationsattheNavier-Stokesorder,aGauss-Hermitequadraturewithatleastthefthdegreeofprecisionisre-quirednoticethatfeq,2T0itselfcontainssecond-ordertermsof.
Forheq,2T0,ontheotherhand,aGauss-Hermitequadraturewithatleastthefourthdegreeofpreci-sionisrequiredbecausethesecond-ordermomentneedstobeevaluatedexactlyinthederivationoftheenergyequation.
Therefore,aGauss-HermitequadraturewiththefthdegreeofprecisioncanbechosentodeterminethediscretevelocitysetforthekineticmodelwiththemodiedEDFs36.
Forthetwo-dimensionalcase,wecanchoosethenine-pointfth-degreeGauss-Hermitequadrature,whichleadstothefollowingdiscretevelocitiesD2Q9model:ci=0,0i=0,ccosi12,sini12,i=1,2,3,4,ccosi922,sini922,i=5,6,7,8,wherec=3RT0.
Theweightcoefcientscorrespondingtothesevelocitiesarew0=4/9,w1=w2=w3=w4=1/9,andw5=w6=w7=w8=1/36.
Similarly,forthree-dimensionalcasewecanobtainthe15-velocityD3Q15and19-velocityD3Q19models33.
OncetheabscissaoftheGauss-Hermitequadratureischosen,theintegralofafunctionof,say,canbeevalu-atedas,T0d=i=1bwici,40wherecii=1,2,…,bistheabscissaandwiisthequadra-tureweight.
Therefore,ifwedenefix,t=wifx,ci,tci,T0,hix,t=wihx,ci,tci,T0,wecanevaluatetheintegralsinEq.
20orEq.
28usingthequadratureanddeterminetheuidvariablesas=ifi,u=icifi,E=ihi.
41Theevolutionequationsforthereduceddistributionfunc-tionsfiandhicanbeeasilyderivedfromthekineticequa-tions18forfandh,whichleadtothefollowingdiscretevelocitymodel;GUOetal.
PHYSICALREVIEWE75,0367042007036704-6tfi+ci·fi=1ffifieq+Fi,42athi+ci·hi=1hhihieq+Zihffifieq+qi,42bwhereZi=ci·uu2/2,andfieq=wi1+ci·uRT0+12ci·uRT02u22RT0,43ahieq=wip0ci·uRT0+ci·uRT02u2RT0+12ci2RT0D+Efieq,43bFiandqiaretwotermsrelatedtotheexternalforce:Fi=wici·aRT0+ci·aci·uRT02a·uRT0,44aqi=wiEci·aRT0+fici·a.
44bThroughtheChapman-Enskoganalysis,thethermohydro-dynamicequationscorrespondingtotheDVM42canbederivedattheNavier-StokeslevelasseeAppendixB:t+·u=0,45atu+·uu=p0+·+a,45btE+·p0+Eu=·T+··u+u·a,45cwherep0=RT0,=S,=fp0,and=cvhp0.
AlthoughEqs.
45looksimilartothosederivedfromtheoriginalkineticmodel,i.
e.
,Eqs.
21,thefollowingdiffer-encesbetweenthemshouldbenoticed.
First,theequationofstateandthetransportcoefcientsinEqs.
45dependonlyonthereferencetemperatureT0,whilethoseinEqs.
21dependonthelocaltemperature.
Inotherwords,thethermo-hydrodynamicequations21arefullycoupled,whileinEqs.
45theenergyequationisdecoupledfromthemomentumequationsinceitcanbesolvedindependentlyoncethersttwoequationsaresolved.
Inthissense,theDVM42canbetermedadecoupledDVM.
Theseconddifferenceliesintheviscousstresstensor.
InEqs.
45thetraceofisnonzero,whichmeansthatthebulkviscosityisalsononzeroandequalto2/D.
Onthecontrary,in21theviscousstressistracelessandthebulkviscosityiszero.
Thenaldifferenceappearsinthethermalconductivity:inEqs.
21itisgivenby=cphp,butinEqs.
45itisgivenby=cvhp.
There-fore,fortheDVMthetworelaxationtimesarerelatedtothePrandtlnumberasPr=cp/=f/h,whichisdifferentfromtheresultofthecontinuouskineticmodel.
IV.
THERMALLATTICEBOLTZMANNMODELA.
LatticeBoltzmannequationsBaseontheDVMpresentedintheabovesection,wecanconstructathermalLBEmodelforlowMachowsbydis-cretizingthetemporalandspatialderivativesfollowingsomestandardprocedures.
First,thetimediscretizationforEq.
42acanbemadebyintegratingtheequationalongthechar-acteristicline,whichleadstofix+cit,t+tfix,t=0tfx+cit,t+t+Fix+cit,t+tdt,46wheretisthetimestepandf=fieqfi/f.
AsarguedinRef.
9,theintegralontheright-handsidemustbeevalu-atedwithatleastsecond-orderaccuracy.
Thetrapezoidalrulecanservethispurposeandleadstothefollowingtime-discretescheme:fix+cit,t+tfix,t=t2fx+cit,t+t+Fix+cit,t+t+t2fx,t+Fix,t.
47AssuggestedbyHeetal.
9,theimplicitnessoftheaboveschemecanbeeliminatedbyintroducingthefollowingdis-tributionfunction:fi=fit2f+Fi,48fromwhichonecanobtainfifieq=1+t2f1fifieq+t2Fi49andfifi=t2f+Fi.
50WiththeaidsofEqs.
49and50,Eq.
47canberewrittenasfix+cit,t+tfix,t=ffix,tfieqx,t+t1f2Fi,51wheref=2t/2f+t.
FromEq.
48,itcanbeeasilyveri-edthatthedensityandvelocityoftheuidcanbecom-putedfromthenewVDFas=ifi,u=icifi+t2a.
52ItisnotedthatthetreatmentoftheforcingterminEq.
51isjustthesameasthatproposedin34.
THERMALLATTICEBOLTZMANNEQUATIONFORLOW…PHYSICALREVIEWE75,0367042007036704-7Thediscretekineticequation42bcanbediscretizedus-ingasimilarapproach.
Specically,theuseofcharacteristicdiscretizationandtrapezoidalquadratureleadstothefollow-ingimplicitscheme:hix+cit,t+thix,t=0thx+cit,t+t+Six+cit,t+tdt,53whereSi=Zif/hf+qi.
TheimplicitnesscanbeeliminatedbyintroducingthefollowingenergyDF:hi=hit2h+Si,54fromwhichonecanobtainhihieq=1+t2h1hihieq+t2Si55andhihi=t2h+Si.
56Therefore,theLBEforthisenergyDFisnowhix+cit,t+thix,t=th+Si=hhix,thieqx,t+t1h2Si=hhix,thieqx,t+t1h2qi+t1h2Zihffifieq,57wheref=2t/2f+t.
TheunderlinedtermintheaboveequationcanbefurtherrewrittenintermsoftheVDFfi.
Infact,bynoticingthatf=tf10.
5andh=th10.
5,wecanobtainthat1hf1+t2f1=hft1h/2.
58Therefore,aftersubstitutingtheexpression49,theunder-linedtermbecomeshfZififieq+t2Fi.
59Asaresult,thenalformulationofthetime-discreteschemefortheenergyequationcanbewrittenashix+cit,t+thix,t=hhix,thieqx,t+t1h2qi+hfZififieq+t2Fi,60andthetotalenergycannowbedeterminedbyE=ihi+t2u·a.
61Afterthetimeisdiscretized,wenowarriveatanotherkeypointforconstructingtheLBE,i.
e.
,thediscretizationofthespace.
InastandardLBEmodel,thespatialspaceisdis-cretizedintoaregularlatticeLwithaspacingsuchthat,foranodexL,thepointx+citisthenearestlatticepointofxalongthedirectionci.
Thismeansthat,ifaparticlemovingwithadiscretevelocityislocatedonthelattice,itwilljumptothenearestneighboronthelatticeinthenextstep.
Forinstance,intheD2Q9modelthelatticespacingischosentobex=ctwithc=3RT0.
Therefore,thelatticeiscloselydependentonthediscretevelocityset,meaningthatthespa-tialdiscretizationiscoupledwiththevelocitydiscretization.
Withthetimeintegrationandthespacediscretization,thetwodiscretekineticequations51and60becomefullydiscretenow,andtheyconstituteourdecouplingthermalLBEmodelforlowMachnumberows.
ThepresentLBEmodelcanalsobeeasilyextendedtopolyatomicgasesbysimplychangingthespecicheatcvfromDR/2toD+KR/2.
ItcanbeeasilyshownthatthethermohydrodynamicequationsderivedfromtheLBEmodelarejustthoseoftheDVM,i.
e.
,45.
Therefore,theoreticallythepresentLBEmodelisapplicableonlytoBoussinesqowswherethesoundspeedandthetransportcoefcientsareindependentoftemperature.
Itisalsobecauseofthisfactthatwecallthemodeladecouplingone.
However,fromthecomputationalpointofview,therestrictiononthetransportcoefcientscanbereleasedtosomeextent.
Thatis,ifweknowtherelationbetweenthesecoefcientsandthelocaltemperatureinad-vance,wecanmakefandhfunctionsoftemperature,f=Tp0andh=Tcvp0.
Withsuchmodications,theLBEmodelcanalsobeappliedtoowswithtemperature-dependenttransportcoefcients.
ItisalsointerestingtomakeacomparisonbetweenthepresentLBEmodelandtheHCDmodel9.
First,itisnotedthatbothmodelssharemanysimilarfeatures:thelowMachnumberlimit,thereplacementoflocaltemperaturewithareferenceoneintheEDFs,andthedecouplingbetweenthemomentumandenergyequations.
Ontheotherhand,thedifferencesbetweenthetwomodelsarealsoapparent:theHCDmodelusestheinternalenergydistribution,whilethepresentmodelemploysthetotalenergydistribution.
Asaresultofthisdifferentchoice,theHCDmodelcontainsacomplicateddifferentialtermthatneedsspecialtreatment,whilethepresentmodelisabletoavoidsuchdifculty.
B.
Somespecialcases1.
FlowswithnegligiblecompressionworkandviscousdissipationInmanypracticalapplications,compressionworkandvis-cousheatcanbeneglected.
ThepresentLBEmodelcanbeGUOetal.
PHYSICALREVIEWE75,0367042007036704-8easilyappliedtosuchproblemsbysimplytakingcp=cv→inthemodel.
Thiscanbeseenmoreclearlyifwere-writetheenergyequation45cinnondimensionalformast+·u=·PrReEcp0·u+EcReS:u,62where=TT0/T,withT0andTbeingthecharacteris-tictemperatureandtemperaturevariation.
TheparametersPr=cp/,Re=Lu0/,andEc=u02/cpTarethePrandtlnumber,Reynoldsnumber,andEckertnumber,respectively,withLbeingthecharacteristiclengthandu0thecharacteris-ticvelocity.
=cp/cvistheratioofthespecicheats,whichcanbetakentobeunityforincompressibleuids.
Ascpislargeenough,theEckertnumberwillbecomesufcientsmallsothatthecompressionworkandtheviscousheatdissipationtermsarenegligible.
Itisnotedthat,fortheHCDmodel,theevolutionequa-tionsshouldalsobemodiedifcompressionworkandvis-cousheatareneglected9.
However,despitetheneglectofthetwofactors,acomplicatedgradienttermsimilartothatintheoriginalHCDmodelstillexistsinthemodiedHCDmodel.
2.
FlowswithbuoyancyforceInnaturalconvectionandmixed-convectionows,thebuoyancyforceshouldbeconsidered.
AsthetemperaturedifferenceTissmallincomparisonwiththeaveragetem-peratureT0,theBoussinesqassumptioncanbeinvoked.
Thatis,theuidpropertiesareassumedtobeindependentofthetemperature,exceptthatthedensityinthegravitationalforceisassumedtobe=01TT0,63where0istheuiddensityattemperatureT0andisthethermalexpansioncoefcient.
Assuch,thegravityforcecanbeexpressedasg=0g0gTT0,64wheregistheaccelerationduetogravity.
Afterabsorbingtheconstantpart0gintothepressure,theeffectiveexternalforcebecomesa=0gTT0.
However,itshouldbenotedthat,withthiseffectiveforce,thepressureeldmodeledbytheLBEisactuallythedy-namicpart,p=p00gy,withgbeingthemagnitudeofthegravityandythecoordinateoppositetothegravityforce.
Withsuchatreatment,theworkdonebythepressurethatenterstheenergyequationcontainsonlythedynamicpart.
Inotherwords,themomentumandenergyequationscorre-spondingtotheLBEareactuallytu+·uu=p+·+a,65atE+·p+Eu=·T+··u+u·a,65bwithp=RT0thedynamicpressureand0thecorre-spondinguiddensity.
Asthecompressionworkisnegli-gible,suchatreatmentcanworkwell.
Ontheotherhand,whenthecompressionworkplaysanimportantrole,thetotalpressureshouldbeusedintheenergyequationforthesakeofthermodynamicconsistency38.
Inordertoaccountforthiseffect,weincludetheworkdonebythestaticpressureintothetermqiintheLBE60,qi=wiERT0+fici·a+iu·g.
66Accordingly,thetotalenergyiscalculatedbyE=ihi+t2u·a+g.
67OnecanshowthatwiththemodiedqitheenergyequationcorrespondingtotheLBEmodeliscvtT+u·T=·Tp·u+·u+u·g,68whichissimilartothatinRef.
38intheincompressiblelimiti.
e.
,0.
C.
BoundaryconditionsInpracticalapplications,theowboundaryconditionsareusuallyspeciedintermsoftheuidvariables.
Inordertotransformthermohydrodynamicboundaryconditionstotheboundaryconditionsforthedistributionfunctions,weem-ploythenonequilibrium-extrapolationapproachinthisworkduetoitssimplicity,second-orderaccuracy,andgoodrobust-ness35.
TheapproachwasoriginallyproposedtorealizeplaneboundariesforisothermalLBEs.
Recently,thisap-proachhasbeenextendedtocurveboundaries36andtoTLBEs11,37.
Thebasicideaofthenonequilibriumextrapolationap-proachistoseparateaDFataboundarynodeintoitsequi-libriumandnonequilibriumparts,wherethehydrodynamicboundaryconditionsareenforcedthroughtheequilibrium,andthenonequilibriumpartisapproximatedbythatoftheDFatthenearestneighbornodeintheuidregion.
Follow-ingthisapproach,thedensityDFfiataboundarynodexbcanbespeciedasfixb=fieqxb,b,ub+fixffieqxf,69wherexfisthenearestuidneighborhood.
Forthevelocityboundaryconditionwherethevelocityubisknown,bisjustaparameter,notnecessarilyequaltothedensityatxb.
Ithasbeendemonstratedthatitisagoodapproximationtosetb=xf11,35,36.
Similarly,forthermalboundarycondi-tionswherethetemperatureattheboundaryisknown,theenergyDFhiisapproximatedashixb=hieqxb,b,Eb+hixfhieqxf,70whereEb=cvTb+ub2/2.
Itisnotedthat,iftheheatuxq˙=T/nisspeciedattheboundary,withnbeingtheunitTHERMALLATTICEBOLTZMANNEQUATIONFORLOW…PHYSICALREVIEWE75,0367042007036704-9vectornormaltotheboundary,theenergyDFcanalsobeapproximatedaccordingtoEq.
70,exceptthatTbisnowgivenbysomenumericalschemesoftheheatux.
V.
NUMERICALTESTSAnumberofsimulations,includingtheplanarthermalPoiseuilleowandthenaturalconvectioninasquarecavity,havebeencarriedouttovalidatethepresentthermalLBEmodel.
Inthesimulations,thetwo-dimensionalnine-speedD2Q9modelisemployed.
A.
PlanarthermalPoiseuilleowThethermalPoiseuilleowinaplanarchannelconsid-eredhereisdrivenbyaconstantforcea,andthetemperatureofthebottomandtopwallsofthechannelarekeptatThandTc,respectively.
Ifthegravityisneglected,thevelocityandthetemperatureprolescanbedescribedasuy=4u0y*1y*,71a=y*+PrEc3112y*4,71bwherey*=y/hhbeingthechannelheight,u0=0ah2/8,and=TTc/ThTc.
ThethermalPoiseuilleowischaracterizedbytheRey-noldsnumberRe=0hu0/,thePrandtlnumberPr=cp/,andtheEckertnumberEc=u02/cpThTc.
WecarriedoutasetofsimulationsfordifferentvaluesofRe,Pr,andEc.
Thespecicheatratioissettobeunitysincetheowcanbeconsideredtobeincompressible.
Inoursimulations,a6464latticeisemployed,andthenonequilibriumextrapola-tionmethodisusedtotreatvelocityandtemperaturebound-aryconditionsforthebottomandtopwalls69and70.
Inthestreamwisexdirection,periodicboundaryconditionsareappliedtotheinletandoutlet.
InoursimulationstheReynoldsnumberistakentobeRe=20andthemaximumvelocityu0issettobe1.
0.
Therelaxationparameterfissettobe0.
8sothatthecomputationalMachnumberu0/3RT0isabout0.
08,whichensuresthelowMachnumberrequire-ment.
Otherparameterscanbedeterminedfromthenondi-mensionalparameters.
ThetemperatureprolesforPr=0.
71asEcvariesfrom0.
1to100arepresentedinFig.
1,whilethetemperatureprolesforaxedEckertnumberEc=10asPrvariesfrom0.
1to4areshowninFig.
2.
Itisclearlyseenthatthenu-mericalresultsareinexcellentagreementwiththeanalyticalsolutions.
TheviscousheateffectsaresuccessfullycapturedbythepresentthermalLBEmodeloverawiderangeoftheproductofPrandEc.
ItisknownthatathermalLBEusuallylosesstabilityoraccuracyasfand/orhapproachto2.
0orwhenthetem-peraturevariationislarge.
Inordertodemonstratetherangeofapplicabilityofthepresentmodel,wetestthemodelbyvaryingthetworelaxationstimesfandh,andthetempera-turedifferenceThTc/Tc.
Inthesimulationsthedrivenforceischosensuchthatthemaximumvelocityu0=0.
1c,withthe6464lattice.
ItisfoundthattheLBEisstillstableandaccurateinbothvelocityandtemperatureevenwhenbothf/tandh/tareassmallas104whenThTc/Tcrangesfrom0to1000.
ItisclearthatapplicabilityrangeandnumericalstabilityofthepresentLBEaresimilartothoseoftheHCDmodel9andaremuchwiderandbetterthanthoseofthemultispeedLBEmodel5.
B.
NaturalconvectioninasquarecavityWenowapplythethermalLBEmodeltothenaturalcon-vectionowinatwo-dimensionalsquarecavity.
ThetwosidewallsleftandrightofthecavityaremaintainedattwodifferenttemperaturesTcandThThTc,respectively,whilethebottomandtopwallsareadiabatic.
TheconvectionowinducedbythetemperaturedifferenceischaracterizedbythePrandtlnumberPrandRayleighnumberRa=02cpgTH3/,where0isthereferencedensity,T=ThTcisthetemperaturedifferencebetweenthehotandcoolwalls,andHistheheightofthecavity.
Werstsimulatedthenaturalconvectionproblemwithnegligiblecompressionworkandheatdissipation.
ThisisachievedbysettingtheEckertnumbertobeassmallas1030.
ThePrandtlnumberissettobe0.
71,andtheRayleighnumberrangesfrom103to106.
Inthecomputationsa128FIG.
1.
Temperaturevariation=TTc/ThTcofthether-malPoiseuilleowatRe=20andPr=0.
71.
a–dEc=0.
1,20,50,and100.
Solidlinesaretheanalyticalsolutions,andthesymbolsarethenumericalresults.
FIG.
2.
Temperaturevariation=TTc/ThTcofthether-malPoiseuilleowatRe=20andEc=10.
0.
a–dPr=0.
1,1.
0,2.
0,and4.
0.
Solidlinesaretheanalyticalsolutions,andthesym-bolsarethenumericalresults.
GUOetal.
PHYSICALREVIEWE75,0367042007036704-10128latticeisemployedandtherelaxationparameterwfischosentobe1.
6forallcases.
Thenonequilibriumextrapo-lationmethodisappliedtospecifytheboundaryconditionsforboththedensityDFfiandtheenergyDFhiatthefoursolidwalls.
StreamlinesandisothermlinespredictedbythepresentTLBEmodelareshowninFigs.
3and4.
ItisseenthatforlowRaacentralvortexappearsasatypicalfeatureoftheow.
ThevortextendstobecomeellipticasRaincreases,andbreaksupintotwovorticesatRa=105.
AsRareaches106,thetwovorticesmovetowardthewallsandathirdvortexappearsinthecenterofthecavity.
TheisothermlinesindicatethechangeofthedominantheattransfermechanismwithRayleighnumber.
ForsmallRa,theheatistransferredmainlybyconductionbetweenthehotandcoldwalls,andtheisothermsarealmostvertical.
AsRaincreases,thedomi-nantheattransfermechanismchangesfromconductiontoconvection,andtheisothermlinesbecomehorizontalinthecenterofthecavity,andareverticalonlyinthethinbound-arylayersnearthehotandcoldwalls.
Alloftheseobserva-tionsareingoodagreementwithresultsreportedinpreviousstudies39,40.
Toquantifytheresults,wecomputedtheNusseltnumberalongthetwosidewallsandthemaximumvelocitiesalongthehorizontalandverticallinesthroughthecavitycenter.
TheresultsarelistedinTableItogetherwiththedatafrompreviousstudies.
Asshown,theTLBEresultsagreewellwiththeavailabledata.
Infact,thedifferencebetweenthepresentLBEresultsandthereferenceonesarewithin1.
0%forthecasesconsidered.
Wenowexaminetheeffectsofcompressionworkandviscousdissipationonnaturalconvection.
ThePrandtlnum-berandtheRayleighnumberaresettobe1.
0and105,re-spectively.
Thesizeofthecomputationalmeshis256256andtherelaxationparameterwfissettobe1.
6.
Thetempera-turesofthecoolandhotsidewallsarekeptat300and310K,respectively.
TwovaluesoftheEckertnumberareusedinoursimulations,i.
e.
,Ec=1030and105,whereintheformercasethepressureworkandviscousdissipationareneglectedwhileinthelattercasetheeffectsareincluded.
InFig.
5,thestreamlinesandisothermallinesarepre-sentedforcomparison.
Itisclearlyseenthattheowandheattransferbehaviorsinthetwocasesarequitedifferent.
Ascompressionworkandviscousdissipationareconsidered,owoccursonlyintheregionsveryclosetothewalls,andtheisothermsareverydenseinthenear-wallregion,andrathersparseintheinteriorregionofthecavity.
TheNusseltnumbersonthecoolandhotwallsinthetwocasesarealsomeasured.
AslistedinTableII,theheattrans-ferisgreatlyenhancedifcompressionworkandviscousdis-sipationareconsidered,whichisconsistentwiththelargertemperaturegradientsinthenear-wallregionsshowninFig.
5.
TheseresultsalsoagreewellquantitativelywiththosereportedinRef.
38forthesamecase.
VI.
SUMMARYInthispaper,wehavedevelopedathermallatticeBoltz-mannequationforlow-speedowsbasedonatwo-TABLEI.
ComparisonsoftheaverageNusseltnumberandthemaximumvelocitycomponentsacrossthecavitycenter.
Thedatainparenthesesarethelocationsofthemaxima.
RaNuumaxyvmaxx103Present1.
11953.
6430.
80473.
69190.
1719Ref.
111.
11683.
65540.
81253.
69850.
1797104Present2.
254516.
12540.
820319.
55770.
1172Ref.
392.
244216.
18020.
826519.
62950.
1193105Present4.
527834.
60330.
851668.
08200.
0703Ref.
394.
521634.
73990.
855868.
63960.
0657106Present8.
774664.
90590.
8516218.
9000.
0391Ref.
398.
825164.
83670.
8505220.
4610.
0390FIG.
3.
StreamlinesforRa=103a,104b,105c,and106dofthenaturalconvectionowinacavity.
FIG.
4.
IsothermlinesforRa=103a,104b,105c,and106dofthenaturalconvectionowinacavity.
FIG.
5.
StreamlinesleftandisothermsrightforRa=105,Pr=1.
Top,without,andbottom,withconsiderationofthepressureworkandviscousdissipation.
THERMALLATTICEBOLTZMANNEQUATIONFORLOW…PHYSICALREVIEWE75,0367042007036704-11relaxation-timekineticmodel.
Theproposedmodeliscon-structedintheDDFframework.
Themostdistinctivefeatureofthepresentmodelisthatanadditionaldistributionfunc-tionisdenedtorepresentthetotalenergy,insteadofrepre-sentingeithertheinternalenergyorthetemperatureinpre-viousstudies.
ThischoicenotonlyenablesthenalTLBEmodeltobesimplebutalsomakestheinclusionofcompres-sionworkandviscousdissipationtobeeasier.
ThenumericaltestsshowthattheresultspredictedbytheTLBEmodelareexcellentagreementwiththeanalyticalsolutionsandnu-mericalresultsreportedinpreviousstudies.
Itshouldbepointedoutthatinthepresentmodeltheenergyequationisdecoupledfromthemomentumequationduetothereplacementofthelocaltemperaturewiththecon-stantreferencetemperatureintheequilibriumdistributionfunctions.
SuchadecouplingcausesthepresentTLBEmodeltobelimitedtoBoussinesqows,inwhichthetemperaturevariationissmall,suchthatthetransportcoefcientsandthesoundspeedbecomealmostindependentoftemperature.
TheextensionofthepresentTLBEmodeltosystemsinwhichthemomentumandenergytransportarecoupledisunderway.
ACKNOWLEDGMENTSThisworkissupportedbytheNationalBasicResearchProgramofChinaGrantNo.
2006CB705804andtheNa-tionalNaturalScienceFoundationofChinaGrantNo.
50606012.
APPENDIXA:CHAPMAN-ENSKOGANALYSISOFTHEKINETICMODELByintroducingthefollowingChapman-Enskogexpan-sions24:t=n=1ntn,=1,a=a1,A1af=n=0nfn,h=n=0nhn,A1bwhereisasmallexpansionparameter,wecanrewritethekineticequations18intheconsecutiveordersoftheparam-eteras0:f0=feq,A2a1:Dt1f0=f1f,A2b2:t2f0+Dt1f1=f2f,A2cand0:h0=heq,A3a1:Dt1h0=h1h+Zf1hf+f0·a1,A3b2:t2h0+Dt1h1=h2h+Zf2hf+f1·a1,A3cwhereDt1=t1+·1+a1·.
EquationsA2aandA3ain-dicatethatfnd=0,fnd=0,hnd=0A4forn0,becauseuE=fdfdhd=feqdfeqdheqd.
A5Furthermore,aftersomestandardalgebrawecanobtainthefollowingresults:0=f0d=p+uu,A6aQ0=h0d=p+Eu,A6bf0d=pu+u+uuuu,A6ch0d=pRT+E+2p+Euu,A6dwhereistheKroneckerdeltawithtwoindices.
FromEqs.
A2bandA3b,wecanobtainthethermohy-drodynamicequationsattherstorder:t1+1·u=0,A7at1u+1·uu+pI=a1,A7bTABLEII.
TheaverageNusseltnumbersathotNuhandcoolNucwallsforPr=1andRa=105.
EcNuhNuc1030Present4.
61284.
6128105Present13.
120113.
2447Ref.
3813.
20013.
198GUOetal.
PHYSICALREVIEWE75,0367042007036704-12t1E+1·p+Eu=u·a1,A7cwherep=RTisthepressure.
Similarly,themomentsofEqs.
A2bandA3bleadtothethermohydrodynamicequationsattheorderof2:t2=0,A8at2u+1·1=0,A8bt2E+1·Q1=0,A8cwhere1=f1dandQ1=h1d.
NotethatfromEq.
A7wecanobtainthatt1p+1·pu=2Dp1·u,A9at1uu+1·uuu=ups+aus,A9bt1pu+1·puu=RT1p+pa12D·upu,A9ct1Eu+1·p+Euu=uu·a1+Ea1Ep+pu·u,A9dwherethesymbol·sdenotesthesymmetricsummationofthebracketedsecond-ordertensor,suchasAs=A+ATthesuperscriptTdenotesthetransposeofthetensor.
WiththeaidsoftheseresultsandEqs.
A6candA6d,wecanobtainfromEqs.
A2bandA3bthat1f1=t10+1·f0d+a1·f0d=t1pI+t1uu+1pus+1·puI+1·uuua1us=pS12D1·uIA10and1hQ1+1hf1·u+0·a=t1Q0+1·h0d+a1·h0d=t1pu+t1Eu+1pRT+E+1·puu+1·p+EuuEa1=0·a1+pR1T2Du1·u+u·1u+1E=0·a1+D+22pR1T+pu·S12D1·uI=0·a1+D+22pR1T1f1·u,orQ1=D+22hpR1T+1·u,A11whereS1=1us.
Combiningtherst-andthesecond-orderresultsA7andA8,togetherwithEqs.
A10andA11,wearrivedatthethermohydrodynamicequationsattheNavier-Stokesor-der,t+·u=0,A12atu+·uu=p+·+a,A12btE+·p+Eu=·T+··u+u·a,A12cwhere=S2/D·uI,and=fpand=D+2R2hparetheviscosityandthermalconductivity,respectively.
APPENDIXB:THEHYDRODYNAMICEQUATIONSOFTHEDISCRETEVELOCITYMODELTheChapman-EnskoganalysisoftheDVM42issimi-lartothatpresentedinAppendixA.
Themaindifferencesliesinthetermsrelevanttothetemperature.
Specically,withthefollowingChapman-Enskogexpansions,t=n=1ntn,=1,a=a1,B1afi=n=0nfin,hi=n=0nhin,B1bwehaveFi=Fi1,qi=qi1+2qi2+,withFi1=wici·a1RT0+ci·a1ci·uRT02a1·uRT0,qi1=wiEci·a1/RT0+fi0ci·a1,qik=fik1ci·a1fork1.
Furthermore,itcanbeeasilyveriedthatifin=0,icifin=0,ihin=0B2forn0,andTHERMALLATTICEBOLTZMANNEQUATIONFORLOW…PHYSICALREVIEWE75,0367042007036704-130=icicifi0=p0+uu,B3aQ0=icihi0=p0+Eu,B3bicicicifi0=p0u+u+ua,B3cicicihi0=p0RT0+E+2p0+Euu,B3dwherefi0=fieqandhi0=hieq.
ItisalsonotedthatiFi=0,iciFi=a,iciciFi=aus,B4aiqi=u·a,iciqi=EI+·a.
B4bwhere=icicifi.
Therefore,therst-orderequationsintheexpansionofthediscretekineticequations42,Dit1f0=fi1f+Fi1,B5aDit1h0=hi1h+Zifi1hf+qi1,B5bwhereDit1=t1+ci·1canleadtothefollowingrst-orderthermohydrodynamicequations:t1+1·u=0,B6at1u+1·uu+p0I=a1,B6bt1E+1·p0+Eu=u·a1.
B6cWiththeseresults,wecanfurtherobtainthat1icicifi1andQ1icihi1havethefollowingexpres-sions:1f1=t10+1·icicicifi0a1us=t1p0I+t1uu+1p0us+1·p0uIa1us=p0S1B7and1hQ1=t1Q0+1icicihi01hf1uiciqi1=t1p0+Eu+1p0RT0+E+2p0+Euu1hf1up0+Ea+u·a1u=t1p0u+1p0RT0ta+p0uu+t1Eu+1p0E+p0+Euu+fhfp0S1up0+Ea+u·a1u=RT0t1u+1p0+1uu+ut1E+1p0+Eu+Et1u+u1u+p0u1u+1p0E+fhfp0S1utap0+Ea+u·a1u=p0a+u·a1u+Ea1p0+p0u1u+p01E+E1p0+fhfp0S1up0+Ea+u·a1u=p0u1u+1E+fuhfS1u=p0cv1T+p0S1u+fhfp0S1u=p0cv1T+fhp0S1u,B8wherewehaveneglectedthetermsoforderMa3intheabovedeductions.
Therefore,fromthesecond-orderequa-tionsoftheDVM,t2fi0+Dit1fi1=fi2f,B9at2hi0+Dit1hi1=hi2h+Zifi2hf+qi2,B9bwecaneasilyobtainthethermohydrodynamicequationsatthesecondorderof2:GUOetal.
PHYSICALREVIEWE75,0367042007036704-14t2=0,B10at2u+1·1=0,B10bt2E+1·Q1=0.
B10cFinally,basedontheresultsattheordersofand2,weobtainthefollowingthermohydrodynamicequationsattheNavier-Stokeslevel:t+·u=0,B11atu+·uu=p0+·+a,B11btE+·p0+Eu=·T+··u+u·a,B11cwhere=Swith=fp0,and=cvhp0withcv=DR/2.
1P.
LallemandandL.
-S.
Luo,Phys.
Rev.
E68,0367062003.
2S.
Succi,TheLatticeBoltzmannEquationforFluidDynamicsandBeyondOxfordUniversityPress,NewYork,2001.
3F.
J.
Alexander,S.
Chen,andJ.
D.
Sterling,Phys.
Rev.
E47,R22491993.
4Y.
Qian,J.
Sci.
Comput.
8,2311993.
5Y.
Chen,H.
Ohashi,andM.
Akiyama,Phys.
Rev.
E50,27761994.
6H.
Chen,C.
Teixeira,andK.
Molving,Int.
J.
Mod.
Phys.
C8,6751997;M.
Soe,G.
Vahala,P.
Pavlo,L.
Vahala,andH.
Chen,Phys.
Rev.
E57,42271998.
7G.
J.
M.
EggleandJ.
A.
Somers,Int.
J.
HeatFluidFlow16,3571995.
8X.
Shan,Phys.
Rev.
E55,27801997.
9X.
He,S.
Chen,andG.
D.
Doolen,J.
Comput.
Phys.
146,2821998.
10B.
J.
PalmerandD.
R.
Rector,J.
Comput.
Phys.
161,12000.
11Z.
L.
Guo,B.
Shi,andC.
Zheng,Int.
J.
Numer.
MethodsFluids39,3252002.
12Y.
Peng,C.
Shu,andY.
T.
Chew,Phys.
Rev.
E68,0267012003.
13Y.
Shi,T.
S.
Zhao,andZ.
L.
Guo,Phys.
Rev.
E70,0663102004.
14G.
Barrios,R.
Rechtman,J.
Rojas,andR.
Tovar,J.
FluidMech.
522,92005.
15T.
Watanabe,Phys.
Lett.
A321,1852004.
16T.
Watanabe,Phys.
Fluids16,9722004.
17A.
D'OrazioandS.
Succi,FGCS,FutureGener.
Comput.
Syst.
20,9352004.
18A.
D'Orazio,M.
Corcione,andG.
P.
Celata,Int.
J.
Therm.
Sci.
43,5752004.
19C.
Shu,X.
Niu,andY.
T.
Chew,J.
Stat.
Phys.
121,2392005.
20C.
-H.
WangandR.
Yang,Appl.
Math.
Comput.
173,12462006.
21Y.
Shi,T.
S.
Zhao,andZ.
L.
Guo,Comput.
Fluids35,12006.
22H.
N.
DixitandV.
Babu,Int.
J.
HeatMassTransfer49,7272006.
23D.
ChatterjeeandS.
Chakraborty,Phys.
Lett.
A351,3592006.
24S.
ChapmanandT.
G.
Cowling,TheMathematicalTheoryofNon-UniformGases,3rded.
CambridgeUniversityPress,Cambridge,U.
K.
,1970.
25P.
L.
Bhatnagar,E.
P.
Gross,andM.
Krook,Phys.
Rev.
94,5111954.
26L.
H.
Holway,Phys.
Fluids9,16581966.
27C.
Cercignani,TheBoltzmannEquationandItsApplicationsSpringer,Berlin,1988.
28L.
C.
Woods,AnIntroductiontotheKineticTheoryofGasesandMagnetoplasmasOxfordUniversityPress,Oxford,1993.
29K.
H.
PrendergastandK.
Xu,J.
Comput.
Phys.
107,531993.
30T.
Abe,J.
Comput.
Phys.
131,2411997.
31X.
HeandL.
-S.
Luo,Phys.
Rev.
E55,R63331997;56,68111997.
32X.
ShanandX.
He,Phys.
Rev.
Lett.
80,651998;X.
Shan,X.
-F.
Yuan,andH.
Chen,J.
FluidMech.
550,4132006.
33Y.
H.
Qian,D.
d'Humieres,andP.
Lallemand,Europhys.
Lett.
17,4791992.
34Z.
Guo,C.
Zheng,andB.
Shi,Phys.
Rev.
E65,0463082002.
35Z.
L.
Guo,C.
Zheng,andB.
Shi,Chin.
Phys.
11,03662002.
36Z.
L.
Guo,C.
Zheng,andB.
Shi,Phys.
Fluids14,20072002.
37G.
H.
Tang,W.
Q.
Tao,andY.
L.
He,Phys.
Rev.
E72,0167032005.
38V.
A.
F.
Costa,Int.
J.
Numer.
MethodsFluids48,23332005.
39M.
HortmannMandM.
Peri,Int.
J.
Numer.
MethodsFluids11,1891990.
40G.
Barakos,E.
Mitsoulis,andD.
Assimacopoulos,Int.
J.
Nu-mer.
MethodsFluids18,6951994.
THERMALLATTICEBOLTZMANNEQUATIONFORLOW…PHYSICALREVIEWE75,0367042007036704-15