亚洲免费av电影一区二区三区,日韩爱爱视频,51精品视频一区二区三区,91视频爱爱,日韩欧美在线播放视频,中文字幕少妇AV,亚洲电影中文字幕,久久久久亚洲av成人网址,久久综合视频网站,国产在线不卡免费播放

        ?

        計算流體力學(xué)的時空觀:模型的時空關(guān)聯(lián)性及算法的時空耦合性

        2021-06-23 14:50:42李杰權(quán)
        空氣動力學(xué)學(xué)報 2021年1期
        關(guān)鍵詞:黎曼通量時空

        李杰權(quán)

        (1. 北京應(yīng)用物理與計算數(shù)學(xué)研究所 計算物理重點實驗室, 北京 100088; 2. 北京大學(xué) 應(yīng)用物理與技術(shù)中心, 北京 100871)

        0 引 言

        1959年,華羅庚先生在《大哉數(shù)學(xué)之為用》(數(shù)與量)[1]一文中關(guān)于“時空”有以下論述:四維空間聽來好象有些神秘,其實早已有之,即以“宇宙”二字來說,“往古來今謂之宙,四方上下謂之宇”(《淮南子·齊俗訓(xùn)》),就是宇是東西、南北、上下三維擴展的空間,而宙是一維的時間。

        接著他引用相對論和生活的常識敘述了時空既獨立又統(tǒng)一的性質(zhì):愛因斯坦不再把“宇”、“宙”分開來看,也就是時間也在進行著。每一瞬間三維空間中的物質(zhì)都占有它一定的位置。他根據(jù)麥克斯韋——洛倫茲的光速不變假定,并繼承牛頓的相對性原理,提出了狹義相對論。狹義相對論中的洛倫茲變換把時空聯(lián)系在一起,但并不消滅時空特點。如向東走三里再向西走三里就回到原處,但時間則不然,共用了走六里的時間。時間是一去不復(fù)返地流逝著。

        華先生意在說明數(shù)學(xué)在描述“宇宙之大”中的作用。本文的引用,旨在強調(diào)時空之間的關(guān)系深刻地植根于計算流體力學(xué)的模型理解和算法構(gòu)造。離開了時空演化,一切歸于“穩(wěn)態(tài)”。正因為時空演化,流體的世界才色彩斑瀾、波瀾壯闊。

        定性甚至定量描述流體世界的時空演化并不容易,數(shù)值模擬同樣是貌似簡單,實則不易。下面舉一個眾所周知的例子:

        ?tu(x,t)+a?xu(x,t)=0

        (1)

        其中a是一個常數(shù),x表示空間坐標,t是時間變量。此模型表示了變量u(x,t)的空間變化(擾動、波)與時間變化(傳播)的關(guān)系,描述了一個時空關(guān)聯(lián)的性質(zhì)。

        u(x,t+Δt)=u(x-aΔt,t)

        (2)

        式(2)事實上比式(1)更加根本和直接,并且不牽涉函數(shù)u(x,t)本身更多性質(zhì),比如正則性。即使u(x,t)是一個方波或脈沖, 式(2)仍然有效。用平移算子e-aΔt?x來表示式(2),則有:

        u(x,t+Δt)=u(x-aΔt,t)=e-aΔt?xu(x,t)

        (3)

        從這里出發(fā),如果想精確表達式(1)的輸運性質(zhì),則需要進行展開(如泰勒展開):

        u(x,t+Δt)=u(x-aΔt,t)=e-aΔt?xu(x,t)

        (4)

        為了設(shè)計實際工程中的計算方法(簡稱算法),需進行必要操作:一是截斷處理,涉及算法與控制方程的相容性,以及精度的概念;二是對空間變化?xu(x,t)的處理(如有限差分),不同的處理得到不同的算法。這些都與函數(shù)u(x,t)的正則性和模型(1)的內(nèi)在性質(zhì)息息相關(guān):正則性越好,精度越高; 模型(1)的時空關(guān)聯(lián)性需要通過式(2)精確描述。這些處理充滿技巧和無限可能,怎樣“令人信服”地設(shè)計時空耦合算法且精準描述相應(yīng)模型的內(nèi)在時空關(guān)聯(lián)性,是一個需要深入思考的本質(zhì)問題。

        近年來,計算流體力學(xué)蓬勃發(fā)展,各種算法和相應(yīng)的應(yīng)用軟件目不暇接。有些算法具有堅實的理論基礎(chǔ),如基于變分原理的有限元方法等。而流體力學(xué)中關(guān)于可壓縮流動的研究,由于流場結(jié)構(gòu)非常復(fù)雜(激波、物質(zhì)界面、渦團和湍流等),人們對流動(解的)性質(zhì)缺乏很好理解,模型的適定性和相關(guān)數(shù)值模擬中算法的探討相當工程化,形式化的表達往往導(dǎo)致似是而非的結(jié)果,常常是模型具有很好的性質(zhì),而相應(yīng)數(shù)值模擬的離散模型失去了該保持的性質(zhì),比如模型的時空關(guān)聯(lián)性和算法的時空耦合性等。不同的出發(fā)點導(dǎo)出的數(shù)值算法往往截然不同。本文試圖在這方面做一點嘗試: 基于連續(xù)介質(zhì)假定探討流體力學(xué)模型時空變化的內(nèi)在關(guān)聯(lián)性,闡述相應(yīng)算法應(yīng)保持的時空耦合性質(zhì)。

        從流體力學(xué)微團法的直接建模開始思考模型的時空關(guān)聯(lián)性。通過對通量的研究,發(fā)現(xiàn)經(jīng)典的瞬時Cauchy通量無法反映時空關(guān)聯(lián)性,提出了某一時間段內(nèi)時間區(qū)間通量,在任意特定時間段內(nèi)反映流體的動力學(xué)過程[2-3],進而提出有限體積格式的基本原理:在連續(xù)介質(zhì)假定性下,任一時間段內(nèi)流體通過控制體邊界的時間區(qū)間通量關(guān)于邊界的擾動是Lipschitz連續(xù)的。這一Lipschitz連續(xù)性恰好體現(xiàn)流場時空關(guān)聯(lián)性的本質(zhì)特征,實現(xiàn)了數(shù)學(xué)上流體力學(xué)方程組弱解和物理上積分型平衡律之間的統(tǒng)一,后者其實就是全離散有限體積格式的物理起源。通過對流體力學(xué)方程組時空關(guān)聯(lián)性質(zhì)的研究,實現(xiàn)有限體積框架下的時空耦合。在算法的實現(xiàn)過程中,物理量(質(zhì)量、動量、能量等)以及它的變化率這一對量起著重要的作用。一個量的變化率通過它的空間變化關(guān)系來反映。變化率不是一類常規(guī)的數(shù)學(xué)演算,而是反映了內(nèi)在的物理規(guī)律。例如,加速度是速度的變化率,它等于控制體表面受到的力,即應(yīng)力,這是牛頓第二定律。在流體力學(xué)計算方法中使用這樣的量是自然的也是必須的,這反映了時空變化的聯(lián)系。在文獻[4]中,筆者倡導(dǎo)的時空耦合高精度算法正是這一思想的具體表現(xiàn)。本文再次描述實現(xiàn)時空耦合算法的一種基本流程,并通過算例展示時空耦合算法的數(shù)值表現(xiàn)。

        從偏微分方程理論上的Cauchy-Kowalevski方法到數(shù)值解中的Lax-Wendroff方法,都是某種意義下的時空耦合方法?,F(xiàn)代計算流體力學(xué)的算法時空觀常有體現(xiàn),如迎風(fēng)格式、廣義黎曼解法器[5-7]、守恒元/解元(CE/SE)方法[8]、氣體動理學(xué)解法器[9-10]等。最近幾年,筆者致力于思考高精度時空耦合算法的必要性與理論基礎(chǔ),但遠遠不夠深入。隨著計算技術(shù)的提高和工程需求的精細化,這方面的研究應(yīng)該得到重視。應(yīng)當建立相應(yīng)的時空觀,這對算法的設(shè)計與實現(xiàn)很有幫助。故作此文,拋磚引玉。

        1 流體力學(xué)方程組的時空關(guān)聯(lián)性

        眾所周知,流體力學(xué)建模過程常用微團法,即在某一時間段(t,t+δt),研究控制體Ω(t)=(x,x+δx)內(nèi)流體質(zhì)團的運動。時間間隔δt非常重要,給出了微團運動的動力學(xué)過程響應(yīng)時間。為了描述方便,基于連續(xù)介質(zhì)假定,流體力學(xué)方程組寫成如下形式:

        給定一個控制體Ω(t),流體運動滿足下列的積分平衡律(略去外力場等效應(yīng)):

        這里u是密度函數(shù),n是?Ω(t)的單位外法向。相應(yīng)地,稱:

        為控制體Ω(t)上的質(zhì)量,右邊:

        為Cauchy通量,f為通量密度函數(shù),簡稱為通量函數(shù),它是u以及空間變化量?u等的函數(shù)f=f(u,?u,…)。方程(5)反映了質(zhì)量時間變化率與通過控制體邊界?Ω(t)的通量之間的瞬時關(guān)系。

        Cauchy在特定連續(xù)介質(zhì)假定下,論述了Cauchy通量的性質(zhì)[11-12],這是他對連續(xù)介質(zhì)力學(xué)最重要的貢獻[11]。要深入理解這一關(guān)系并不容易,這與Gauss-Green定理密切相關(guān)。基于這一定理,可以將式(5)寫成散度型偏微分方程(PDE)形式:

        ?tu+?·F(u,?u,…)=0.

        (8)

        這里,略去雷諾輸運定理等細致討論,并設(shè)Ω不隨流場運動,即歐拉型控制體。相應(yīng)地,式(1)中的控制體為拉格朗日型控制體。

        方程組(8)對光滑流動是有效的,可由偏微分方程的知識建立時空的關(guān)聯(lián)性。注意這里的通量函數(shù)F(u,?u,…)與式(7)略有不同,它涉及雷諾輸運定理等有關(guān)變換,這里不做詳細討論。以后歐拉控制體邊界上的通量用大寫字母表示,否則用小寫字母。對于實際流動來說,這一轉(zhuǎn)化是非常粗糙的,微妙之處在于Cauchy通量C(?Ω;t)的數(shù)學(xué)性質(zhì),即正則性??蓧嚎s流場蘊涵豐富的現(xiàn)象,如沖擊波、物質(zhì)界面、湍流等效應(yīng),流場(或稱方程的解)的正則性非常弱,人們對它的認識很少,故而該定理的應(yīng)用并不顯然。歷史上有很多研究[13-16],直到最近,這方面的研究仍是如火如荼[17],遺憾的是所取得的結(jié)果與實際的期望仍然相差甚遠。正如最近的專著[18]的1.3節(jié)有這樣一段論述,“Regarding the right-hand side of (1.3) one needs to keep in mind the following comment concerning the identification of the boundary flux: the drawback of this functional analytic, demonstration is that it does not provide any clues on how the qDmay be computed from A” ,其中qD對應(yīng)于這里的f(u)·n,A對應(yīng)于一個給定的應(yīng)力張量。也就是說:人們還不知道如何從應(yīng)力張量中計算出通量密度函數(shù)。

        仔細審視可以看出Cauchy通量C(?Ω;t)是表示一個特定時刻t流動的瞬時行為,方程(5)只是表示一個非常弱的“時-空”關(guān)聯(lián)性質(zhì)?,F(xiàn)在常用的一個觀點是下面的“弱解”概念,它來自于偏微分方程(8),利用了分布理論。

        定義1.1設(shè)Ω是n中的一個有界區(qū)域,n=1,2,3, 0≤t1≤t2。 對任意試驗函數(shù)如果函數(shù)u(x,t)滿足:

        (9)

        則稱u(x,t)是方程(8)的弱解。

        這個表達式比式(5)前進了一大步:如果解u(x,t)是光滑的,式(8)和式(9)是等價的。近年來偏微分方程理論得到極大的發(fā)展,人們對于式(8)的認識比直接對式(5)的認識要深刻得多,且偏微分方程(8)的時空關(guān)聯(lián)性一目了然:流場的空間攝動可以通過式(8)反映到時間的變化中, 線性波的傳播式(1)正是這種時空關(guān)聯(lián)性的特例。另外,當一個間斷面在無限薄的假定下,可以導(dǎo)出Rankine-Hugoniot間斷關(guān)系,即間斷面(線)兩側(cè)解的“跡”的關(guān)系。

        定義1.1并沒有對函數(shù)u(x,t)做過多假設(shè),只要式(9)成立即可。一旦流場出現(xiàn)復(fù)雜間斷或其它流場結(jié)構(gòu)時,這個定義并不能給出太多的信息。我們回避不了的事實是:在可壓縮流動中,流場即使初始性質(zhì)非?!昂谩?,但在不久的將來由于復(fù)雜的非線性相互作用,可能變得非?!霸愀狻?。這可以從最簡單的Burgers方程看到[19-20]。因此,人們不得不面對復(fù)雜的情形,深入思考為什么應(yīng)用Gauss-Green定理理解Cauchy通量非常困難。

        任何流動都有一個動力學(xué)過程,這也許是個常識, 但仍需要嚴格論證。最近的研究表明[3-4],流動的時空關(guān)系既相互獨立又彼此關(guān)聯(lián)。下面的原理深刻地反映這一事實。

        有限體積方法基本原理。設(shè)u(x,t)是定義1.1下方程(8)的弱解,Ω是n中任一緊集。假定:

        (i)u(x,t)是局部有界;

        (ii) 質(zhì)量m(t)是時間t的連續(xù)函數(shù)。

        那么有下列結(jié)論:

        (10)

        (ii) 對任意δ>0,時間段(t,t+δt)內(nèi)流過Ω邊界?Ω的時間區(qū)間通量:

        關(guān)于Ω的邊界?Ω擾動是Lipschitz連續(xù)的。

        由此,我們給出方程(8)的弱解另一種表述,它恰恰是物理上常見的積分型平衡律。

        定義1.2設(shè)Ω是n中任一緊集,δt>0。如果它滿足下列條件:

        (i)u(x,t)局部有界,且質(zhì)量m(t)是時間t的連續(xù)函數(shù);

        (ii) 整體通量(4)有意義且關(guān)于Ω的邊界?Ω擾動是連續(xù)的;

        (iii) 積分平衡律成立:

        其中n是?Ω的單位外法向, 則稱函數(shù)u(x,t)是偏微分方程(8)的弱解。

        事實上,已經(jīng)證明定義1.1和定義1.2是等價的[3-4],并且給出了通量(11)的正則性質(zhì)。這個定義說明后面討論的有限體積方法就是計算流體力學(xué)的基本格式,和物理最原始的建模一致。流體力學(xué)方程組(8)可理解為:如果流場光滑,用偏微分方程組(8)刻畫;但當流場失去正則性時,解需要滿足積分平衡律式(12)。流體力學(xué)有限體積格式的構(gòu)造過程就是基于式(12)的直接建模過程:在偏微分方程(8)誘導(dǎo)出的時空關(guān)聯(lián)解輔助下,構(gòu)造出和式(12)相容的離散模型,精準反映真實的物理過程。

        2 有限體積方法及其時空耦合性

        2.1 半離散有限體積方法

        在Ωj上,對式(8)的空間散度項應(yīng)用Gauss-Green公式,可得半離散有限體積方程:

        它直接對應(yīng)于流體力學(xué)方程組(5)。在初始時刻t=0, 給定u(x,t)的分布:

        將之帶入(13)即得:

        (19)

        半離散有限體積方法屬于線方法,從某種意義上來說是時空解耦方法。在光滑流場中,流體力學(xué)方程組的偏微分方程關(guān)系可以直接反映時空關(guān)聯(lián)性。但對間斷解來說,式(19)中的局部誤差的累積會給數(shù)值模擬帶來巨大傷害。

        2.2 全離散有限體積方法

        本文將把式(12)應(yīng)用到時空控制體Ωj×[tn,tn+1],tn+1=tn+Δtn,Δtn>0為時間步長,得到有限體積格式的全離散形式:

        也可以把式(20)看作是式(13)在時間段(tn,tn+1)上的積分。與半離散情形下求解瞬時通量(16)不同,這里需要利用偏微分方程(8)的時空關(guān)聯(lián)性質(zhì)逼近[tn,tn+1]上時間區(qū)間通量:

        (21)

        將之代入(20),得到有限體積公式:

        這里記逼近解在控制體Ωj上的平均值為:

        由有限體積方法基本原理,可以得到:

        =O(Δtn)q+1|Γjl|

        (24)

        其中q>0。值得注意的是,式(24)中的誤差和式(19)中的誤差非常不同,式(19)中的誤差是用解的局部跳躍(變差)來測量,而式(24)中的誤差直接用時間步長(等價于網(wǎng)格尺度)來測量。當流場相對光滑時,這兩者是等價的,但當流場中有劇烈間斷、脈動時,兩者差別是明顯的。即使網(wǎng)格加密,式(19)也得不到收斂解(在標量方程的研究中,網(wǎng)格加密是可以得到收斂解的,原因在于全局變差有界性質(zhì)。到目前為止還沒有例證可證明該性質(zhì)在流體力學(xué)方程組成立[7])。后文算例4.1可以看出數(shù)值通量的重要性。

        2.3 有限體積方法的相容性

        所謂相容性,描述的是離散格式與背景方程之間的關(guān)系。傳統(tǒng)上常常使用Taylor展開的方式研究數(shù)值格式與偏微分方程(如式(8))之間的相容性,這樣的運算只對光滑流場成立。當流場比較復(fù)雜時,目前大部分的數(shù)值分析某種意義上只是啟發(fā)性的,比如以標量方程為模型建立相關(guān)的理論[21]。

        對半離散格式(18)來說,在每個固定時刻給出的通量誤差至多為:

        =O((Δu)jl)|Γjl|

        =O((Δu)j)|Γj|dj

        (26)

        其中(Δu)j表示在控制體附近解的局部跳躍,|Γj|是Ωj邊的最大長度,dj=diam(Ωj)。并且在每個時間層,隨著Runge-Kutta步的增加,誤差也會累積。在式(26)中,dj來源于不同方向上通量差的獲利,在3.5節(jié)可以進一步看到。

        對于全離散方法(25),我們討論每個時間層數(shù)值通量的逼近,并期望有下列的估計:

        =O(Δtn)q+2|Γj|.

        (27)

        其中q>0。比較式(26)和式(27),它們有本質(zhì)的差別,即(Δu)j與Δtn的差別。對光滑解來說它們是等價的。這樣我們給出下列的定義。

        定義2.1 (有限體積格式的相容性)。設(shè)Ωj是任一控制體,j=1,…,N, 如果對于某個q>0,式(27)成立,則有限體積格式(25)相容于平衡律(12),并具有q階相容性。

        如上所述,相容性概念常常相對與偏微分方程所言。如果式(8)是雙曲守恒律,Lax和Wendroff 給出了有限差分方法的相容性,常稱為Lax相容性[23]。由于Lax相容性概念影響深遠,有必要進行回顧。

        定義2.2 (Lax相容性[23])??紤]雙曲守恒律方程組:

        ?tu+?xf(u)=0

        (28)

        其2p+1點守恒性差分格式為:

        g(u,…,u)=f(u)

        (30)

        則稱差分方法(29)與雙曲守恒律(28)相容。

        基于這個定義,建立了影響深遠的Lax-Wendroff收斂定理,直至今日,該定理仍被廣泛應(yīng)用。隨著高階精度格式的發(fā)展以及工程應(yīng)用的需求,這個定理常常被誤用,典型表現(xiàn)為:

        (ii) 在此基礎(chǔ)上建立的Lax-Wendroff定理只適用于一致網(wǎng)格或結(jié)構(gòu)網(wǎng)格,對非結(jié)構(gòu)網(wǎng)格并不適用[24-25]。數(shù)值驗證過程中使用的網(wǎng)格加密方法測試收斂階,需要倍加小心。

        實際上,定義2.1第一次給出高精度有限體積數(shù)值格式與積分平衡律之間的相容性[3]。有限體積方法與網(wǎng)格的結(jié)構(gòu)無關(guān),因此適用于無結(jié)構(gòu)網(wǎng)格。特別地,收斂階測試時一般針對光滑流進行,數(shù)據(jù)重構(gòu)部分能夠達到應(yīng)有的精度,通量逼近階事實上就是收斂階。但是,當流場含有間斷或其他復(fù)雜結(jié)構(gòu)時,數(shù)據(jù)重構(gòu)仍存在諸多爭議,通量的逼近效果往往被忽略。通過以上分析可以看到,如果通量相容性(27)不成立,逼近解不可能收斂。也就是說,式(27)是有限體積格式相容性的一個必要條件。

        2.4 再訪Godunov方法

        談到流體力學(xué)中的有限體積方法,需要回顧Godunov方法[26]。經(jīng)過半個多世紀的發(fā)展和檢驗,它已成為現(xiàn)代計算流體力學(xué)基石??紤]一維可壓縮歐拉方程組:

        ?tu+?xF(u)=0

        (31)

        即積分平衡律(12)的形式。對應(yīng)于式(22)有限體積格式為:

        ?tu+?xF(u)=0,x∈,t>tn,

        由于

        可知式(35)和式(12)完全一致。因此,從積分平衡律的逼近來說,基于分片常數(shù)的初值逼近,式(35)與式(12)完全相容或通量有無窮階相容性。從整體逼近式(25)來說,所有的誤差來自于投影過程。因此習(xí)慣上稱Godunov格式為一階格式。

        如果用逼近的黎曼解法解,界面上的值一般可以寫為:

        這個誤差在強間斷的情形下對計算結(jié)果有巨大的傷害。如式(19)所述,這個誤差并不隨著網(wǎng)格加密而減小。

        在多維情形下,例如二維雙曲守恒律情形:

        ?tu+?xF(u)+?yG(u)=0

        (38)

        (i) 相對于界面的法向黎曼問題

        由于流體力學(xué)方程組的伽利略不變性,通過旋轉(zhuǎn)變換,總可以設(shè)x-方向為Γjl的外法向,Γjl兩側(cè)單元的值記為uL和uR,法向黎曼問題定義為:

        ?tu+?xF(u)=0, (x,y)∈2,t>0

        它的求解和上面的一維黎曼問題(34)沒有本質(zhì)差異。顯然,此解只反映了沿Γjl法向流場的變化情況。

        (ii) 頂點處二維黎曼問題

        這是真正的多維問題,可以寫成如下形式[27]:

        ?tu+?xF(u)+?yG(u)=0,

        (x,y)∈2,t>0

        u(x,y,0)=u, (x,y)∈Θ

        (40)

        這里Θ,=1,…,K,是從原點出發(fā)的角形區(qū)域,見圖1。由于有限傳播性質(zhì),從中心點發(fā)出的復(fù)雜波結(jié)構(gòu)只會影響有限的區(qū)域。僅從通量的構(gòu)造來說,頂點處解對界面通量的貢獻占比很小,可以適當限制Courant數(shù),減少頂點周邊流場對數(shù)值通量的影響,從而可以忽略此問題的求解。在移動網(wǎng)格方法中,特別是拉格朗日方法中,需要用頂點解確定頂點運動速度,于是頂點處二維黎曼問題(40)的解非常重要。但由于解的結(jié)構(gòu)過于復(fù)雜,倒不如通過頂點近似黎曼解法器來處理更為方便可靠[28]。

        圖1 二維黎曼問題的初始數(shù)據(jù)分布

        黎曼問題的(逼近)解被廣泛用到雙曲守恒律,并延展到更一般的流體力學(xué)方程組半離散數(shù)值方法中。與下面廣義黎曼問題的解進行比較可以發(fā)現(xiàn),這個解不能充分反映出流體的動力學(xué)過程;即使從微觀角度,歐拉方程反映了粒子無窮碰撞的結(jié)果。再者,在每個控制單元上及初始時間層,流場處于常狀態(tài),空間的變化依賴相鄰單元之間的間斷來實現(xiàn)。因此,Godunov格式的解不能充分反映瞬時行為。對長時間、多物理以及多尺度現(xiàn)象,數(shù)值模擬結(jié)果常常出現(xiàn)似是而非的現(xiàn)象。

        到目前為止黎曼問題只限于對雙曲守恒律進行研究, 相關(guān)的拓展可以從下面的廣義黎曼問題研究中看出。

        2.5 高階數(shù)值通量與廣義黎曼問題

        一般地,方程組(8)的廣義黎曼問題可以寫成如下形式:

        ?tu+?·F(u,?u,…)=0,x∈n,t>0

        其中Φ(x)=0是定向的n-1維光滑曲面,經(jīng)過適當?shù)淖鴺俗儞Q,可記該曲面為x=0 (記空間坐標為x=(x,y,z)),即式(41)可化為如下形式:

        根據(jù)有限體積方法基本原理,需要直接逼近:

        這里Aε={(0,y,z);y∈(y0-ε,y0+ε),z∈(z0-ε,z0+ε)}是面x=0上的一部分。在一維情況下,(43)即為:

        根據(jù)已有流體力學(xué)方程組相關(guān)的偏微分方程理論與應(yīng)用研究,可以有效地逼近或求解式(42),滿足與式(27)對應(yīng)的相容性要求。具體地,式(44)有:

        這里需要被積函數(shù)F(u(0;t))滿足關(guān)于時間t的正則性。對于式(42)中的初始數(shù)據(jù),這一點可以得到保障,從而可以對式(44)進行相容性逼近。

        定義2.3 (有限體積方法的時空耦合性)。考慮時空關(guān)聯(lián)模型(8)的有限體積方法(25)。如果方程(8)的解可以用來逼近數(shù)值通量,并使式(25)在定義2.1的意義下相容, 數(shù)據(jù)重構(gòu)也充分利用式(8)的時空關(guān)聯(lián)信息,則算法(25)是時空耦合的。

        注意文獻[22]給出的Hermite型數(shù)據(jù)重構(gòu)方法是時空耦合的,該方法已用在兩步四階方法[4,29]中。文獻中數(shù)據(jù)重構(gòu)的時空耦合性研究還不充分,上述定義中“數(shù)據(jù)重構(gòu)也充分利用式(8)的時空關(guān)聯(lián)信息”這句話比較模糊,需要做更深入的研究。

        2.6 幾個時空耦合算法的例子

        下面給出幾個時空耦合算法的例子。

        2.6.1 線性輸運方程

        對于線性方程(1), 其有限體積格式可以寫成:

        從而由(46)可得:

        這就是迎風(fēng)格式。數(shù)值通量的構(gòu)造完全使用了解的時空關(guān)聯(lián)信息,因此迎風(fēng)格式(49)是唯一的一階時空耦合格式。眾所周知,迎風(fēng)格式在所有一階穩(wěn)定格式中具有“最優(yōu)”性質(zhì)。

        如果初始數(shù)據(jù)是分片多項式,特別是分片線性函數(shù)時:

        則式(1)的解式(2)可寫為:

        t∈(tn,tn+1)

        (51)

        直接計算可得:

        以及

        更一般地,我們可以利用文獻[29,22]的方式來構(gòu)造高階時空耦合方法,或者Lax-Wendroff型的單步高階方法,只是推廣到實際工程問題時,真正非線性和多維性會給單步方法帶來實質(zhì)性的困難。

        2.6.2 線性對流擴散方程

        考慮下列方程:

        這里物理黏性系數(shù)μ>0。把它寫成散度形式有:

        ?tu+?x(au-μ?xu)=0,μ>0

        (55)

        在時空控制體上,把它表示為式(32)的形式,進行通量逼近。文獻[8]中的守恒元/解元方法展示其時空耦合的屬性。由于其構(gòu)造需要一點篇幅展示,有興趣的讀者可以查閱文獻[8]及其后來的一系列文章。

        相對來說, Navier-Stokes方程組及其相關(guān)時空關(guān)聯(lián)模型的時空耦合算法研究較少. 盡管形式上可以進行簡單的時空對換處理,但對耗散項等高階空間導(dǎo)數(shù)項的處理仍需要嚴格數(shù)學(xué)論證。當然,GKS方法間接提供了Navier-Stokes方程的數(shù)值算法,具有時空耦合性[9,30]。

        2.6.3 線性松弛系統(tǒng)

        考慮下列簡單的松弛系統(tǒng):

        u(x,0)=u0(x)

        (56)

        其中τ>0是松弛參數(shù),a為常數(shù),v也是常數(shù)。這個方程的有限體積形式為:

        為了構(gòu)造時空耦合算法,可以用式(56)的解表達式:

        (59)

        很顯然,格式(59)~式(61)是時空耦合的,并具有時空二階精度。

        3 基于廣義黎曼解法器 (GRP solver)的時空耦合算法

        數(shù)值通量的構(gòu)造是執(zhí)行有限體積格式的兩個核心步驟之一。對于非線性問題我們一般無法像上述線性方程那樣,得到解的表達式。下面將要利用廣義黎曼問題(41)或(42)解的局部正則性質(zhì),發(fā)展廣義黎曼解法器。其主要思想可以概括為:

        為了方便敘述,這一節(jié)統(tǒng)一把界面放在x=0, 把廣義黎曼問題的求解點平移到坐標原點(0,0), 初始數(shù)據(jù)表示成式(42)的形式。記:

        廣義黎曼問題(GRP)解法器:給定控制方程以及形如式(42)的分片光滑函數(shù)作為初始數(shù)據(jù),求解式(42) 并得到形如式(62)的極限值。

        一旦有了極限值(62),受益于解u(x,t)在界面上關(guān)于時間t的連續(xù)性,我們有:

        u(0,δt)=u*+(?tu)*δt+O(δt2)

        (63)

        由于u*只表示了一種瞬時行為,其解由相應(yīng)的黎曼解給定:

        u*=R(0;uL(0),uR(0))

        (64)

        接下來的關(guān)鍵任務(wù)是求值(?tu)*。“非常不嚴格”地由控制方程(42)直接得到:

        然后求得極限值。這樣解u(x,t)的變化率可以通過空間的變化反映,就是經(jīng)典的Lax-Wendroff方法[23],或Cauchy-Kowalevski方法[31]。通過這種途徑把模型的時空關(guān)聯(lián)性質(zhì)嵌入到數(shù)值格式中,所得的算法具有時空耦合性質(zhì)。

        3.1 線性GRP解法器

        對于線性系統(tǒng):

        ?tu+A?xu=Cu+D,t>0

        其中A、C、D是常數(shù)矩陣,A可以對角化,Λ=R-1AR, 其中Λ是由A的特征值λi組成的矩陣,Λ=diag{λi}。記w=R-1u, 有:

        ?tw+Λ?xw=R-1Cu+R-1D,t>0

        (I+R-1CuL+I-R-1CuR)+R-1D

        (68)

        (RI+R-1CuL+RI-R-1CuR)+D

        (69)

        特別當uL=uR=u*時有:

        從中可以看出如何應(yīng)用Lax-Wendroff解法器和間斷追蹤計算(?tu)*。

        多維情形是類似的。例如考慮二維線性方程組:

        ?tu+A?xu+B?yu=Cu+D,t>0

        這里A、B、C、D是常矩陣,且A、B可對角化。切向變化量B?yu可看作相對法向的一個擾動, 將線性方程組寫成如下形式:

        ?tu+A?xu=-B?yu+Cu+D,t>0

        與上面類似方法,我們可以得到:

        (?tu)*=-[A+(?xu)L+A-(?xu)R]-

        RI+R-1B(?yu)L+RI-R-1B(?yu)R]+

        (RI+R-1CuL+RI-R-1CuR)+D

        (73)

        特別強調(diào),通過GRP解法器把切向變化自然地蘊含到數(shù)值通量中,這是法向(一維)黎曼解法器做不到的[38]。換句話說,如果只使用黎曼問題解法器,導(dǎo)致格式失去多維性,該格式應(yīng)用到多維問題模擬時自然會有數(shù)值缺陷,不得不用其他方式進行補償。再者,GRP解法器是保證數(shù)值格式真正多維性的有效手段,在合理的數(shù)據(jù)投影基礎(chǔ)上,算法具有渦量保持等性質(zhì)。

        3.2 聲波近似 ——線性化GRP解法器

        對于非線性系統(tǒng)來說,當流場是光滑的或者波的強度不大時,使用近似解法器進行通量的逼近就可以取得不錯的效果,即聲波近似或線性化GRP解法器。Toro等的ADER解法器就是GRP的聲波近似版本[35]。

        先看下面的一維雙曲平衡律:

        ?tu+?xF(u)=H(x,u),x∈,t>0

        ?tv+A(u*)?xv=H(x,u*),u=u*+v

        (75)

        這是一個線性系統(tǒng)。類似于式(66),從中可以計算(?tu)=(?tv)*,

        這里A±(u*)=R(u*)Λ±(u*)R-1(u*),Λ±(u*)是由Λ(u*)的“±”部分組成Λ+=diag{max(λi),0)},Λ-=diag{min(λi),0)}。從而

        當uL(0-0)≠uR(0+0)時,聲波近似策略如下:以局部黎曼解u*=R(0;uL(0),uR(0))為背景解, 線性化式(74)得到式(75),從而可得線性化GRP解法器式(77),具體見文獻[32]。

        對于多維系統(tǒng)(仍然以二維為例):

        ?tu+?xF(u)+?yG(u)=H(x,y,u),t>0

        采用聲波近似的策略,線性化這個方程組得到:

        ?tv+A(u*)?xv=-B(u*)?yv+H(x,y,u*),t>0

        從而得到類于式(73)的解法器, 這里B(u*)=?uG(u)。

        3.3 真正非線性GRP解法器

        在式(74)中,如果初始數(shù)據(jù)在原點(相對網(wǎng)格尺寸)有“很大”跳躍,就會出現(xiàn)強非線性波。線性化GRP解法器不足以分辨強非線性波,需要使用真正非線性GRP解法器,否則就會出現(xiàn)明顯的數(shù)值缺陷[7,34-35]。一般的界定標準為:

        ‖uL(0)-uR(0)‖?αΔx

        (80)

        參數(shù)α>0是一個重要的經(jīng)驗參數(shù)。求解GRP解法器的基本思想是:解析強稀疏波,追蹤強間斷。特別需要指出,在強間斷的情形下,熱力學(xué)效應(yīng)起著重要作用,真正非線性的GRP解法器反映了熱力學(xué)相容的性質(zhì)[7]。

        本文不具體給出真正非線性GRP解法器。對于可壓縮歐拉方程組,可參閱文獻[5,36]。后期發(fā)展的不依賴于坐標系統(tǒng)的GRP解法器,詳見文獻[6]。對一般的雙曲方程組,可參見文獻[37]。

        對于多維情形,仍然可以把切向的變化和間斷面的變化看作擾動,考慮有下列形式問題:

        ?tu+?xF(u)=-?yG(u)+H(x,y,u),t>0

        從而利用3.1和3.2節(jié)中法向解法器求解式(72)和式(73),這主要源自于一個關(guān)鍵的事實:切向擾動在法向的傳播是連續(xù)的。由此切向的變化在通量中得到充分反映,得到的數(shù)值方法具有真正的多維性[38]。

        正如前面所討論的那樣,除非涉及(依賴于網(wǎng)格移動)中心拉格朗日型數(shù)值方法,這里不關(guān)注頂點多維黎曼解法器。由于真正多維黎曼問題的理論基礎(chǔ)很不成熟[27,39-41],真正多維頂點黎曼解法器常常帶來不穩(wěn)定的數(shù)值結(jié)果[42]。在實踐中,人們更傾向使用健壯的逼近GRP解法器, 例如,Maire等利用新的框架并結(jié)合拉格朗日型GRP解法器[43],發(fā)展了穩(wěn)定的中心型高階拉氏方法。

        3.4 一般時空關(guān)聯(lián)模型的GRP解法器

        對于一般的時空關(guān)聯(lián)模型,如多相流、多介質(zhì)模型[44]和Navier-Stokes方程等,廣義黎曼問題解法器可以進行類似的構(gòu)造,簡述如下。

        對于多相流、多介質(zhì)模型,廣義黎曼問題解法器可以歸結(jié)為一般的非線性平衡律的范疇,模型的時空關(guān)聯(lián)性質(zhì)在GRP解法器中可以得到充分體現(xiàn)[45-47]。困難在于介質(zhì)組份以及混合引起的數(shù)值振蕩、物質(zhì)分數(shù)為負等非物理解,涉及物理建模本身的缺陷以及有限體積格式的投影(平均)過程。由于熱力學(xué)關(guān)系的高度非線性,投影過程不能反映精確的物理過程, 例如內(nèi)能的平均過程導(dǎo)致物質(zhì)界面附近的壓力振蕩. 因此,在投影步實現(xiàn)時空耦合,充分利用解的信息也許可以提高相關(guān)數(shù)值質(zhì)量。

        對于Navier-Stokes方程組等宏觀模型,基本的想法是類似的。對于對流項(歐拉部分),使用上述標準的廣義黎曼解法器。需要指出的是:初始數(shù)據(jù)式(41)需要進行適當?shù)钠ヅ?。也就是說該初始數(shù)據(jù)至少是二階以上分片多項式。另外可使用的策略為,對于對流占優(yōu)問題,采用松弛策略,把相應(yīng)模型雙曲化[48-49]。這樣的策略是基于能量原理——在對流占優(yōu)情形下,能量集中于初始能量傳播影響區(qū)域內(nèi)(雙曲性質(zhì))。

        從Boltzmann方程直接出發(fā)的動理學(xué)方法是設(shè)計流體力學(xué)數(shù)值方法的一條有效途徑。一般來說,很難寫出Boltzmann方程的解,但在特定的近似假定下,如BGK模型[50],可以類似于式(58)那樣隱式得出解的表達式,并將之用于數(shù)值通量的構(gòu)造,得出的數(shù)值方法如氣體動理學(xué)格式(GKS)[9]、統(tǒng)一氣體動理學(xué)格式UGKS[10]等。特別地,GKS可以用來模擬Navier-Stokes方程描述的流動,可以作為Navier-Stokes的解法器來使用。按照這樣的定義,在通量的構(gòu)造部分,GKS和UGKS解法器顯然是時空耦合的。

        3.5 高精度方法中黎曼解法器和GRP解法器的比較

        在高精度數(shù)值方法中,黎曼解法器和GRP解法器都可用來構(gòu)造數(shù)值通量。前者在高階線方法中常用,如Runge-Kutta型高階方法,后者用在兩步四階方法中。下面仍以一維雙曲守恒律方程組(31)來比較兩種解法器的差別。

        由解關(guān)于時間t的正則性得到:

        所以對于間斷解來說,通量的截斷誤差為:

        誤差階q=0。不過,對于光滑解來說,式(84)中的差賦予了一個“階數(shù)”,

        從而誤差階為q=1。

        與上面討論類似,對于間斷解,GRP解法器有一階精度q=1,而對于光滑解有二階精度q=2。

        3.6 時空耦合數(shù)值邊界條件

        邊界條件的數(shù)值處理是計算流體力學(xué)中的一個核心問題,大部分的數(shù)值處理基于對虛擬區(qū)域的延拓(外插技巧)或特征傳播理論[51]。最近逆Lax-Wendroff方法[52]用來數(shù)值處理邊界條件,這種方法蘊含了時空耦合性。我們認識到,流體力學(xué)方程組的數(shù)值邊界條件等價于單邊廣義黎曼問題的數(shù)值求解[53]。事實上,先前的工作已經(jīng)認識到單邊黎曼問題在數(shù)值邊界條件處理上的重要性[54-55]。這些工作與直接的逆Lax-Wendroff方法相比,有以下特點:

        (ⅰ) 在有限體積框架下,計算區(qū)域的邊界自然為控制體的邊界,所以數(shù)值邊界條件處理等價于邊界上的通量數(shù)值逼近,它歸結(jié)為單邊黎曼解法器的構(gòu)造。

        (ⅱ) 單邊黎曼解法器中的(?tu)*其實反映了邊界條件上的受力情況,即牛頓第二定律的數(shù)學(xué)表現(xiàn),這是時空耦合算法的體現(xiàn)。在文獻[53]中,這一事實是單邊黎曼解法器的基石。

        (ⅲ) 在技術(shù)上,如此處理數(shù)值邊界條件可以盡量少使用虛擬網(wǎng)格。在時空二階格式中無需使用虛擬網(wǎng)格; 相比較其它更高階方法,至多使用一半的虛擬網(wǎng)格。

        (ⅳ) 無需使用更高階導(dǎo)數(shù),避免了人為的復(fù)雜性和虛假物理性質(zhì)等[52]。

        單邊黎曼解法器的使用將極大簡化邊界條件的數(shù)值處理,特別是解決無結(jié)構(gòu)網(wǎng)格下虛擬區(qū)域的插值問題[56]。

        3.7 高階時空耦合算法

        高階精度數(shù)值方法對小尺度問題的數(shù)值模擬非常重要,如湍流等?;谟邢摅w積格式的框架式(25),高階方法需要在數(shù)值通量和數(shù)據(jù)重構(gòu)兩部分具有高階逼近性質(zhì)。在通量的逼近部分,可以采用Lax-Wendroff方法,但流體力學(xué)方程組(8)的非線性性質(zhì)導(dǎo)致高階展開過于復(fù)雜,缺乏實際工程價值。更糟糕的是,解的間斷性質(zhì)意味著高階展開的運算沒有意義, 因此需要尋求其他方式。在過去的幾十年中,各種空間重構(gòu)技術(shù)以及Runge-Kutta型時間推進方法取得了巨大的進步,可以容易地查閱到相關(guān)文獻。

        最近,文獻[29]發(fā)展了兩步四階時間推進方法,成功地融合了Runge-Kutta型方法簡單性優(yōu)點和基于Lar-Wendroff型解法器的時空耦合性質(zhì)。文獻[4]詳述了該類方法的特性: (i) 計算效率。同等條件下其計算效率是Runge-Kutta方法的一半(二維情形為例)[57]。(ii) 緊致性。由于時間推進步驟減少一半,模版就會減少一半;時空耦合的HWENO型重構(gòu)[22, 58]可以使計算格式更加緊致。(iii) 穩(wěn)定性。已經(jīng)證明其穩(wěn)定區(qū)域比Runge-Kutta大[59]。(iv) 兼容性。該方法很容易和其它方法相結(jié)合,如GKS解法器[60]、多矩方法[61]以及CRP重構(gòu)技術(shù)[62]等。

        目前這類方法還局限在“顯式”框架內(nèi),相關(guān)的研究處于起步階段。為了應(yīng)用的需要,亟須發(fā)展“隱式”兩步四階方法,以適應(yīng)“強剛性”等多尺度問題的求解。

        3.8 投影過程時空耦合性的一點注釋

        從而有:

        在最近的兩步四階方法中,我們也使用了這一策略[22],得到的數(shù)據(jù)更加緊致,從而格式也具有緊致性。

        4 數(shù)值算例

        在可壓縮流體的算法設(shè)計及數(shù)值模擬中,有大量基準算例。隨著算法進步和工程需求的提高,基準算例應(yīng)該與時俱進,以面對相當極端的環(huán)境。在文獻[65]中,一系列基準算例值得用新的數(shù)值算法進行測試。下面幾個算例,可以看出時空耦合性的重要性。

        4.1 大壓力比(大密度比)問題

        這個問題又稱為LeBlanc問題,它是經(jīng)典的激波管問題的極端情形,從中可以看出數(shù)值方法對強稀疏波的分辨能力以及它對激波產(chǎn)生的影響??刂品匠虨橐痪S歐拉方程(31),初始數(shù)據(jù)具有如下形式:

        (90)

        多方指數(shù)取1.4,計算的終止時間t=0.12。文獻[66]利用這個例子檢測了多個數(shù)值格式的性能。這里我們同樣用這個算例針對性討論本文涉及的一些觀點,數(shù)值結(jié)果來源于文獻[7,29]。

        首先在圖2中展現(xiàn)使用不同解法器的一階格式??梢钥闯鲆浑AGodunov方法隨著網(wǎng)格加密,可以收斂到精確解;使用逼近解法器,收斂相對較慢,但仍然具有收斂性。圖3使用了空間二階重構(gòu),而時間離散使用一階向前歐拉方法,解法器分別使用黎曼解法器和逼近的黎曼解法器(HLLC,Roe), 數(shù)值表現(xiàn)很差,不具有收斂性。正如在3.5節(jié)所述,當空間具有高階精度,使用一階黎曼解法器構(gòu)造數(shù)值通量等,算法不具有相容性;類似的情形反映在圖4, 即使時間離散使用二階Runge-Kutta方法。

        (a) 200網(wǎng)格點

        (a) 1000網(wǎng)格點

        (a) 1000網(wǎng)格點

        圖5考察了在這種極端情況下使用線性化方法的數(shù)值表現(xiàn),黎曼解使用聲波近似解法器。其實文獻[34-35]已經(jīng)指出線性化解法器的數(shù)值缺陷。圖6使用了非線性GRP解法器,數(shù)值結(jié)果立即改善,說明了強非線性出現(xiàn)后真正非線性GRP解法器的重要性。圖7展示的結(jié)果是基于GRP解法器的兩步四階方法[29],從中看到用100網(wǎng)格點可以很好分辨間斷,300網(wǎng)格點可以得到完美效果。這是許多數(shù)值方法很難達到的,從而說明了時空耦合算法的必要性。

        (a) 1000網(wǎng)格點

        (a) 1000網(wǎng)格點

        (a) m=50

        4.2 激波和熵波相互作用的問題

        這個算例是Shu-Osher算例的推廣,考慮了更極端的情形,又稱為Titarev-Toro問題[67]??刂品匠桃廊粸闅W拉方程,初始數(shù)據(jù)為:

        這里使用了GKS解法器[9]得到圖8的數(shù)值結(jié)果,詳細說明見文獻[60]。

        圖8 Titarev-Toro問題. 這里使用1000網(wǎng)格點數(shù),空間重構(gòu)采用了不同的重構(gòu)策略,計算終止時間t=5

        4.3 激波和懸浮剛體的相互作用

        圖9 激波和懸浮剛體的相互作用后壓力的分布情況.

        4.4 多介質(zhì)激波和氣泡的相互作用

        這個例子展示了激波和氣泡相互作用的情況(圖10),該問題已經(jīng)成為多相流領(lǐng)域一個標準算例[68]。下面的結(jié)果取自文獻[45],分別用能量分裂的Godunov方法和GRP方法模擬所得結(jié)果(圖11)。顯然,GRP很好地提高了流場的分辨率。

        圖10 激波撞擊氦氣泡的裝置示意圖. 初始時刻激波的馬赫數(shù)Ms=1.22,計算網(wǎng)格為2500×890, 其他參數(shù)詳見文獻[45]的算例D

        圖11 激波撞擊氣泡不同時刻的陰影圖,每個子圖的左邊是Godunov格式的結(jié)果,右邊是GRP的結(jié)果. (a) t=32 μs; (b) t=62 μs; (c) t=72 μs; (d) t=102 μs; (e) t=427 μs; (f) t=674 μs

        4.5 各向同性可壓縮湍流的模擬

        圖12 各向同性可壓縮湍流的模擬[69]. 其中湍流馬赫數(shù)Mat=1.2, 初始泰勒微尺度雷諾數(shù)Reλ=72, 速度梯度張量第二不變量的等值面Q=25.

        5 討論與展望

        流體力學(xué)的時空關(guān)聯(lián)模型刻畫了物理量空間變化在時間方向上的演化的傳播,如各種波的傳播等。當進行數(shù)值模擬時,這種性質(zhì)理應(yīng)得到繼承,相應(yīng)地所用算法應(yīng)該具備時空耦合特性。實踐過程中這個看似簡單課題理應(yīng)得到關(guān)注。遺憾的是,相對于時空解耦方法,時空耦合算法需要更好理解模型的物理性質(zhì)或數(shù)學(xué)理論,人們自然“避難就易”。一個直白的問題是——即使物理問題的數(shù)學(xué)模型十分完美,當相應(yīng)的算法缺乏相應(yīng)完美性質(zhì)時,其數(shù)值結(jié)果怎么令人信服呢?

        本文首先嚴格論述有限體積方法。正如引言所解釋的原因:(i) 無論從物理定律的形成還是數(shù)值算法的構(gòu)造,有限體積方法是解流體力學(xué)問題最自然的框架,本文總結(jié)了這個框架形成的數(shù)學(xué)理論和內(nèi)在涵義;(ii) 流體現(xiàn)象的復(fù)雜性(如間斷等)客觀地阻礙了有限體積方法嚴格數(shù)學(xué)理論的形成,對于這一框架的認識常出現(xiàn)諸多似是而非的爭論; (iii) 對于許多實際工程問題,基于無結(jié)構(gòu)網(wǎng)格的算法設(shè)計是一個自然要求,在此之上許多方法相互沖突[24],有必要從最基本的地方出發(fā)建立一個根本準則。幸運地是,許多工程軟件,如Fluent, 自然選擇在有限體積框架下生成;許多工程人員當遇到難以跨越的困難時,往往借用有限體積框架度過難關(guān)。從論證過程可以看到時空耦合是算法的根本屬性。

        有限體積格式的步驟很簡單:投影過程和數(shù)值通量的構(gòu)造。目前CFD算法的大部分研究者將注意力集中于前者,基本上在空間逼近論的范疇內(nèi)進行探索。盡管黎曼不變量等物理量以及其他技術(shù)用于重構(gòu)過程,但是在隨機選取方法中重要的時空耦合性質(zhì)相當缺乏。 黎曼問題的解被用在隨機選取中,在某種意義下,它可以看成是一種時空耦合投影方法。本文側(cè)重于后者的討論,給出了有限體積格式與積分平衡律(12)之間的相容性準則。特別對于數(shù)值通量介紹了黎曼解法器和GRP解法器,并在3.5節(jié)中進行了對比。簡單地總結(jié)如下,具體內(nèi)容可以到文獻[2]中查閱。

        (i)關(guān)于Riemann解法器。對于雙曲守恒律(歐拉方程),當初始數(shù)據(jù)是分片常數(shù)時,由Riemann解法器產(chǎn)生的數(shù)值通量具有無窮階相容性,即Godunov格式就是積分平衡律;當初始數(shù)據(jù)是非常數(shù)的分片光滑函數(shù)時,Riemann解法器產(chǎn)生的數(shù)值通量對光滑解是一階相容的,但對間斷解是不相容的(見3.5節(jié))。也就是說對于高階空間重構(gòu),如果不能有匹配的數(shù)值通量,產(chǎn)生的數(shù)值格式是不相容的。值得注意的是,如果控制方程包含外力項時,Godunov格式永遠是不相容的。

        (ii)關(guān)于GRP解法器。GRP解法器給出時空耦合的數(shù)值通量。對于光滑流場,GRP解法器得到二階相容的數(shù)值通量;但對于間斷解只有一階時間精度。GRP解法器是保證有限體積格式的逼近解收斂的一個必要條件。

        這里之所以再次強調(diào)這點并細致給出第算例4.1節(jié)的, 原因在于為了看清它與傳統(tǒng)理解上的差異。事實上,算法時空耦合性也體現(xiàn)在其他方面,比如最近研究氣體動理學(xué)的漸近性質(zhì)時,提出了統(tǒng)一保持性質(zhì) (unified preserving property, UP) 的概念[71]。對于一個動理學(xué)方法,不僅應(yīng)該具有歐拉層面的漸近保持性質(zhì)(Asymptotic preserving property, AP),還應(yīng)該根據(jù)需要具有Navier-Stokes或更深層面的UP性質(zhì),這個過程實際上是通過時空耦合的思想實現(xiàn)的。該文分別用時空耦合DUGKS方法以及時空解耦CLR格式,對二維不可壓縮Taylor渦進行數(shù)值模擬,見圖13,具體見文獻[71]。特別指出的是,這一問題本質(zhì)上是多尺度問題。在多尺度問題及其相關(guān)的研究中,時空耦合策略是一條有效途徑。

        圖13 關(guān)于二維不可壓縮Taylor渦的DUGKS(時空耦合)和CLR(時空解耦)模擬的比較

        由于篇幅限制,本文只給出了有限體積方法基本原理的相容性準則以及通過GRP解法器實現(xiàn)相容性的技術(shù)細節(jié),對很多根本性質(zhì)還沒有涉及,如時空相容格式的熵穩(wěn)定性等。對于可壓縮流體力學(xué)來說,熵穩(wěn)定性對應(yīng)于熱力學(xué)相容性[7],這部分留在將來的文章中探討。

        客觀地說,時空耦合算法的研究非常不成熟,目前只在相對比較“純粹”的模型上進行分析和驗證,但是這一思想應(yīng)該加以推廣與應(yīng)用,比如如何將時空耦合算法的思想應(yīng)用于一般的時空關(guān)聯(lián)湍流模型[72]、粒子建模和時空耦合算法等。就算法本身來說,時空耦合的隱式格式研究還沒有開展,這一領(lǐng)域的研究應(yīng)該是下一階段的重點。

        在工程應(yīng)用中時空耦合的程序顯得更少,一方面是受限算法模塊的歷史延續(xù)性制約,另一方面是工程領(lǐng)域內(nèi)的“理性”力學(xué)越來越少。加強工程算法的科學(xué)性是需要目前亟待解決的觀念問題。

        后記。作為一個數(shù)學(xué)工作者,很忐忑地接受邀請在力學(xué)的專業(yè)期刊《空氣動力學(xué)學(xué)報》奉獻此類稿件,畢竟隔行如隔山。但想到數(shù)學(xué)、力學(xué)本就是“一家”,向力學(xué)家“學(xué)習(xí)”本就是“數(shù)學(xué)”的一大研究動機,心里稍微坦然點。

        致謝:此文的想法是在與下列合作者合作過程中逐步形成的,在此一并致謝,他們是Matania Ben-Artzi、齊進、汪玥、杜知方、雷昕、徐昆、郭照立、李啟兵,等。特別感謝曹貴瑜博士提供的關(guān)于各向同性可壓縮湍流模擬算例4.5, 感謝徐昆教授、郭照立教授和汪玥副研究員提出許多根本性改進意見,感謝陳海波給予了文字上的潤色。

        猜你喜歡
        黎曼通量時空
        非齊次二維Burgers方程的非自相似黎曼解的奇性結(jié)構(gòu)
        冬小麥田N2O通量研究
        跨越時空的相遇
        緊黎曼面上代數(shù)曲線的第二基本定理
        鏡中的時空穿梭
        數(shù)學(xué)奇才黎曼
        少兒科技(2019年4期)2019-01-19 09:01:15
        非等熵 Chaplygin氣體極限黎曼解關(guān)于擾動的依賴性
        玩一次時空大“穿越”
        時空之門
        緩釋型固體二氧化氯的制備及其釋放通量的影響因素
        化工進展(2015年6期)2015-11-13 00:26:29
        国产激情一区二区三区在线蜜臀| 影音先锋女人aa鲁色资源| 98久9在线 | 免费| 91超碰在线观看免费| 亚洲中字永久一区二区三区| 午夜秒播久久精品麻豆| 少妇久久久久久被弄到高潮| 91视频免费国产成人| 国产乱老熟视频乱老熟女1| 极品粉嫩小仙女高潮喷水操av| 久久精品国产亚洲av电影网| 欧美国产小视频| 亚洲二区三区四区太九| 国产人妖网站在线视频| 少妇私密会所按摩到高潮呻吟| 久久综合五月天| 少妇高潮太爽了免费网站| 日韩a级精品一区二区| 日韩亚洲av无码一区二区三区| 国产肉体XXXX裸体784大胆| 天堂麻豆精品在线观看| 亚洲色精品三区二区一区| 欧美v亚洲v日韩v最新在线| 人人爽亚洲aⅴ人人爽av人人片| 草逼视频污的网站免费| 日韩人妻无码一区二区三区久久| 男女扒开双腿猛进入免费看污| 亚洲美女性生活一级片| 蜜臀av一区二区三区免费观看| 中文字幕久久精品一二三区| 毛片无码高潮喷白浆视频| 不卡av一区二区在线| 亚洲国产精品无码久久| 亚洲乱码视频在线观看| 亚洲一区二区三区国产精品| 欧美大片va欧美在线播放| 中文字幕乱码人妻一区二区三区 | 在线欧美精品二区三区| 精品在线亚洲一区二区三区| 未发育成型小奶头毛片av| 亚洲av无码一区二区乱子伦as|