The method pro posed here exploits physics unique to smoke in order to design a numerical method that is both fast and ef64257cient on the relatively coarse grids traditionally used in computer graphics applications as compared to the much 64257ner ID: 22882
Download Pdf The PPT/PDF document "Visual Simulation of Smoke Ronald Fedkiw..." is the property of its rightful owner. Permission is granted to download and print the materials on this web site for personal, non-commercial use only, and to display it on your personal computer provided you do not modify the materials and that you retain all copyright notices contained in the materials. By downloading content from our website, you accept the terms of this agreement.
VisualSimulationofSmokeRonaldFedkiwStanfordUniversityJosStam wavefrontHenrikWannJensenStanfordUniversityInthispaper,weproposeanewapproachtonumericalsmokesimulationforcomputergraphicsapplications.Themethodpro-posedhereexploitsphysicsuniquetosmokeinordertodesignanumericalmethodthatisbothfastandefÞcientontherelativelycoarsegridstraditionallyusedincomputergraphicsapplications(ascomparedtothemuchÞnergridsusedinthecomputational StanfordUniversity,GatesComputerScienceBldg.,Stanford,CA wavefront,1218ThirdAve,8thFloor,Seattle,WA98101,U.S.A.StanfordUniversity,GatesComputerScienceBldg.,Stanford,CAtheanimationofgasessuchassmoke.Weproposeamodelwhichisstable,rapidanddoesnÕtsufferfromexcessivenumericaldis-sipation.Thisallowsustoproduceanimationsofcomplexrollingsmokeevenonrelativelycoarsegrids(ascomparedtotheonesusedinCFD). toreducethenumericaldissipationinherentinsemi-Lagrangianschemes.WeachievethisbyusingatechniquefromtheCFDlit-eratureknownasÓvorticityconÞnementÓ[20].Thebasicideaistoinjecttheenergylostduetonumericaldissipationbackintotheßuidusingaforcingterm.ThisforceisdesignedspeciÞcallytoin-creasethevorticityoftheßow.Visuallythiskeepsthesmokealiveovertime.ThisforcingtermiscompletelyconsistentwiththeEu-lerequationsinthesensethatitdisappearsasthenumberofgridcellsisincreased.InCFDthistechniquewasappliedtothenumeri-calcomputationofcomplexturbulentßowÞeldsaroundhelicopterswhereitisnotpossibletoaddenoughgridpointstoaccuratelyre-solvetheßowÞeld.Thecomputationoftheforceonlyaddsasmallcomputationaloverhead.ConsequentlyoursimulationsarealmostasfastastheoneÕsobtainedfromthebasicStableFluidsalgorithm[17].Ourmodelremainsstableaslongasthemagnitudeoftheforcingtermiskeptbelowacertainthreshold.Withinthisrange,ourtimestepsarestillordersofmagnitudehigherthantheonesusedinexplicitschemes.Semi-Lagrangianschemesareverypopularintheatmosphericsciencescommunityformodelinglargescaleßowsdominatedbyconstantadvectionwherelargetimestepsaredesired,seee.g.[19]forareview.Weborrowfromthisliteratureahigherorderinter-polationtechniquethatfurtherincreasesthequalityoftheßows.ThistechniqueisespeciallyeffectivewhenmovingdensitiesandtemperaturesthroughthevelocityÞeld.Finally,ourmodel,likeFosterandMetaxasÕ[6],isabletohan-dleboundariesinsidethecomputationaldomain.Therefore,weareabletosimulatesmokeswirlingaroundobjectssuchasavirtualactor.Therestofthepaperisorganizedasfollows.Inthenextsectionwederiveourmodelfromtheequationsofßuidßow,andinsection3wediscussvorticityconÞnement.Insection4,weoutlineourimplementation.Insection5,wepresentbothaninteractiveandahighqualityphotonmapbasedrenderertodepictoursmokesimu-lations.Subsequently,insection6,wepresentsomeresults,whilesection7concludesanddiscussesfuturework.2TheEquationsofFluidFlowAttheoutset,weassumethatourgasescanbemodeledasinviscid,incompressible,constantdensityßuids.Theeffectsofviscosityarenegligibleingasesespeciallyoncoarsegridswherenumericaldissipationdominatesphysicalviscosityandmoleculardiffusion.WhenthesmokeÕsvelocityiswellbelowthespeedofsoundthecompressibilityeffectsarenegligibleaswell,andtheassumptionofincompressibilitygreatlysimpliÞesthenumericalmethods.Conse-quently,theequationsthatmodelthesmokeÕsvelocity,denotedbyu;v;w,aregivenbytheincompressibleEulerequations[14] Thesetwoequationsstatethatthevelocityshouldconservebothmass(Equation1)andmomentum(Equation2).Thequantitythepressureofthegasandaccountsforexternalforces.Alsowehavearbitrarilysettheconstantdensityoftheßuidtoone.Asin[7,6,17]wesolvetheseequationsintwosteps.FirstwecomputeanintermediatevelocityÞeldbysolvingEquation2overatimestepwithoutthepressureterm AfterthisstepweforcetheÞeldtobeincompressibleusingaprojectionmethod[3].ThisisequivalenttocomputingthepressurefromthefollowingPoissonequation withpureNeumannboundarycondition,i.e., atabound-arypointwithnormal.(Notethatitisalsostraightforwardtoim-poseDirichletboundaryconditionswherethepressureisspeciÞeddirectlyasopposedtospecifyingitsnormalderivative.)Theinter-mediatevelocityÞeldisthenmadeincompressiblebysubtractingthegradientofthepressurefromitWealsoneedequationsfortheevolutionofboththetempera-andthesmokeÕsdensity.Weassumethatthesetwoscalarquantitiesaresimplymoved(advected)alongthesmokeÕsvelocity (u BoththedensityandthetemperatureaffecttheßuidÕsvelocity.Heavysmoketendstofalldownwardsduetogravitywhilehotgasestendtoriseduetobuoyancy.WeuseasimplemodeltoaccountfortheseeffectsbydeÞningexternalforcesthataredirectlypropor-tionaltothedensityandthetemperature=(0pointsintheupwardverticaldirection,istheambienttemperatureoftheairandaretwopositiveconstantswithappropriateunitssuchthatEquation8isphysicallymeaningful.Notethatwhen,thisforceisEquations2,6and7allcontaintheadvectionoperatorAsin[17],wesolvethistermusingasemi-Lagrangianmethod[19].WesolvethePoissonequation(Equation4)forthepressureusinganiterativesolver.WeshowinSection4howthesesolverscanalsohandlebodiesimmersedintheßuid.3VorticityConÞnementUsuallysmokeandairmixturescontainvelocityÞeldswithlargespatialdeviationsaccompaniedbyasigniÞcantamountofrotationalandturbulentstructureonavarietyofscales.Nonphysicalnu-mericaldissipationdampsouttheseinterestingßowfeatures,andthegoalofournewapproachistoaddthembackonthecoarsegrid.Onewayofaddingthembackwouldbetocreatearandomorpseudo-randomsmallscaleperturbationoftheßowÞeldusingei-theraheuristicorphysicallybasedmodel.Forexample,onecouldgenerateadivergencefreevelocityÞeldusingaKolmogorovspec-trumandaddthistothecomputedßowÞeldtorepresentthemiss-ingsmallscalestructure(see[18]forsomeCGapplicationsoftheKolmogorovspectrum).Whilethisprovidessmallscaledetailtotheßow,itdoesnotplacethesmallscaledetailsinthephysicallycorrectlocationswithintheßowÞeldwherethesmallscaledetailsaremissing.Instead,thedetailsareaddedinahaphazardfashionandthesmokecanappeartobeÒaliveÓ,rollingandcurlinginanonphysicalfashion.ThekeytorealisticanimationofsmokeistomakeitlooklikeapassivenaturalphenomenaasopposedtoaÒlivingÓcreaturemadeoutofsmoke.OurmethodlooksforthelocationswithintheßowÞeldwheresmallscalefeaturesshouldbegeneratedandaddsthesmallscalefeaturesintheselocationsinaphysicallybasedfashionthatpro-motesthepassiverollingofsmokethatgivesittherealisticturbu-lentlookonacoarseCGgrid.Withunlimitedcomputingpower, Figure1:Discretizationofthecomputationaldomainintoidenticalvoxels(left).ThecomponentsofthevelocityaredeÞnedonthefacesofeachvoxel(right).anyconsistentnumericalmethodcouldbeusedtoobtainaccept-ableresultssimplybyincreasingthenumberofgridpointsuntilthedesiredlimitingbehaviorisobserved.However,inpractice,computationalresourcesarelimited,gridsarefairlycoarse(evencoarserinCGthaninCFD),andthediscretedifferenceequationsmaynotbeasymptoticallycloseenoughtothecontinuousequa-tionsforaparticularsimulationtobehaveinthedesiredphysicallycorrectfashion.Ourkeyideaistodesignaconsistentnumericalmethodthatbehavesinaninterestingandphysicallyplausiblefash-iononacoarsegrid.Ingeneral,thisisverydifÞculttodo,butluckilyavorticityconÞnementmethodwasrecentlyinventedbySteinhoff,seee.g.[20],forthenumericalcomputationofcomplexturbulentßowÞeldsaroundhelicopterswhereitisnotpossibletoaddenoughgridpointstoaccuratelyresolvetheßow.TheÞrststepingeneratingthesmallscaledetailistoidentifywhereitcomesfrom.Inincompressibleßow,thevorticityprovidesthesmallscalestructure.EachsmallpieceofvorticitycanbethoughtofasapaddlewheeltryingtospintheßowÞeldinaparticulardirection.ArtiÞcialnumericaldissipationdampsouttheeffectofthesepaddlewheels,andthekeyideaistosimplyadditback.Firstnormalizedvorticitylocationvectors thatpointfromlowervorticityconcentrationstohighervorticityconcentrationsarecomputed.ThenthemagnitudeanddirectionofthepaddlewheelforceiscomputedasisusedtocontroltheamountofsmallscaledetailaddedbackintotheßowÞeldandthedependenceonthespatialguaranteesthatasthemeshisreÞnedthephysicallycorrectsolutionisstillobtained.ThistechniquewasinventedbySteinhoffabout10yearsagowithaformsimilartoEquation11withoutthedependenceonseeforexample[20].ThismethodhasbeenusedsuccessfullyasanengineeringmodelforverycomplexßowÞelds,suchasthoseasso-ciatedwithrotorcraft,whereonecannotcomputationallyaffordtoaddenoughgridpointstoresolvetheimportantsmallscalefeaturesoftheßow.4ImplementationWeuseaÞnitevolumespatialdiscretizationtonumericallysolvetheequationsofßuidßow.AsshowninFigure1,wediceupthe Figure2:Semi-LagrangianpathsthatendupinaboundaryvoxelareclippedagainsttheboundariesÕface.computationaldomainintoidenticalvoxels.Thetemperature,thesmokeÕsdensityandtheexternalforcesaredeÞnedatthecenterofeachvoxelwhilethevelocityisdeÞnedontheappropriatevoxelfaces(seeFigure1,right).NoticethatthisarrangementisidenticaltothatofFosterandMetaxas[6]butdiffersfromtheoneusedbyStam[17]wherethevelocitywasdeÞnedatthevoxelcentersaswell.OurstaggeredgridarrangementofthevelocityÞeldgivesim-provedresultsfornumericalmethodswithlessartiÞcialdissipation.SeeappendixAformoredetailsonourdiscretization.Tohandleboundariesimmersedintheßuidwetagallvoxelsthatintersectanobjectasbeingoccupied.Alloccupiedvoxelcellfaceshavetheirvelocitysettothatoftheobject.Similarly,thetemperatureatthecenteroftheoccupiedvoxelsissettotheobjectÕstemperature.Consequentlyananimatorcancreatemanyinterestingeffectsbysimplymovingorheatingupanobject.ThesmokeÕsdensityisofcourseequaltozeroinsidetheobject.However,toavoidasuddendrop-offofthedensityneartheobjectÕsboundary,wesetthedensityatboundaryvoxelsequaltothedensityoftheclosestunoccupiedvoxel.Oursolverrequirestwovoxelgridsforallphysicalquantities.WeadvanceoursimulationbyupdatingonegridfromtheotheroveraÞxedtimestep.Attheendofeachtimestepweswapthesegrids.Thegridmayinitiallycontainsomeuserprovideddata,butinmostcasesthegridsaresimplyempty.WeÞrstupdatethevelocitycomponentsoftheßuid.Thisisdoneinthreesteps.First,weaddtheforceÞeldstothevelocitygrid.TheforcesincludeusersuppliedÞelds,thebuoyancyforcedeÞnedbyEquation8andthenewconÞnementforcedeÞnedbyEquation11.Thisisdonebysimplymultiplyingeachforcebythetimestepandaddingittothevelocity(seeappendixA).NextwesolvefortheadvectionterminEquation3.Wedothisusingasemi-Lagrangianscheme,see[19]forareviewand[17]foritsÞrstapplicationincomputergraphics.Thesemi-LagrangianalgorithmbuildsanewgridofvelocitiesfromtheonesalreadycomputedbytracingthemidpointsofeachvoxelfacethroughthevelocityÞeld.Newvelocitiesarethenin-terpolatedatthesepointsandtheirvaluesaretransferredtothefacecellstheyoriginatedfrom.Itispossiblethatthepointendsupinoneoftheoccupiedvoxels.InthiscasewesimplyclipthepathagainstthevoxelboundaryasshowninFigure2.Thisguaranteesthatthepointalwaysliesintheunoccupiedßuid.Simplelinearinterpola-tioniseasytoimplementandcombinedwithournewconÞnementforcegivessatisfactoryresults.Itisalsounconditionallystable.Higherorderinterpolationschemesare,however,desirableinsomecasesforhighqualityanimations.Thetrickypartwithhigheror-derschemesisthattheyusuallyovershootthedatawhichresultsin instabilities.InappendixBweprovideacubicinterpolatorwhichdoesnotovershootthedata.FinallyweforcethevelocityÞeldtoconservemass.Asalreadystatedinsection2,thisinvolvesthesolutionofaPoissonequationforthepressure(Equation4).Thediscretizationofthisequationresultsinasparselinearsystemofequations.WeimposefreeNeu-mannboundaryconditionsattheoccupiedvoxelsbysettingthenor-malpressuregradientequaltozeroattheoccupiedboundaryfaces.Thesystemofequationsissymmetric,andthemostnaturallinearsolverinthiscaseistheconjugategradientmethod.Thismethodiseasytoimplementandhasmuchbetterconvergencepropertiesthansimplerelaxationmethods.ToimprovetheconvergenceweusedanincompleteCholeskipreconditioner.Thesetechniquesareallquitestandardandwereferthereadertothestandardtext[11]formoredetails.Inpracticewefoundthatonlyabout20iterationsofthissolvergaveusvisuallyacceptableresults.Afterthepres-sureiscomputed,wesubtractitsgradientfromthevelocity.SeeappendixAfortheexactdiscretizationoftheoperatorsinvolved.Afterthevelocityisupdated,weadvectboththetemperatureandthesmokeÕsdensity.Wesolvetheseequationsusingagainasemi-Lagrangianscheme.Inthiscase,however,wetracebackthecentersofeachvoxel.Theinterpolationschemeissimilartothevelocitycase.5RenderingForeverytimestep,oursimulatoroutputsagridthatcontainsthesmokeÕsdensity.Inthissectionwepresentalgorithmstorealisti-callyrenderthesmokeundervariouslightingconditions.Wehaveimplementedbotharapidhardwarebasedrendererasin[17]andahighqualityglobalilluminationrendererbasedonthephotonmap[12].ThehardwarebasedrendererprovidesrapidfeedbackandallowsananimatortogetthesmoketoÒlookrightÓ.Themoreex-pensivephysics-basedrendererisusedattheendoftheanimationpipelinetogetproductionqualityanimationsofsmoke.WeÞrstbrießyrecalltheadditionalphysicalquantitiesneededtocharacterizetheinteractionoflightwithsmoke.Theamountofinteractionismodeledbytheinverseofthemeanfreepathofaphotonbeforeitcollideswiththesmokeandiscalledtheextinc-tioncoefÞcient.TheextinctioncoefÞcientisdirectlyrelatedtothedensityofthesmokethroughanextinctioncross-section.Ateachinteractionwiththesmokeaphotoniseitherscatteredorabsorbed.Theprobabilityofscatteringiscalledthe.Avalueofthealbedonearzerocorrespondstoverydarksmoke,whileavaluenearunitymodelsbrightgasessuchassteamandclouds.Ingeneralthescatteringoflightinsmokeismostlyfocusedintheforwarddirection.Thedistributionofscatteredlightismodeledthroughaphasefunctionwhichgivestheprobabilitythatanincidentphotonisdeßectedbyanangle.AconvenientmodelforthephasefunctionistheHenyey-Greensteinfunction (1+wherethedimensionlessparametermodelstheanisotropyofthescattering.Valuesnearunityofthisparametercorrespondtogaseswhichscattermostlyintheforwarddirection.Wementionthatthisphasefunctionisquitearbitraryandthatotherchoicesarepossible[1].5.1Hardware-BasedRendererInourimplementationofthehardware-basedrenderer,wefollowthealgorithmoutlinedin[17].InaÞrstpass,wecomputetheamountoflightthatdirectlyreacheseachvoxelofthegrid.ThisisachievedusingafastBresenhamlinedrawingvoxeltraversalal-gorithm[8].Initiallythetransparenciesofeachrayaresettoone).Then,eachtimeavoxelishit,thetransparencyiscom-putedfromthevoxelÕsdensity:=exp(,wherethegridspacing.ThenthevoxelÕsradianceissettowhilethetransparencyoftherayissimplymultipliedbythevoxelÕstransparency:.SincethetransparencyoftheraydiminishesasittraversesthesmokeÕsdensitythispasscorrectlymimicstheeffectsofself-shadowing.Inasecondpasswerenderthevoxelgridfromfronttoback.Wedecomposethevoxelgridintoasetoftwo-dimensionalgrid-slicesalongthecoordinateaxismostalignedwiththeviewingdirection.Theverticesofthisgrid-slicecorrespondtothevoxelcenters.Eachsliceisthenrenderedasasetoftransparentquads.Thecolorandopacityateachvertexofaquadcorrespondtotheradianceandtheopacity,respectively,ofthecorrespondingvoxel.Theblendingbetweenthedifferentgridsliceswhenrenderedfromfronttobackishandledbythegraphicshardware.5.2PhotonMapRendererRealisticrenderingofsmokewithahighalbedo(suchaswaterva-por)requiresafullsimulationofmultiplescatteringoflightinsidethesmoke.Thisinvolvessolvingthefullvolumerenderingequa-tion[2]describingthesteady-stateoflightinthepresenceofpar-ticipatingmedia.Forthispurposeweusethephotonmappingal-gorithmforparticipatingmediaasintroducedin[12].ThisisatwopassalgorithminwhichtheÞrstpassconsistsofbuildingavolumephotonmapbyemittingphotonstowardsthemediumandstoringtheseastheyinteractwiththemedium.Weonlystorethephotonscorrespondingtoindirectillumination.Intherenderingpassweuseaforwardraymarchingalgorithm.Wehavefoundthistobesuperiortothebackwardraymarchingalgorithmproposedin[12].TheforwardraymarchingalgorithmallowsforamoreefÞcientcullingofcomputationsinsmokethatisobscuredbyothersmoke.InadditionitenablesamoreefÞcientuseofthephotonmapbyallowingustouselessphotonsinthequeryastheraymarchergetsdeeperintothesmoke.Ourforwardraymarcherhastheformistheopticaldepth,isthefractionoftheinscatteredradiancethatisscatteredindirectionthesizeofthethstep,isarandomlychosenlocationinthethsegment.Thefactorcanbeconsideredtheweightofthethsegment,andweusethisvaluetoadjusttherequiredaccuracyofthecomputation.Thecontributionduetoin-scatteredradiance,,isgivenbyx;~)= x;~Wesplittheinscatteredradianceintoasinglescatteringterm,andamultiplescatteringterm,.Thesinglescatteringtermiscomputedusingstandardraytracing,andthemultiplescatteringtermiscomputedusingthevolumeradianceestimatefromthepho-tonmapbylocatingthenearestphotons.Thisgives:x;~x;~ 4 isthepowerofthethphotonandisthesmallestsphereenclosingthe Figure6:Twostillsfromtherotoranimation.Aboxisrotatinginsidethesmokecloudcausingittodisperse.Noticehowthesmokeissuckedinverticallytowardstheboxasitispushedoutwardshorizontally.Thesimulationtimefora120x60x120gridwasroughly60seconds/frame. Figure3:Risingsmoke.Noticehowthevorticiesarepreservedinthesmoke.Thesimulationtimefora100x100x40gridwasroughly30seconds/frame. Figure4:Lowalbedosmokepassingthroughseveralobjects.Eachobjectinteractswiththesmokeandcauseslocalturbulenceandvor-ticity.Thesimulationtimefora160x80x80gridwasroughly75 Figure5:Risingsmokeswirlingaroundasphere.Noticehowthesmokecorrectlymovesaroundthesphere.Thesimulationtimefora90x135x90gridwasroughly75seconds/frame. Figure7:Sixframesrenderedusingourinteractivehardwareren-dererofthesmoke.Thesimulationtimefora40x40x40gridwasroughly1second/frame. Figure8:Comparisonoflinearinterpolation(top)andournewmonotoniccubicinterpolation(bottom).Thesimulationtimefora20x20x40gridwasroughly0.1second/frame(linear)and1.8sec-onds/frame(thirdorder).6ResultsThissectioncontainsseveralexamplesofsmokesimulations.Wehaverunmostofthesimulationsincludingtherenderingonadual-Pentium3-800orcomparablemachine.TheimagesinFigures3-6havebeenrenderedatawidthof1024pixelsusing4samplesperpixel.Thesephotonmaprenderingsweredoneusing1-2millionphotonsinthevolumephotonmapandtherenderingtimesforallthephotonmapimagesare20-45minutes.Figure3isasimpledemonstrationofsmokerising.Theonlyexternalforceonthesmokeisthenaturalboyancyofthesmokecausingittorise.Noticehoweventhissimplecaseisenoughtocreatearealisticandswirlyapperanceofthesmoke.Figures4and5demonstratethatoursolvercorrectlyhandlestheinteractionwithobjectsimmersedinthesmoke.Theseobjectsneednotbeatrest.Figure6showstwostillsfromananimationwherearotatingcubeisinsideasmokecloud.Therotationofthecubecausesthesmoketobepushedouthorizontallyandsuckedinvertically.ThegridresolutionsandthecostofeachtimesteparereportedintheÞgureFigure7showssixframesofananimationrenderedusingourinteractiverenderer.TherenderingtimeforeachframewaslessthanasecondonanVidiaQuadrographicscard.Thespeed,whilenotreal-time,allowedananimatortointeractivelyplacedensitiesandheatsourcesinthesceneandwatchthesmokeraiseandbillow.Finally,Figure8demonstratesthebeneÞtsofusingahigheror-derinterpolantinthesemi-Lagrangianscheme.Thethreepicturesonthetopshowtheappearanceoffallingsmokeusingalinearin-terpolant,whilethepicturesonthebottomshowthesamesmokeusingournewmonotoniccubicinterpolant.Clearlythenewinter-polationreducestheamountofnumericaldissipationandproducessmokesimulationswithmoreÞnedetail.7ConclusionsInthispaperweproposedanewsmokemodelwhichisbothstableanddoesnotsufferfromnumericaldissipation.Weachievedthisthroughtheuseofanewforcingtermthataddsthelostenergybackexactlywhereitisneeded.Wealsoincludedtheinteractionofobjectswithoursmoke.WebelievethatourmodelisidealforCGapplicationswherevisualdetailandspeedarecrucial.WethinkthatvorticityconÞnementisaveryelegantandpower-fultechnique.Weareinvestigatingvariantsofthistechniquecus-tomtailoredforotherphenomenasuchasÞre.Wearealsoinvesti-gatingtechniquestoimprovetheinteractionofthesmokewithob-jects.Inourcurrentmodelobjectsmaysometimesbetoocoarselysampledonthegrid.8AcknowledgementsWewouldliketothankJohnSteinhoff(FlowAnalysisInc.andUTSI)andPatHanrahan(StanfordUniversity)formanyhelpfuldiscussions.TheworkoftheÞrstauthorwassupportedinpartbyONRN00014-97-1-0027.Theworkofthelastauthorwassup-portedbyNSF/ITR(IIS-0085864)andDARPA(DABT63-95-C-ADiscretizationWeassumeauniformdiscretizationofspaceintovoxelswithuniformspacing.ThetemperatureandthesmokeÕsdensityarebothdeÞnedatthevoxelcentersanddenotedbyi;j;k;N;respectively.ThevelocityontheotherhandisdeÞnedatthecellfaces.ItisusualintheCFDliteraturetousehalf-wayindexnotationforthis;N;j;k;N;;N;i;k;N;;N;i;j;N:UsingthesenotationswecannowdeÞnesomediscreteoperators.ThedivergenceisdeÞnedaswhilethediscretegradientsare(noteThediscreteLaplacianissimplythecombinationofthedivergenceandthegradientoperators.ThediscreteversionofthevorticityisdeÞnedasfollows.Firstwecomputethecell-centeredvelocitiesthroughaveraging Figure9:StandardcubicHermiteinterpolation(left)producesovershootswhileourmodiÞedinterpolationscheme(right)guar-anteesthatnoovershootsoccur.AllofourforceÞeldsaredeÞnedatthecenterofthegridvoxels.Togetvaluesatthefaceswesimplyaverageagain.IftheforceÞeld,thenthevelocityisupdatedasBMonotonicCubicInterpolationInthisappendixwepresentacubicinterpolationschemewhichdoesnotovershootthedata.Sinceourvoxelgridsareregularthethree-dimensionalinterpolationcanbebrokendownintoase-quenceofone-dimensionalinterpolationsalongeachcoordinateaxis.Therefore,itissufÞcienttodescribetheone-dimensionalcaseonly.ThedataconsistsofasetofvaluesdeÞnedatthelocations.Avalueatapointpointtk;tk+1]canbeinterpolatedusingaHermiteinterpolantasfollows[8]However,thisinterpolantusuallyovershootsthedataasweshowonthelefthandsideofFigure9.Wewanttoavoidthis,sincemonotoneinterpolationguaranteesstability.Onesolutionistosim-plycliptheinterpolationagainstthedata,butthisresultsinsharpdiscontinuities.Anotherremedyistoforcetheinterpolanttobemonotonicovereachintervalaltk;tk+1].AnecessaryconditionforthistobethecaseisthatsignsignsignInourimplementationweÞrstcomputeandthensettheslopestozerowhenevertheyhaveasigndifferentfrom.OntherighthandsideofFigure9,weshowtheournewinterpolantappliedtothesamedata.ClearlytheovershootingproblemisÞxed.References[1]P.Blasi,B.LeSaec,andC.Schlick.ARenderingAlgo-rithmforDiscreteVolumeDensityObjects.ComputerGraph-icsForum(EUROGRAPHICS93ConferenceProceedings)12(3):201Ð210,1993.[2]S.Chandrasekhar.RadiativeTransfer.Dover,NewYork,[3]A.Chorin.ANumericalMethodforSolvingIncompressibleViscousFlowProblems.JournalofComputationalPhysics2:12Ð26,1967.[4]Y.Dobashi,K.Kaneda,T.Okita,andT.Nishita.ASimple,EfÞcientMethodforRealisticAnimationofClouds.InGRAPH2000ConferenceProceedings,AnnualConference,pages19Ð28,July2000.[5]D.S.EbertandR.E.Parent.RenderingandAnimationofGaseousPhenomenabyCombiningFastVolumeandScanlineA-bufferTechniques.ComputerGraphics(SIGGRAPH90ConferenceProceedings),24(4):357Ð366,August1990.[6]N.FosterandD.Metaxas.RealisticAnimationofLiq-GraphicalModelsandImageProcessing,58(5):471Ð483,1996.[7]N.FosterandD.Metaxas.ModelingtheMotionofaHot,TurbulentGas.InSIGGRAPH97ConferenceProceedings,AnnualConferenceSeries,pages181Ð188,August1997.[8]J.D.Fowley,A.vanDam,S.K.Feiner,andJ.F.Hughes.ComputerGraphics:PrinciplesandPractice.SecondEdi-.Addison-Wesley,Reading,MA,1990.[9]M.N.Gamito,P.F.Lopes,andM.R.Gomes.Two-dimensionalSimulationofGaseousPhenomenaUsingVor-texParticles.InProceedingsofthe6thEurographicsWork-shoponComputerAnimationandSimulation,pages3Ð15.Springer-Verlag,1995.[10]G.Y.Gardner.VisualSimulationofClouds.puterGraphics(SIGGRAPH85ConferenceProceedings)19(3):297Ð384,July1985.[11]G.GolubandC.VanLoan.MatrixComputations.TheJohnHopkinsUniversityPress,Baltimore,1989.[12]H.W.JensenandP.H.Christensen.EfÞcientSimulationofLightTransportinSceneswithParticipatingMediausingPhotonMaps.InSIGGRAPH98ConferenceProceedings,AnnualConferenceSeries,pages311Ð320,July1998.[13]J.T.KajiyaandB.P.vonHerzen.RayTracingVolumeDen-ComputerGraphics(SIGGRAPH84ConferencePro-,18(3):165Ð174,July1984.[14]L.D.LandauandE.M.Lifshitz.FluidMechanics,2ndedi-.Butterworth-Heinemann,Oxford,1998.[15]K.Perlin.AnImageSynthesizer.ComputerGraphics(SIG-GRAPH85ConferenceProceedings),19(3):287Ð296,July[16]G.Sakas.FastRenderingofArbitraryDistributedVolumeDensities.InF.H.PostandW.Barth,editors,ProceedingsofEUROGRAPHICS90,pages519Ð530.ElsevierSciencePublishersB.V.(North-Holland),September1990. [17]J.Stam.StableFluids.InSIGGRAPH99ConferencePro-ceedings,AnnualConferenceSeries,pages121Ð128,August[18]J.StamandE.Fiume.TurbulentWindFieldsforGaseousPhenomena.InSIGGRAPH93ConferenceProceedings,An-nualConferenceSeries,pages369Ð376,August1993.[19]A.StaniforthandJ.Cote.Semi-lagrangianintegrationschemesforatmosphericmodels:Areview.MonthlyWeatherReview,119:2206Ð2223,1991.[20]J.SteinhoffandD.Underhill.ModiÞcationoftheeulerequa-tionsforÒvorticityconÞnementÓ:Applicationtothecomputa-tionofinteractingvortexrings.PhysicsofFluids,6(8):2738Ð2744,1994.[21]L.YaegerandC.Upson.CombiningPhysicalandVisualSim-ulation.CreationofthePlanetJupiterfortheFilm2010.puterGraphics(SIGGRAPH86ConferenceProceedings)20(4):85Ð93,August1986.[22]G.Yngve,J.OÕBrien,andJ.Hodgins.Animatingexplosions.SIGGRAPH2000ConferenceProceedings,AnnualCon-ferenceSeries,pages29Ð36,July2000.