基于navier-stokes方程的煙霧湍流效果模擬導(dǎo)師版

結(jié)合kkSmoketurbulenceeffectSimulationBasedonNavier-StokesEquationMajor:ComputerApplicationsPostgraduate:Xiang Tutor:ZhangTurbulenceisthemostcommonphenomenoninthenature.Itpresentsabundantfluiddetails.Thesimulationofturbulencefluidiswidelyusedinsystemsimulation,computergameeffects,movieeffectsandsomeotherfields.Itisverydifficulttosimulateturbulencefluidbasedonphysicalmethods.Althoughitcanproducerealeffects,contradictionbetweenreal-timeandrealityisfarbeyondsolved.Thispapertriestoimprovesmokerealityinconditionofkeeitseffectswhensimulatingsmoketurbulence.Thispapermakesastudyofturbulencealgorithms,yzesitsgeneratingconditionsandfactors,establishesasimulationalgorithmcombinedgridsystemwithparticlesystemaftersufficientresearchesonrelatednoise,turbulencemodels,turbulenceenergyspectrumandturbulencefluid,andreducescouplingbetweenparticlesystemandgridsystemwhencomputingturbulence.Therefore,itmeetsbothrealityandefficiencyofturbulenceatthesametime.ThecorealgorithmofthisisthattosimulateturbulencedetailsofsmokefluidaccuraybykmodelcombinedwithKolmogorovlaw，efficientparticlesystemandwaveletnoisetechnology.Comparedwithexistingalgorithm,theadvantagesofthispaperareasfollows:itreducedcouplingbetweenturbulenceparticlesystemandgridsystem,ensuredrealityandreal-timeofturbulenceflowatthesametimeandmadeitsreal-timeapplicationoffluidturbulencepossible.ThemaincontentofthispaperDiscardingPerlinNoise,itusesWavletnoiseinstead.SoitcouldedetaillossandaliasingofPerlinNoisewhenrenderingscene.

k

simulation,itsynthesisturbulencewithenergyspectrumdefinedbyKolmogorov’slaw.Finally,itcombinesgridsystemandparticlesystemtogether.Particlesk

turbulencemodelonlyneedtosamplevelocitiesofgridsystemtheircorrespondingpositions.Inthisway,itreducescouplingbetweengridsystemandparticlesystemandgeneratesturbulenceinlowergridresolutionsuccessfully.Andthenitcancontrolturbulencescaleartificiallybydecidingnumberofparticles.Whilemaintainingsmokeefficiencyitimprovessmokerealitybyaddingturbulenceeffect.TheexperimentalresultsshowthatitcanmeetapplicationrequirementbasedusingNavier-Stokeequationcombinedwithk

modelparticles.Whenresolutionchanging,turbulenceeffectisstillstable.Framerateisalsoimprovedcomparedwithhighresolutiongridmethodatthesameeffect.:Turbulence,particle,k

model,grid,reality, 歷史發(fā) 湍流方 進 湍流模擬的難 Helmholtz-Hodge分 計算平流 計算粘度擴散 初始化和邊界條 優(yōu)缺 湍流粘度恒定模 混合尺寸模 單方程模 正交 上采 小波和下采 構(gòu)建噪聲 小數(shù)部分的縮放與平 在空間中的噪 投影噪 噪聲方 噪聲求 湍流k-ε模 能量頻 湍流大小尺 細節(jié)湍流合成的重要性 湍流速度場的合 穩(wěn)定 初始狀 湍流相干 小波噪 后處理方法流動效 對比實驗 對比實驗二 對比實驗 圖1-1現(xiàn)實天空中的湍 圖2-1JosStam的網(wǎng)格半拉格朗日平流 圖2-2粒子的平 圖2-3求解NS方程示意 圖3-1小波噪聲四個步 圖3-2可改善的均勻二次B樣條曲 圖3-3一維上的小波噪 圖3-4-2≤j≤3時一維影響函數(shù)的衰減過 圖4-1能量的流 圖4-2頻譜、噪聲與速度場的結(jié) 圖4-3湍流的合成示意 圖4-4算法流程示意 圖4-5煙霧后處理方 圖4-6紋理坐標t1或t2變形示意 圖5-1兩種噪聲的比 圖5-2后處理方法下的煙霧效果 圖5-3后處理方法下粒子對應(yīng)時刻的速度 圖5-4湍流運 圖5-5整個煙霧場景 圖5-6后處理方法第100幀煙霧不同視 圖5-732*32*32普通網(wǎng)格方法求解N-S方程三維煙 圖5-864*64*64普通網(wǎng)格方法求解N-S方程三維煙 圖5-9不同分辨率下的煙 圖5-10不同分辨率下的煙霧湍流效率走勢 圖5-11三種煙霧不同條件下幀率對 無論在自然界或者是中，煙霧、水和火的效果隨處可見。它們不可思1-1所示，天空中火山或者是大型的煙氣羽流中，可以看到不同規(guī)模大小從數(shù)米到數(shù)百米的的漩渦，它們的運動極其復(fù)雜。這些復(fù)雜和炫目的運動使得它們成為了和電腦游戲中最具性的部分。1-1動時產(chǎn)生明顯的變化，所以湍流流動是難以控制和預(yù)測的，這是流體模擬雖然的畫面可以展示出相當真實的動畫場景，但是流體模擬中實時性和真實性之間的還遠遠沒有解決。為了獲得高質(zhì)量的視覺效果，動畫的運用，但是流體實時性和真實性之間的這個基本問題卻依然沒有解決。結(jié)流體模擬在計算機圖形學(xué)和計算流體力學(xué)中的發(fā)展歷史和進展，并且主[6]（FVM部分，但是它的缺點在于需要前期的大量計算和分離網(wǎng)格這個多余的步驟。并結(jié)合可壓縮的內(nèi)核來模擬了流體力學(xué)[8]。Müller等人用這個方法解決了自由SPH在計算機圖形學(xué)領(lǐng)域得以流行[9]。但是為了保證粒子在合，流體流線變得模糊，直到流線徹底，流動中形成了許多小漩渦，此時[13][14]的情況下，此模型用一系列微分方程來模擬湍流的動態(tài)狀態(tài)。在航空領(lǐng)域中廣t-maras模型15直接用一個偏微分方程（）來模擬湍流粘（twequationmdelk-ε模型16k-ω模型17到現(xiàn)在仍然是在工程上和計算流體力學(xué)中使用最為廣泛的湍流模型。但是，這esChen199418它依賴于雷諾茲平均數(shù)（Reynoldsaveraging）來分離平均流動和湍流，（coherenteddies）的時候。鑒于此，Smagorinsky在湍流模型[19]中引入了大渦模擬（LargeEddySimulations，LES。LES模型用更精確的頻譜濾波器來分離平均流動和湍流，它可以直接模擬相干漩渦。LES模型最初用于氣象學(xué)中描述空氣的內(nèi)部流動，但是現(xiàn)在經(jīng)過改進可以用于計算流體力學(xué)的諸多應(yīng)用中了[20]。LES模型的一個改進是分離了渦模擬（DetachedEddySimulations(DES)）[21]，這個模型更經(jīng)典模型所提供的信息。當湍流粘度足以描述平均流動的虛擬擴散時，（Reynoldsrestraortmodels）解決了上面的這個問題。它不僅模擬湍流粘度，而且還（turblentstensor2223lLangeinEquation）24]這樣的密度方程直接使用粒子方法描述了湍流的高階統(tǒng)計學(xué)性質(zhì)。而橢圓松弛方法（ElliticRelaxtion）25卻在湍流模型中用全局優(yōu)化法代替了流動局部對流擴散的處理方法。雖然，所有的這些方法比經(jīng)典的湍流模型在理論上擁有更高階的預(yù)測精度，但是它們的計流體力學(xué)的研究領(lǐng)域和高品質(zhì)的模擬仿真中。隨著近些年來計算機計算能力的不斷增長，S和雷諾茲應(yīng)力傳輸模型的應(yīng)用逐漸變得廣泛起來。7080年代，計算流體力學(xué)中的湍流模型的基本原理陸續(xù)在湍流比較成領(lǐng)域。近些年來，比較流行的湍流話題有可壓縮流體中的湍現(xiàn)多尺寸湍流模擬技術(shù)[29]，以及用自適應(yīng)尺寸來模擬湍流（scaleadaptivesimulation(SAS)）[30]?？聽柲觐l譜（Kolmogorovspectrum）產(chǎn)生的無散度湍流引入了湍流合成中[31]這個方法曾經(jīng)被用來模擬的[32]和火焰的燃燒[33]Bridson提出用Kim用小波分解來確定局部湍流的強度[35]，并且使用與小波頻譜匹配的旋度噪湍流會存在問題。Schechter[36]Narain[37]用傳輸模型來導(dǎo)出湍流參數(shù)。這雖然條件假設(shè)，比如：用不同尺度混合模型來近它們的單方程（one-equation）能ThreKim提出了一種把網(wǎng)格和粒子以及兩方程模型結(jié)合的算法[38]。在下面主要討論湍流合成的進展。Navier-Stokes方程（Reynolds-averagedNavier-Stokesequation，RANS）湍流pdf方法在粒子系統(tǒng)上解出了郎之萬方程[40]。最后，更新的速度直速度作為實體粒子的運動速度時，合成的運動既不是無散度的又不是空間相干的。這產(chǎn)生真實漩渦或者是湍流的必備條件。另外，Zho提出了作用于上采樣網(wǎng)格上的隨機外力來合成湍流的方法41]。他可以通過預(yù)計算一個給定的能量頻譜來得到無散度的外力場，這種方法不需種方法用于尺寸稍大一些的湍流中，這種暫時的不一致將會導(dǎo)致湍流中出現(xiàn)人有的思路。Th remm發(fā)現(xiàn)當流動被自由ClosesttMethod（M）[42]技術(shù)成功地跨過了這個生成液體表面時所的巨大。CM運算符不是天然的三維生成形態(tài)，因此，可以利用逆阿貝爾變換（ivereltranform減少巨大的模糊性。最后，利用兩個湍流耦合方法來生成高分辨率的液體表面波浪。從上文湍流方法介紹和實際中所遇到的可以看出，湍流模擬的難點主湍流噪聲的選取。柏林噪聲（Perlinnoise）是計算機圖形學(xué)領(lǐng)域中最為湍流流動模型和能量譜的選擇。眾所周知，湍流流動模型有許多種類，選擇模型對湍流進行模擬直接關(guān)系到湍流流動的最終效果，而能量譜的正確選擇能夠準確反映湍流的能量變化。因此，需要選擇正確的湍流流動模型和能量譜；并不能夠?qū)崟r和真地模擬出湍流的效果為此后來又引入了兩方程模型但實性之間的卻沒有有效地得到解決。絕大多數(shù)時候，效果的實現(xiàn)者不得不Navier-StokesJosStam求Navier-Stokes方程的經(jīng)典歐拉網(wǎng)格法方法為基礎(chǔ)，并且把流體的不可壓縮性Navier-StokesGPU上求解出加入了限Navier-Stokes方程，這是本文的基礎(chǔ)。運用小波噪聲（waveletnoise）技術(shù)替代了傳統(tǒng)流體模擬中所使用的在流體模擬中引入了k湍流模型這個兩方程模型，它使用了與湍流粘度類似的定義，即混合長度模型，但是它卻用了湍流勢能項k和湍流擴散項來定義。k湍流模型可以有效地模擬小規(guī)模的湍流，這些湍流無法直在結(jié)合前面（1）、（2）、（3）GPU強大的并行第六章將對本文工作做一個全面的總結(jié)，并本文所用方法存在的不足其中，流體的速度是最為重要的一個變量，因為它決定了流體本身以及它所攜2-1所示。起始階 2-1JosStam

（2-（2-Navier-Stokes方程包括了流體動力學(xué)的三個基本控制方程：連續(xù)性方程、Claude-LouisNavierGeorgeGabrielStokes1821年和1845Navier-Stokes的連續(xù)性方程，如（2-1）所示，它是描述流動物質(zhì)運動的方程。其也十分重要；表示的是梯度算子；表示的是散度算子；2斯算子。（2-2）是流體的不可壓縮的約束條件2-22-2其它很多液體更加黏稠因此蜂蜜流動得十分緩慢而自來水或者卻流動JosStam1999MAC離散方法，引入了帶有一階壓力投影的半拉格和基礎(chǔ)。JosStamNavier-Stokes方程應(yīng)用于計算機圖形學(xué)流體模擬中，首先Helmholtz-Hodge分解成為（2-3）的形式：wu

（2-其中，u0，而且與D平行，即在Du?n=0（關(guān)于它的證明Stokes0的速度：uw

（2-根 3， （2-3）兩邊，可以得到w(up)u

（2-因為（2-2）中u0，所以（2-5）可以簡化為2p （2-程。這說明已知w，可以通過（2-6）得到p。所以，如（2-6）所示，還需要找到一個計算w的方法。眾所周知，有mnmnmn的一個單位矢mn方向上所得到的分量。因此，Helmholtz-HodgeP，它的意思是把一個矢wuP帶入方程（2-3）PPwPuuP(p)0Navier-Stokes方程了。把投影算符帶入（2-1）得：PuP((u)u1pv2uF

（2-u0，公式（2-7）左邊的導(dǎo)數(shù)也是無散度的，所以，有8 （2-w，再對它進行投影運算得到一個新的無散度的速度場u。只需要從（2-6）中解出p，然后再帶入（2-ADFP。

（2-2-3在計算平流時，JosStam把每個網(wǎng)格單元看成是一個粒子。所以，計算平流就一個時間段tr

（2-平流會造成流體在t比較大的時候不能夠保持其本身的穩(wěn)定性，如果u(t)t因此，JosStam使用了逆向和隱式積分法。他不是通過粒子“移動”來計再先前位置的量（這個量可以是任何需要求解的物理量。假設(shè)q是需要

（2-JosStam2-1uv （2-

（2-在方程（2-13）中，2是方程（2-3）拉斯算子的不連續(xù)形式。如同計算平流的顯式歐拉方法一樣，對于t和v的值比較大時就很不穩(wěn)定。根據(jù)JosStam的方法，對（2-13）使用隱式解法得：

（2-x(k

x(k

x(k

x(k

x(k1)i i i,j i,j ii

（2-其中，a和GPU-tokes方程每個時間步長中，都需要求解速度項和壓力項。本文使用了最簡單的邊界條件，把速度項和壓力項在邊界上都0，造成任何影響。Navier-Stokes2-3NSNavier-StokesNavier-Stokes方程中某些2-3求解NS如羅納德均勻Navier-Stokes方程等平均的模擬方法來說，雷諾茲應(yīng)力張量（Reynoldsstresstensor）是其中一個非常重要的變量。因此研究湍流的很和上的假設(shè)。這些方法可以按照建模變量的數(shù)量和它們的完整度來分類。VT為恒定，但該模型只能夠模擬簡單的流動。一種改進的方法是將粘度VT表示為多個不同尺寸大小的漩渦模型系數(shù)的混Smagorinsky、Baldwin和Lomax共同混合尺寸模型是模擬湍流中很vL2

tensor編碼，它的精確性取決于準確的尺寸值。Lm的解析式只能適用于幾種簡單類型Spalart-Allmaras單方程模型是湍流模型中比較完整使用非常廣泛的一種模 k

[(vT

j

ij

D

vTk Spalart-Allmaras對于大量的分離流動，無剪切流以及衰變湍流效果不好。章節(jié)中引入的兩方程k湍流模型做了必要的鋪墊?，F(xiàn)今計算機圖形學(xué)中使用最廣泛的噪聲是柏林噪聲[44]。它是由KenPerlinWiener插值法（Wienerinterpolationmethod）來構(gòu)建一個包含隨機頻譜的紋丟失和產(chǎn)生鋸齒之間的。 (b)下采樣圖像 (d)噪聲帶圖像N=R-3-1接下來，本文將就小波噪聲的算法展開詳細的討論小波噪聲的算法包括了以下的四個步驟，如圖3-1所示下采樣（downsampling）RR↓（3-上采樣（upsampling）RR↓↑（3-RRN

bmaxwN(2bx bb

（3-b整個頻帶值大小的典型的渲染使用濾波核（filterkernel）把一個場景方程轉(zhuǎn)換成像素。因此，i個像素的值為：Pixel(i)S(x)K(xi

（3-iK(x)參通常是非負的。場景方程S(x包含了一些典型的場景信息，比如位置信息和幾何片元的辨率噪聲的縮放和平移定位估計。如果Si個像素的縮放，k表示偏移，假設(shè)S(x)M(2s(xk))i個像素的值為： b

bj N(2jxl)K(xijbb

（3-為了標記的方便，在結(jié)果中，用Sbj，用2Sbk替換l。在數(shù)學(xué)上，j1上沒j0j0

（3-這就是小波噪聲N(x與濾波核K(x)的正交性。本文小波噪聲的實現(xiàn)就例如：B樣條曲線（B-splinecurves）B樣條基礎(chǔ)方程累加起來假設(shè)基礎(chǔ)方程(x)x=0i就可F(x)fi(xii其中，fi表示基礎(chǔ)函數(shù)(xi)0x20S0S0(x){F(x)|F(x f(xi

（3-（3-G(x)gi(2xii

（3-小波分析中最的思想是所有在S0空間中的函數(shù)同樣都在S1空間中。換言之，S0S1，S1完全包含S0(x)和(2xk)之間存在一個系數(shù)Pk滿足以下的條件(x)pk(2xk （3-（refinable均勻B樣條基礎(chǔ)方程正弦方程daubechies基礎(chǔ)方程而函（Gaussians）110 - 3-2B3-2BB數(shù)p，如（3-5）所示。那么存在一個在空間S1中代表方程的系數(shù)pi有iF(x)f(2xiii

（3-8ifpi2kfkik

（3-其中，序列fi長度是序列fk長度的兩倍。因此， S0空間的基礎(chǔ)方程轉(zhuǎn)換到分辨率S1Algorithm3-1：BUpsample(*org,*des,n,*p,Coeffs[4]={0.25,0.75,0.75,p=foreachiin[0,n]des[i*stride]=foreachkin[0,i/2+1]des[i*stride]+=p[i-2*k]*return所有在S1中的方程可以完全用S0中的方程來表示，通常這里就會造成一些細節(jié)的損失。但是，可以把細節(jié)的損失最小化。如（3-11）所示：G(x)G(x)D(x

（3-其中G(x)是最接近G(x)D(x)是G(x)與G(x)D(x)越小G(x)與G(x)D(xl)(xi)dx

（3-由于D(x)在所有整數(shù)位置上與S0上所有的函數(shù)都是正交的。因此，細節(jié)函數(shù)D(x)組成的矢量空間被稱為小波空間W0。設(shè)已知g是函數(shù)G(x基礎(chǔ)函數(shù)的系數(shù)，而g是G(x 那么，一定存在一個分析系數(shù)ag igak2igkik

（3-gg 數(shù)a2B樣條曲線來說，序列 因為(x)(2jx)jN。（3-j的值不會對圖像結(jié)果造成任何影響。Algorithm3-2：B#defineNUMDownsample(*org,*des,n,*a,Coeffs[2*NUM]= 0.000334,-0.001528,0.000410, -0.000938,-0.008233,0.002172, -0.005040,-0.044412,0.011655, -0.025936,-0.243780,0.033979, 0.655340,0.033979,-0.243780,- 0.019120,0.002172,-0.008233,- 0.003546,0.000410,-0.001528,a=&Coeffs[NUM]foreachkin[0,i/2+1]foreachiin[0,n/2]des[i*stride]=foreachkin[2*i-NUM,2*i+NUM]des[i*stride]+=p[k-2*i]*現(xiàn)在構(gòu)建噪聲帶的工作已經(jīng)準備好了N(x)。因為小波空間W0是由濾波曲線基礎(chǔ)函數(shù)B(x)。原因如下：B它次數(shù)低、依賴小使得對N(x用隨機數(shù)創(chuàng)建R(x)：用一個隨機數(shù)創(chuàng)建一個系數(shù)序R(...,r,...)S1中定義了一個函數(shù)R(x)B(x)，i （3-7）變成了R(x)riB(2xi

把R(x)下采樣得到R(x)：正如 部分：一部分是R(x)在S0中，另一部分是細節(jié)部分N(x)，它在W0中。R(x)R(x)N(xN(x)R(x)R(x帶入（3-5）得iiir其中，系數(shù)r可以 （3-13）中系數(shù)計算得到ri

（3-（3-（3-irak2iik

（3-3-3（b）（把R(x)計算上采樣得到R(x)：把3-15）和（3-19）帶入（17）

（3-二項求和是在S0中使用基礎(chǔ)函數(shù)B(xi)要一個通用的基礎(chǔ)函數(shù)。這是可以實現(xiàn)的，通過上采樣可以在S1中精確地表示S0中的R(x)：iii其中，系數(shù)r可以通過上采樣系數(shù)r r r

（3- R(x)3-3（c）

i2k

（3-從R(x)中減去R(x)得到N(x

（3-niB(2xii

（3- nrr得到的。噪聲帶N(x)如圖3- B樣條曲線得到。3-3N(2jx)j0來說在S0中都是正交的。因此，不會對圖像造成任何影響。另外，噪聲帶彼此Algorithm3-3Bif *temp1,*temp2,foreachiin[0,sz]noise[i]=foreachiyin[0,n]foreachizin[0,n]i=iy*n+Downsample(&noise[i],&temp1[i],n,Upsample(&temp1[i],&temp2[i],n,foreachixin[0,n]foreachizin[0,n]i=ix*n+Downsample(&temp2[i],&temp1[i],n,nUpsample(&temp1[i],&temp2[i],n,nforeachixin[0,n]foreachiyin[0,n]i=ix+Downsample(&temp2[i],&temp1[i],n,n*nUpsample(&temp1[i],&temp2[i],n,n*nforeachiin[0,n*n*n]foreachixin[0,n]foreachiyin[0,n]foreachizin[0,n]temp1[i++]=noise[Mod(ix+offset,n)++Mod(iz+offset,n)*n*nforeachiin[0,n*n*n]噪聲帶N(x已經(jīng)被構(gòu)建出來了，因此當l是任意整數(shù)，j是非負整數(shù)N(2jxl對最后的像素值沒有任何影響。但是，實際上lj不總是整數(shù)，這j0C(j,l)N(2jxl)B(x（3-是一個關(guān)于ljjl1/201233-42≤j≤3j0在空間中的噪二維噪聲函數(shù)N(x,y)是由二維的系數(shù)矩陣 構(gòu)建成的，因此均勻B

ix

ix

nix

B(2xix)B(2yiy

（3-Nix,j

是從隨機數(shù)矩陣Rix,j

,...)首先，對Rix,j

,...)進行下采樣濾波，再對每一列R N(2jx,2jy)B(xi)B(yi)dxdy0,j

（3- N(2jx,2jy,2jz)B(xi)B(yi)B(xi)

（3-舉個例子，傅里葉切片理論（FourierSliceTheorem）闡明了如果噪聲帶的2jzzN(2jx,2jy,z)在(xy) (2jx,2jy)N(2jx,2jy,2jz)B(zz

（3- 28， (2jx,2jy)B(xi)B(yi)dxdy

（3-，，

N（x,y,z）而建立起來的。因此，頻帶的正交性也被保存了。只要ij那么Nz0(2ix,2iy)和Nz0(2jx,2jy)jj3-4從隨量ri到噪聲系數(shù)ni的每一步都是線性的，這表明噪聲系數(shù)ni是ri的線性組合。從統(tǒng)計學(xué)的角度來說，這樣獨立的線性組合與隨機分布變是等價的。令x1,…,xn是隨機分布變量，它的方差是2xnn

ywixii1

（3-2

nw

（3- x i1ni的方差2ni數(shù)值計算得到2 差與方差2nQPixel(i)BqN(xqQq

（3-N(xq)B(x-xq)x=xqniN(xq)。WNoiseVec inti,f[3],c[3],mid[3],floatw[3][3],t,result=foreachiin[0,3) mid[i]=ceil(p[i]-0.5));t=mid[i]–(p[i]- foreachf[2]in[-1,1]do{foreachf[1]in[-1,1]do{foreachf[0]in[-1,1] foreachf[0]in[-1,1]do{c[i]=Mod(mid[i]÷f[i],n);weightresult+=weight*return但是如果Q的采樣數(shù)量太少時，N(xq)替代N(xq)來獲得更精確的像素值N(xq)是N(x)在xq附近圈中的平均值當Q=1時那么有Pixel(i)N(i)， （3-2）的積分結(jié)果等價。計算N(xq)的最佳方法受到渲染器正交計算細節(jié)B(x)的影響。例如，當Q=1時，用B((12j)(xxq))，來對噪聲系數(shù)進j=06j=-14.5個噪聲系數(shù)普遍收縮的B(xxq)Algoritm3-4所示。湍流kεk（4-4這并不同于單方程模型，因為單方程模型沒有額外依賴的假設(shè)條件。兩方程模k和k表示表示的湍流結(jié)構(gòu)的耗散。k模型的偏微分方程表示如下：(kUkvTk)P(1 （4-1(UvT)(CPC(2 2

（4-（由4-1）和（4-2）可以看出，兩個 含了一個在平均流動中的平流項。兩 的等號右邊依次包含了粘度擴散（diffusion（dissipationU和生成項中的平均流動張力SPT

2（4-22vC2 u

（4-其中，S0.5*(U/xU/x)，C0.09 VT是一個虛擬的擴散項，這可以通過在標準k有粒子的k和值造成影響，還能顯著地提高程序的效率。DkPD(CPC

（4-

（4-為了獲取粒子系統(tǒng)中的細節(jié)運動以及基礎(chǔ)流動的速度，將會在粒子上直接接下來，首先將介紹能量模型的頻譜，然后再描述計算各項同性湍流的方法。前面介紹了湍流勢能的空間分布。下面將要說明湍流的頻譜分布。分析Navier-Stokes方程在頻域中的性質(zhì)有很多方法。這些頻譜方法有一些理想的性LogLog（波數(shù)Log(能量4-14-1，根據(jù)能量的流動過程可以把能量譜分成三個部分：描述是十分的。在平均流動應(yīng)力過程中，湍流的產(chǎn)生主要發(fā)生模型所依賴柯爾莫戈（Kolmogorov）定度表示出來?？聽柲昙僭O(shè)進一步解釋道，在慣性子區(qū)間的能量譜可以（4-E(k)C2/3k

（4-數(shù)，是能量的耗散值。它也叫做柯爾莫戈5-3定律或者叫K41理論。2，。01之間的隨機數(shù)疊加成N(x)E(k)eikxi2(k （4-Ei表示，那么Nf(x)EiNi(xi

（4- N Nu

asasN

Ny(r) -頻率頻 噪聲紋-4-2

kaaakAlgorithm4-1計算張力場forxU和S2計算湍流粘度：vC22 u2PT

更新kkt(|P更新 (CPC 通過 和 的限制使k和都在湍流穩(wěn)定的范圍 uU+2(ak)1/2N

ixx更新qqt(u(t)+(1-)u(t)) U13Navier-Stokes得到速度和張力。從第4行到第15行就是本文的，其中第4行到111214行，在把計算得到的值寫入粒子中。OBjectFBO4-5圍8個點進行采樣，再取平均值，然后得到后處理后的值，如圖4-5Algorithm4-2FragmentShader4-5vec2_win_size=vec2(1.0/win_x,foreachiin[-radium,radium]foreachjin[-radium,radium]color+=texture2D(SmokeTex,gl_TexCoord[0].st+ color.xyz/= gl_FragColor=vec4(color.xyz,本文需要在GPU上同時實現(xiàn)基于歐拉方法的流體模擬和基于粒子方法的半拉格朗日平流法而對于拉格朗日湍流模型每個粒子了它自身的位置速度以及湍流特有的系數(shù)k、、k和q。這些變量的計算由（4-5Ak湍流勢能k的正常范圍可以由k1.5*U2I20動力學(xué)領(lǐng)域的，得到湍流密度的最小值為Imin=1.0，而自然的湍流密同樣地，我們也可以定義一個然最低值）來獲得的最小值。根據(jù)（4-3）湍流粘度的定義，可以得

k

可以得到

。其中，Lmin0.1 k和多數(shù)情況下，在低湍流區(qū)域加入湍流，因此可以以k和初始值k和 可以估計最典型的k和值，這類似于前一小節(jié)描述的估計最大值的方法。因此，最小長度尺寸Linlet成為模型中另外一個很重要的參數(shù)，它用來對產(chǎn)生湍流的區(qū)域進行調(diào)節(jié)。還有一種可能是用k和的下限值k和 q

Dq/DtuUqG將會作為局部相關(guān)聯(lián)粒子的紋理坐標的參考系。相對于導(dǎo)向粒子的位置來計算， （4-q=x-x （4-工跳變現(xiàn)象，如圖4-6所示。因此，本文使用（4-9）來解決這個問題。u'=u(t)+(1-)u(t

（4-t1t2的線性方程組合而成。其中，起始階 4-6t1t2本章首先介紹湍流的k度場和湍流細節(jié)速度場最終結(jié)合的方法5-15-1In(R)Core(TM)i3-2100CPU@3.10GHz3.10GHzNVIDIAGeGTSWindow764VisualStudio2008OpenGL3.0，動的粒子系統(tǒng)，運用小波噪聲的生成和k湍流模型，對煙霧的湍流流動做了用。為了方便計算，生成的小波噪聲是在一維數(shù)組中的，如圖5-1（b）所像素有4個顏色通道，每個通道可以一個顏色值）的紋理中，但是小(a)128*128*128圖5-1兩種噪聲的比較小波噪聲數(shù)據(jù)保存在圖5-1（b）所示的扁平化的小波噪聲紋理中，第二次運行程序是就直接從紋理 小波噪聲就可以了，這可以提高程序的運行效率對于本文所追求的低分辨率的網(wǎng)格粒子方法，128*128*128的小波噪聲完全可5-1（a）所示，柏林噪聲無論在形態(tài)上還是在實際運用中相比于小波噪聲都要遜色不少。從柏林噪聲中可以看到明(a)第24 (b)第46(c)第68 (d)第1005-25-2所示，其中，分辨X=10Y=25Z=10D=44以每幀在有煙霧密度的網(wǎng)格隨機生成1-100粒子數(shù)目，使用后處理方法所生成(a)第24 (b)第48(c)第68 (d)第1005-3始做湍流運動，按k模型的湍流統(tǒng)計學(xué)性質(zhì)，使粒子按湍流軌5-45-532*32*32的網(wǎng)格分辨率下，只保留湍流粒子所得到的細（a）（d 5-632*32*32的網(wǎng)格下，對三維煙霧從六個方向(a)左視 (b)右視(c)前視 (d)后視(e)俯視 (f)仰視5-61005-75-232*32*32這樣比較Navier-Stokes方程所生成的煙霧流動效果。把兩者對32*32*32這樣比較粗糙的分辨率下，看不見煙霧 5-732*32*32N-S 5-864*64*64N-S5-85-264*64*64這樣分辨率比較Navier-Stokes方程所生成的煙霧流動效果。在本文所使用的電腦1幀/秒，雖然煙霧的細節(jié)有所體現(xiàn)，但是還是能夠看見煙5-95-25-10以本方法為基礎(chǔ)，分別從數(shù)據(jù)抽象方向和折線圖直觀方向?qū)Σ煌W(wǎng)格分辨 5-92328653435-105-11程序的幀率開始慢慢降。最后由于生成和的粒子達到平衡狀態(tài)，所以7.335幀/42幀/5幀//秒，42幀/5幀/10585520k湍3點ray-castingray-marching序效率，但是仍然一些問題。雖然本文使用OpenGL在CPU和GPU下進GPU并行計算來提高CPUOpenMPGPU的CUDA并行技術(shù)。否則，本文可以加入數(shù)量級更大的粒子使細節(jié)更加完善，效106-107級別，從而使湍流結(jié)果更加完善和出色。由于本方法是基于的統(tǒng)計學(xué)的k湍流模型，為了提高渲染效率，并本文主要對以Navier-Stokes方程為基礎(chǔ)的煙霧湍流效果進行了模擬。Navier-Stokes方程進行了研究和分析，k湍帶入k湍流模型中，求得的湍流勢能k和湍流結(jié)構(gòu)的耗散值；再次，把卷Navier-Stokes方程的煙霧湍流效果模擬的優(yōu)點在于：用的是kNavier-Stokes的生成放置于GPU中，使用GPU加速結(jié)構(gòu)，例如：CUDA來大大地提ray-castingray-marching等，來增強煙霧本身的視覺AndrewSelle,RonaldFedkiw,ByungMoonKim,YingjieLiu,andJarekRossignac.AnunconditionallystableMacCormackmethod.JournalofScienti?cComputing,2008.D.C.Wilcox.TurbulencemodellingforCFD.DCWIndustries,StephenB.Pope.TurbulentFlows.CambridgeUniversityPress,RobertBridson.FluidSimulationforComputerGraphics.AKPeters,J.T.OdenandL.C.Wellford.ysisofviscous?owbythe?niteelementmethod.AIAAJ.,10:1590,1972.C.Hsu.Acurviliniear-coordinatemethodformomentum,heatandmasstransferinsofirregulargeometry.PhDthesis,UniversityofMinnesota,1981.F.HarlowandE.Welch.Numericalcalculationoftime-dependentviscouspressible?owof?uidswith surface.PhysicsofFluids,8,1966.R.A.GingoldandJ.J.Monaghan.Smoothedparticlehydrodynamics:theoryandapplicationtonon-sphericalstars.Mon.Not.R.Astron.Soc.,pages375–389,M.Müller,B.Solenthaler,R.Keiser,andM.Gross.Particle-based?uid-?uidinction.ACMSIGGRAPH/EGSymposiumonComputerAnimation,2005.[10]J.Hardy,O.DePazzis,andJ.Pomeau.Moleculardynamicsofaclassicallatticegas:Transportpropertiesandtimecorrelationfunctions.PhysicalReviewA,13:1949–1960,1976.[11]N.Thuerey,K.Iglberger,andU.Rüde.SurfaceFlowswithMovingandDeformingObjectsforLBM.ProceedingsofVision,ModelingandVisualization2006,pages193–200,Nov2006.[12]L.Prandtl.¨ubereinneuesformelsystemf¨urdieausgebildeteNachr.Akad.Wiss.G¨ottingenK1,pages6–10,[13]J.Smagorinsky.Generalcirculationexperimentswiththeprimitiveequations.i.thebasicexperiment.MonthlyWeatherReview,1963.[14]B.S.BaldwinandH.Lomax.ThinLayerApproximationandAlgebraicModelforSeperatedTurbulentFlows.AmericanInstituteofAeronauticsandAstronauticsJournal,1978.[15]P.R.SpalartandS.R.Allmaras.Aone-equationturbulencemodelforaerodynamic?ows.AIAAPaper,92:0439,1994.[16]B.E.LaunderandD.B.Sharma.Applicationsoftheenergy-dissipationmodelofturbulencetothecalculationof?ownearaspinningdisc.Lett.HeatMassTransf.,1:1031–138,1974.[17]D.C.Wilcox.TurbulencemodellingforCFD.DCWIndustries,[18]M.JonesandM.Chen.Anewapproachtotheconstructionofsurfacescontourdata.ComputerGraphicsForum,13(3):pp.75–84,[19]J.Smagorinsky.Generalcirculationexperimentswiththeprimitiveequations.i.thebasicexperiment.MonthlyWeatherReview,1963.[20]D.C.HaworthandK.Jansen.Large-eddysimulationonunstructureddeformingmeshes:towardsreciprocatingicengines.ComputersandFluids,29:493–524,[21]P.R.Spalart.Detachededdysimulation.AnnualReviewofFluidMechanics,41:181–202,2009.[22]B.E.Launder,G.J.Reece,andW.Rodi.Progressinthedevelopmentofareynolds-stressturbulenceclosure.J.FluidMech.,68:537–566,1975.[23]M.K.ChungandS.K.Kim.Anonlinearreturn-to-isotropymodelwithturbulent?uctuations.Phys.Fluids,7:1425–1436,1995.[24]S.B.Pope.Alagrangiantwo-timeprobabilitydensityfunctionequationforinhomogeneousturbulent?ows.Phys.Fluids,26:3448–3450,1983.[25]P.A.Durbin.Areynoldsstressmodelfornear-wallturbulence.J.FluidMech.,249:465–498,1993.[26]B.Aupoix.Modelingofcompressibiliyeffectsinmixinglayers.JournalofTurbulence,5,2004.[27]R.B.Langtry,F.R.Menter,S.R.Likki,andY.B.Suzen.Acorrelation-basedtransitionmodelusinglocalvariables.J.Turbomach.,128:123–143,2006.[28]J.Dans,E.Garnier,andP.Sagaut.Numericalsimulationofactiveseparationcontrolbyasyntheticjet.J.FluidMech.,574:25–58,2007.[29]B.ChaoatandR.Schies.Fromsingle-scaleturbulencemodelstomultiple-scaleandsubgrid-scalemodelsbyfouriertransform.Theor.andComp.FluidDyn.,21(3):201–229,2007.[30]F.MenterandM.Kuntz.Ascale-adaptivesimulationmodelusingtwo-equationmodels.AIAApaper05-1095,2005.JosStamandEugeneFiume.Turbulentwind?eldsforgaseousphenomena.InProceedingsofACMSIGGRAPH,1993.NickRasmussen,DucQuangNguyen,WilliGeiger,andRonaldFedkiw.Smokesimulationforlargescalephenomena.InProceedingsofACMSIGGRAPH,ArnauldLamorletteandNickFoster.Structuralmodelingof?amesforaproductionenvironment.InProceedingsofACMSIGGRAPH,2002.RobertBridson,JimHouriham,andMarcusNordenstam.Curl-noiseforprocedural?uid?ow.ACMSIGGRAPHpapers,26(3):Article46,2007. reKim,NilsThue