可壓縮多介質(zhì)流動問題的高精度數(shù)值模擬方法

何志偉, 田保林,2,*, 李 理, 李海鋒, 張又升, 孟寶清

(1. 北京應(yīng)用物理與計(jì)算數(shù)學(xué)研究所, 北京 100094; 2. 北京大學(xué) 應(yīng)用物理與技術(shù)研究中心, 北京 100871)

0 引 言

可壓縮多介質(zhì)流動問題，是指流場涉及多種物質(zhì)，且物質(zhì)之間存在相互耦合作用的流動現(xiàn)象。它廣泛存在于各種自然現(xiàn)象、工程應(yīng)用中。如超新星爆發(fā)、超燃沖壓發(fā)動機(jī)中的燃料混合、慣性約束聚變和武器物理等[1]，有著非常重要的科學(xué)研究意義和工程應(yīng)用價(jià)值。

在這種流動問題中，由于移動且可變形界面的存在，界面處流體性質(zhì)不連續(xù)給物理建模與數(shù)值方法帶來了極大的困難。目前，文獻(xiàn)中存在多種具有不同復(fù)雜度和適用范圍的物理模型及數(shù)值計(jì)算方法。從界面處理方式來看，它們可以分為兩類[2]：第一類將界面處理為清晰的間斷，即清晰界面模型/方法；第二類將將界面處理為多物質(zhì)相互擴(kuò)散而成的混合區(qū)，即擴(kuò)散界面模型/方法。

在擴(kuò)散界面模型中，界面被處理為多物質(zhì)相互擴(kuò)散而成的混合區(qū)，就像數(shù)值捕獲空氣動力學(xué)中的接觸間斷一樣[2]。雖然這類擴(kuò)散界面是一種人為做成的數(shù)值混合，但這類模型最重要的優(yōu)點(diǎn)就是可以在固定網(wǎng)格上，所有計(jì)算單元都采用相同的數(shù)值方案來求解，而不用考慮流場中存在什么樣類型的波(沖擊波、稀疏波，界面等)以及它們之間的相互作用。這類模型/方法可進(jìn)一步分為兩小類：基于歐拉/Navier-Stokes方程的擴(kuò)散界面模型和基于多相流方程的擴(kuò)散界面模型。

1) 基于單流體歐拉/Navier-Stokes方程的擴(kuò)散界面模型。這類模型是由傳統(tǒng)的歐拉/Navier-Stokes方程，附加每種物質(zhì)的(質(zhì)量)濃度方程構(gòu)成，每個(gè)計(jì)算網(wǎng)格即使含有多種流體組分仍具有相同的運(yùn)動速度。這類方法雖然形式上簡單，但卻存在一些問題：①過界面非物理振蕩。Abgrall通過數(shù)值實(shí)驗(yàn)表明[3]，對于兩種理想氣體的多介質(zhì)問題，直接采用這類模型會在界面處產(chǎn)生壓力振蕩。為了消除這個(gè)問題，他提出了非守恒的四方程模型(擬守恒方法)[3]及等效四方程模型[4]，在恰當(dāng)數(shù)值離散的基礎(chǔ)上，可以使用同種有限體積/差分算法來處理由界面分開的可壓縮多介質(zhì)流動問題[5]。但是這樣的操作最終導(dǎo)致方程組喪失守恒性。② 等溫等壓假設(shè)限制此模型的應(yīng)用范圍。這類模型通過兩種組分的狀態(tài)方程構(gòu)造混合物的狀態(tài)方程。壓強(qiáng)均衡會導(dǎo)致兩相的聲速必須耦合，而且顯然無法滿足某些真實(shí)物理情況。

2) 基于多相流方程的擴(kuò)散界面模型。多介質(zhì)流動是一類含有可區(qū)分大尺度物質(zhì)界面的多相流問題，因此很多學(xué)者也直接將多介質(zhì)流動稱為多相流動。多相流方程的擴(kuò)散界面模型是指將所有相視為連續(xù)介質(zhì)，利用各種平均方法(時(shí)間平均、空間平均、系綜平均)對各相的場方程進(jìn)行平均后得到的一類忽略了流動的一些微觀細(xì)節(jié)的宏觀均質(zhì)或平均混合模型[6]。由于使用平均方法，最終的模型會出現(xiàn)非守恒項(xiàng)(即最終的模型無法寫成散度形式)和較多的附加項(xiàng)(物質(zhì)之間的質(zhì)量、動量和能量傳遞等過程的宏觀?；?。對這些附加項(xiàng)建立不同的唯象模式，就形成了不同的模型(對于兩種介質(zhì)的流動，又稱雙流體模型或七方程模型)，如Baer-Nunziato模型(BN模型)[7]、Saurel-Abgrall模型(SA模型)[8]。這類模型具體如下一些特點(diǎn)：① 不同相之間是完全非平衡的，即每個(gè)相流體都有自己的密度、壓力、速度、溫度等，因此一維情況下共有7個(gè)方程。由于這類模型包含的方程個(gè)數(shù)多、波系復(fù)雜、存在非守恒項(xiàng)[10]，有時(shí)候并不能滿足嚴(yán)格雙曲特性，給數(shù)值求解帶來了很大挑戰(zhàn)。② 雙流體多相流模型一般只適用于稠密多相流情形，即可以離散可近似看做一種連續(xù)的流體(如大量顆粒、液滴等)。對于顆粒相稀疏或者處于稀疏到稠密過渡情形，雙流體多相流模型的適應(yīng)性還存在問題(為克服這一困難，最近我們提出了一種適用于稀疏到稠密跨流態(tài)可壓縮兩相流的模擬方法(CMP-PIC)[9])。③ 對于一些特殊情況，這類模型還可以退化為一系列簡化模型[1]。

在實(shí)際的工程應(yīng)用中，我們會根據(jù)具體問題的特殊性，針對性選擇一些簡化模型來開展相應(yīng)的數(shù)值模擬工作。但不管選擇何種類型的模型，為建立一套高效的高精度數(shù)值模擬算法，我們遇到的問題都是基本類似的。因此，本文內(nèi)容不局限于某種特殊模型，而將從數(shù)值框架的選擇、非守恒方程相容離散、高精度有界格式構(gòu)造、界面壓縮、界面-單介質(zhì)分區(qū)計(jì)算方法等5個(gè)維度出發(fā)，綜述近幾年我們在可壓縮多介質(zhì)流動問題高精度數(shù)值模擬工作中的一些進(jìn)展。

1 高精度數(shù)值模擬方法

1.1 基本要求

為了準(zhǔn)確模擬可壓縮多介質(zhì)流動問題，尤其是其中涉及的界面大變形運(yùn)動，需要探索和建立一套高效的高精度數(shù)值算法。同面向其他可壓縮流動問題的數(shù)值算法一樣，這套算法同樣需要滿足一定的要求，如守恒性、高分辨率、健壯性等。但可壓縮多介質(zhì)流動問題又有其特殊性，因此有一些額外的特殊要求，其中包括：1) 過界面不產(chǎn)生各種非物理的數(shù)值振蕩；2) 在長時(shí)間計(jì)算過程中，保持界面的銳利結(jié)構(gòu)。界面不能隨著計(jì)算時(shí)間的增加而不斷彌散變寬，最終影響流場其他結(jié)構(gòu)；3) 保界性?？蓧嚎s多介質(zhì)流動的物理模型中往往會涉及其他類型的物理量，如質(zhì)量分?jǐn)?shù)、體積分?jǐn)?shù)。這類物理量有著明確的上界和下界。因此，在對可壓縮多介質(zhì)流動問題進(jìn)行模擬時(shí)，不僅僅需要考慮熱力學(xué)量的保正問題，還得同時(shí)考慮如質(zhì)量分?jǐn)?shù)這類物理量的保界問題。從而保證程序的穩(wěn)定運(yùn)行和界面混合封閉模式的精準(zhǔn)執(zhí)行。

為了實(shí)現(xiàn)構(gòu)建滿足上述約束條件的數(shù)值方法，我們做了一系列探索和研究工作。下面詳細(xì)介紹研究的進(jìn)展以及一些具體工作。

1.2 數(shù)值框架的選擇

流體力學(xué)數(shù)值方法通常可以分為拉格朗日方法和歐拉方法兩種類型。拉格朗日方法，計(jì)算網(wǎng)格隨流體一起運(yùn)動，適宜計(jì)算一些變形不大的多介質(zhì)界面問題。歐拉方法，網(wǎng)格固定不動，能夠模擬界面大變形問題。本文主要關(guān)注歐拉數(shù)值方法。目前歐拉數(shù)值離散方法大致分為有限差分方法、有限體積方法、有限元方法三種基本框架。這些數(shù)值框架能否適用于可壓縮多介質(zhì)流動問題，需要做適應(yīng)性分析。而對這類框架開展適用性分析的指導(dǎo)思想是Abgrall提出的界面無振蕩分析方法[3]——過界面速度-壓力需要保持為常數(shù)。

1.2.1 過界面非物理振蕩現(xiàn)象

考慮如下經(jīng)典的滿足氣體狀態(tài)方程的界面運(yùn)動問題(見圖1)。左邊是一種氣體，如空氣，其比熱比為γL=1.4。右邊是另一種氣體，如氦氣，其比熱比為γL=1.667。跨過界面時(shí)，兩種介質(zhì)的壓力和速度均相同。

圖1 界面運(yùn)動問題

對于這個(gè)問題，采用如下模型(一維形式)[3]：

其中，ρ、u、E、p、YL分別為密度、速度、總能、壓力和左邊物質(zhì)的質(zhì)量分?jǐn)?shù)。其中質(zhì)量分?jǐn)?shù)滿足YL+YR=1，且混合物的內(nèi)能和壓力之間滿足：

(2)

其中WL和WR分別是作用物質(zhì)的摩爾質(zhì)量。對于這個(gè)模型，可以采用經(jīng)典的有限體積方法或有限差分方法進(jìn)行數(shù)值離散。但是大量數(shù)值結(jié)果表明，過界面的壓力和速度無法繼續(xù)保持常數(shù)。此即所謂的界面壓力振蕩現(xiàn)象。

目前，此現(xiàn)象的數(shù)值機(jī)理已經(jīng)很明確。即，傳統(tǒng)的數(shù)值方法具有數(shù)值黏性，會產(chǎn)生界面處的數(shù)值擴(kuò)散，導(dǎo)致界面處的不同介質(zhì)進(jìn)行了偽混合。在這種偽混合區(qū)域，數(shù)值方法計(jì)算出的壓力和溫度是不準(zhǔn)確的。此過程可以用圖2更清晰地進(jìn)行闡述。如圖2所示，初始狀態(tài)1和狀態(tài)2在等壓線上，但是由于數(shù)值黏性的存在，界面處會產(chǎn)生新的混合狀態(tài)3，而這個(gè)新狀態(tài)如果不在等壓線上，那么一定會出現(xiàn)壓力振蕩現(xiàn)象。

圖2 數(shù)值混合的不兼容性示意圖

在文獻(xiàn)[11]中，作者證明，當(dāng)介質(zhì)具有p=Aρe+B(其中A和B是常數(shù))這樣類型的狀態(tài)方程時(shí)，采用傳統(tǒng)的高精度差分方法是可以實(shí)現(xiàn)單介質(zhì)接觸間斷處壓力無振蕩屬性。但是當(dāng)介質(zhì)具有非線性狀態(tài)方程，或者多介質(zhì)混合問題(兩種氣體的混合造成混合物的狀態(tài)方程呈現(xiàn)非線性特點(diǎn))，接觸間斷/界面處一定會出現(xiàn)壓力振蕩現(xiàn)象。

1.2.2 數(shù)值方案引起的非物理振蕩機(jī)制分析

除了上述根本原因外，不同的數(shù)值方案會額外導(dǎo)致過界面非物理振蕩。為此，我們以有限差分方法為例[12]，開展了數(shù)值方案引起的非物理數(shù)值振蕩的分析工作。

基于通量分裂的高精度差分方法[13]，由于算法簡單，容易實(shí)現(xiàn)高精度計(jì)算，因此得到了大量應(yīng)用。我們在文獻(xiàn)[11-12,15]中詳細(xì)分析了此套方案應(yīng)用于多介質(zhì)流動問題時(shí)引起的非物理振蕩現(xiàn)象的數(shù)值機(jī)制。

考慮Abgrall提出的擬守恒方法[3]：

通過求解最后一個(gè)非守恒方程，避免了使用混合封閉模式(2)。此方法的本質(zhì)在于將原有的等溫假設(shè)修正為等壓假設(shè)。這種修正對多介質(zhì)界面運(yùn)動問題是合理的。對于擬守恒方法的控制方程，我們可以將其寫為：

其中，

u=(ρ,ρu,ρE,1/(γ-1))T

f=(ρu,ρu2+p,ρEu+pu,0)T

h=diag(0,0,0,u)

(5)

采用傳統(tǒng)的基于通量分裂的高精度差分方法[13]，可以得到如下的半離散形式：

對于逐特征場有限差分方法，需要在通量分裂后，通過局部特征分解后在每個(gè)特征場上應(yīng)用差分格式獲得相應(yīng)的特征場數(shù)值通量，進(jìn)而通過變換最終獲得物理空間的數(shù)值通量。

在文獻(xiàn)[12]，我們系統(tǒng)分析了逐方程有限差分方法在求解多介質(zhì)界面問題時(shí)存在的問題。

首先，假定采用的差分格式為線性格式，通過推導(dǎo)發(fā)現(xiàn)，傳統(tǒng)的Steger-Warming通量分裂方法[14]無法保持界面處壓力不振蕩。只有Lax-Friedrichs類型通量分裂方法可以滿足質(zhì)量-動量-能量方程的分裂得到的正通量之間/負(fù)通量之間保持相容性，其中g(shù)lobal Lax-Friedrichs通量分裂方法[13]是最簡單有效的方法。具體推導(dǎo)參見文獻(xiàn)[12]。

其次，在固定使用global Lax-Friedrichs通量分裂方法的情況下，有限差分方法需要采用非線性WENO差分格式[13]來計(jì)算各個(gè)方程的數(shù)值通量。但是，由于每個(gè)方程的通量的量級大小不同，每個(gè)方程所計(jì)算出的WENO格式的非線性權(quán)系數(shù)不同，進(jìn)而導(dǎo)致質(zhì)量-動量-能量方程所獲得的數(shù)值通量之間不滿足相容性，進(jìn)而造成界面處的壓力振蕩現(xiàn)象。具體推導(dǎo)參見文獻(xiàn)[12]。

進(jìn)一步，為了使方程之間的非線性WENO離散具有兼容性，我們提出了一套“公共權(quán)”方案。即質(zhì)量-動量-能量方程分別計(jì)算得到的WENO格式的非線性權(quán)修正為所有方程采用一套WENO權(quán)系數(shù)，這樣就可以保證質(zhì)量-動量-能量方程所獲得的數(shù)值通量之間的兼容性。在文獻(xiàn)[12]中，我們給出了一種具體的設(shè)計(jì)方法。但是更為魯棒的公共權(quán)方法需要進(jìn)一步研究。

最后，根據(jù)界面輸運(yùn)問題中速度和壓力要保持常數(shù)這一原則，推導(dǎo)出了非守恒方程相容離散形式為：

其中，φ=1/(γ-1)，α為global Lax-Friedrichs通量分裂方法中所取的全局常數(shù)。

對于雙曲系統(tǒng)，逐特征場的有限差分方法更加匹配雙曲系統(tǒng)。當(dāng)使用局部特征分解時(shí)，我們發(fā)現(xiàn)[11,15]：上述的WENO公共權(quán)方案的約束條件可以得到部分放松；質(zhì)量-動量-能量方程之間的相容性會得到更好的保持。在文獻(xiàn)[11]中，詳細(xì)證明了只要非線性特征場和描述界面演化的線性場采用公共權(quán)技術(shù)，就可以實(shí)現(xiàn)過界面壓力無振蕩屬性。但是設(shè)計(jì)滿足這種條件的公共權(quán)技術(shù)仍然是一個(gè)懸而未決的問題。

基于前期的工作，我們在文獻(xiàn)[15]中最終給出了避免界面壓力振蕩的有限差分方法的一般性的處理方案[15]。具體如下：

首先，將擬守恒方程寫成“守恒+源項(xiàng)”的如下具體形式：

其中,

采用上述的方法結(jié)合局部特征分解技術(shù)可以保證界面輸運(yùn)問題壓力-速度保持常數(shù)的屬性，具體證明見文獻(xiàn)[15]。

再其次，提出一個(gè)新的準(zhǔn)則[15]“單介質(zhì)保持準(zhǔn)則”——多介質(zhì)算法應(yīng)能自動保持單介質(zhì)流動屬性。利用這個(gè)性質(zhì)，推導(dǎo)出了源項(xiàng)的兼容離散形式為：

至此，就完成了新算法的推導(dǎo)過程。此推導(dǎo)過程清楚地表明,新算法不需要對非保守項(xiàng)進(jìn)行任何特殊處理。它成功地避免了設(shè)計(jì)WENO格式的公共權(quán)重技術(shù)。此外，它可以直接應(yīng)用于任何差分格式，而不局限于WENO格式。更詳細(xì)的內(nèi)容請參考文獻(xiàn)[15]。這里以如下的界面輸運(yùn)問題為例給出具體效果。此問題的初始條件為：

此問題的具體設(shè)置見文獻(xiàn)[15]。

從圖3中可以看出，采用五階精度的WENO格式[13]的差分方法直接求解守恒型多組分方程時(shí)，速度和壓力都會產(chǎn)生非物理數(shù)值振蕩。相反，當(dāng)使用新提出的算法求解擬守恒方法時(shí)，獲得了與精確解的吻合較好的數(shù)值結(jié)果。

圖3 界面輸運(yùn)問題的數(shù)值結(jié)果[15]

相比于只能采用通量分裂形式的有限差分方法，有限體積框架則具有更多的靈活性。如在重構(gòu)過程中，可以采用原始變量、守恒變量、特征量變等。在計(jì)算通量的過程中，可以采用各種黎曼求解器、矢通量分裂方法等[16]。Abgrall在最初的文章中就使用了基于原始變量重構(gòu)的有限體積方法。Johnsen[17]進(jìn)一步確認(rèn)了采用原始變量重構(gòu)的有限體積方法可以有效避免數(shù)值框架帶來的非物理數(shù)值振蕩。

這里需要特別指出，Deng等提出的基于半點(diǎn)值的緊致型有限差分方法[18-19]，由于可以借助有限體積中的重構(gòu)和通量計(jì)算策略，因此也比較容易拓展應(yīng)用于可壓縮多介質(zhì)流動問題的數(shù)值模擬中。對此問題，Nonomura等做了相關(guān)分析[20]。

1.3 非守恒方程相容離散

在可壓縮多介質(zhì)流動問題中，我們會遇到一些額外類型的方程。典型的如上節(jié)里提到的對流方程。這類方程與傳統(tǒng)的雙曲守恒方程不同，它們無法寫成守恒形式。從上一節(jié)的分析還可以看出：這類方法如果不恰當(dāng)離散，同樣會導(dǎo)致過界面的非物理振蕩。因此對于一個(gè)包含守恒方程和非守恒方程的可壓縮多介質(zhì)流動問題的控制方程，需要進(jìn)一步解決的問題是：當(dāng)守恒方程采用確定的雙曲守恒離散方法(如有限差分方法、有限體積方法)后，這些非守恒方程如何進(jìn)行相容離散。

其實(shí)，在上一節(jié)數(shù)值框架的分析過程中，已經(jīng)給出了一種有限差分方法框架下的離散方案。此外，對于有限體積方法，不同研究者也已經(jīng)提出了一些不同方法[15,17,21]。最近，我們進(jìn)一步拓展了此工作，給出了一套統(tǒng)一的離散方案。新方案自動兼容我們提出的有限差分離散方法[15]及其他部分研究者提出的有限體積離散方法[21]。此外，新方案還可以進(jìn)一步導(dǎo)出一些黎曼求解器對應(yīng)的新離散方案。關(guān)于這部分工作，詳見文獻(xiàn)[22]。

這里需要指出的是，針對一般性非守恒系統(tǒng)， Dal Maso、LeFloch和Murat的提出了DLM理論[49]，其中跨間斷的非守恒乘積項(xiàng)被定義為連接左右狀態(tài)的沿路徑積分?；谏鲜隼碚摪l(fā)展了所謂的路徑守恒離散方案[50-51]。但是，這類方法在路徑的選擇、收斂性等方面還沒有確切結(jié)論[52]，仍需進(jìn)一步發(fā)展。

1.4 高精度有界數(shù)值格式構(gòu)造

不同于單介質(zhì)歐拉方程和N-S方程，可壓縮多介質(zhì)流動的物理模型中會涉及其他類型的物理量，如質(zhì)量分?jǐn)?shù)、體積分?jǐn)?shù)。這類物理量有著明確的上界和下界，直接影響界面混合封閉模式的精準(zhǔn)執(zhí)行。因此，在對可壓縮多介質(zhì)流動問題進(jìn)行模擬時(shí)，我們不僅需要考慮熱力學(xué)量的保正問題，還得同時(shí)考慮如質(zhì)量分?jǐn)?shù)這類物理量的保界問題。因此，除了要求具備高精度、基本無振蕩屬性等基本性質(zhì)外，如何實(shí)現(xiàn)物理量的保界/保正要求是我們選擇和構(gòu)造可壓縮多介質(zhì)流動問題的高精度數(shù)值格式的額外指導(dǎo)原則。

目前，經(jīng)數(shù)值實(shí)驗(yàn)驗(yàn)證比較成功的高精度數(shù)值格式為WENO類型格式[13]。這類格式在實(shí)現(xiàn)高精度計(jì)算的同時(shí)，也實(shí)現(xiàn)了對各種間斷結(jié)構(gòu)的基本無振蕩捕捉。但是WENO格式同樣存在一些缺點(diǎn)：如WENO-JS格式耗散偏大，導(dǎo)致對一些多尺度問題(如激波-湍流邊界層干擾問題、聲學(xué)問題)的計(jì)算效率過低。為此，相關(guān)研究者進(jìn)行了大量的探索。提出了諸如WENO-M[23]、WENO-Z[24]、TENO[25]、雜交方法[26-27]等各種方法。

除了上述經(jīng)典方法外，其實(shí)還有其他一些類型格式。如：高精度保單調(diào)格式(monotonicty-preserving, MP)[28]。這種格式的基本思想為：首先采用任意類型高精度格式預(yù)估一個(gè)初始結(jié)果；其次通過構(gòu)造局地的上界和下界，形成一個(gè)局地的區(qū)間；最終通過探測初始結(jié)果是否落入此區(qū)間內(nèi)，決定是否對預(yù)估的初始結(jié)果進(jìn)行進(jìn)一步處理。這類格式繼承了TVD類型格式的保界思想，并進(jìn)一步拓展至高精度計(jì)算[29-30]。

因?yàn)閃ENO格式不具備保界性，Balsara在2000年提出了MPWENO格式[31]。其基本思想是用WENO格式預(yù)估初始結(jié)果，然后用MP格式構(gòu)造的局地上界和下界來對初始結(jié)果進(jìn)行進(jìn)一步限制，從而實(shí)現(xiàn)保界性。

在2016年，我們發(fā)現(xiàn)[32]，MP格式的局地上界和下界的設(shè)計(jì)方法會導(dǎo)致形成的區(qū)間過大，無法抑制一些數(shù)值振蕩。這些數(shù)值振蕩會持續(xù)傳播，最終導(dǎo)致MP格式的穩(wěn)定性較差。上述發(fā)現(xiàn)也得到了進(jìn)一步工作的證實(shí)[33]：“It was recently found that this limiting procedure is no longer accurate for long - term integrations. As indicated by He et al., the major problem arises mostly due to the use of considerable amounts of room for the MP constraint, resulting in either problem-dependent solutions or abundant room for the growth of spurious oscillations”。結(jié)合TVD有界理論，我們給出了一種改進(jìn)版本的MP格式[32]，即MP-R格式。數(shù)值實(shí)驗(yàn)表明改進(jìn)的局地保界方案可以較為顯著地提高最終方案的健壯性。

由于保界/保正要求對于模擬多介質(zhì)問題至關(guān)重要，我們自然希望計(jì)算法自動實(shí)現(xiàn)這個(gè)要求。因此，上述基于保界思想構(gòu)造的MP格式是一類比較適用于可壓縮多介質(zhì)流動問題的方法。但是，對于非線性問題而言，這種先驗(yàn)性的保界/保正算法不具備普適性。因此，為了進(jìn)一步提升算法的保界性和保正性，可以進(jìn)一步采用最新發(fā)展的保界/保正處理方法，如Postivity-Preserving方法[34]、Multi-dimensional Optimal Order Detection(MOOD)方法[35]。

1.5 界面壓縮

擴(kuò)散界面方法是將界面視為數(shù)值擴(kuò)散的狹窄區(qū)域，比如氣體動力學(xué)中的數(shù)值模擬所捕捉的接觸間斷。因此造成的結(jié)果是——數(shù)值界面會一直彌散變寬。并且隨著計(jì)算時(shí)間的增加，這種彌散變寬現(xiàn)象會越來越嚴(yán)重，最終可能會導(dǎo)致影響流動的其他物理結(jié)構(gòu)[36]。因此，如何避免或抑制界面的持續(xù)彌散變寬問題，一直是可壓縮多介質(zhì)流動問題界面捕捉算法的一個(gè)研究重點(diǎn)。

一個(gè)自然的想法是發(fā)展高精度間斷捕捉數(shù)值格式。在同樣的網(wǎng)格上，高精度間斷捕捉格式可實(shí)現(xiàn)高階精度計(jì)算，對間斷的分辨能力也顯著提高。但是文獻(xiàn)[37]中發(fā)現(xiàn)，采用高精度間斷捕捉格式，雖然可以在一定程度上提高間斷/界面的分辨能力，但依舊無法抑制長時(shí)計(jì)算情況下的界面彌散變寬問題。此外，由于高精度間斷捕捉格式引入的數(shù)值耗散較小，因此結(jié)果具有更多的數(shù)值振蕩。為了抑制這些數(shù)值振蕩，雖然可以采用局部特征分解方法來獲得基本無振蕩結(jié)果，但如果物質(zhì)由非常復(fù)雜的狀態(tài)方程(EOS)所描述，如何構(gòu)造一個(gè)穩(wěn)定且與物理模型兼容的局部特征分解方法仍需要進(jìn)一步探索[37]。

總之，采用高精度格式依舊無法抑制長時(shí)計(jì)算情況下的界面彌散變寬問題。此外，還會帶來可能的數(shù)值振蕩[37]。

因此，即使采用了高精度間斷捕捉格式來進(jìn)行計(jì)算，界面壓縮算法也是必須的[37]。目前，相關(guān)研究者提出了各種界面壓縮/反耗散算法[36-39,53]。這些方法可以大致分為四類：(1) 限制型下游格式方法[53]。這種方法的基本思想是在滿足TVD穩(wěn)定性理論的基礎(chǔ)上盡可能使用下游格式。目前僅適用于Lagrange-Remap方案。(2) 反擴(kuò)散方法[54]。這類方法的基本思想是直接求解反擴(kuò)散方程從而實(shí)現(xiàn)界面壓縮。(3) 基于修改重構(gòu)過程的代數(shù)型VOF方法，如基于雙曲正切函數(shù)發(fā)展的界面捕捉算法(THINC)[38-39]。(4) 界面重整化(又叫擴(kuò)散-銳利化)方法[36]。這類方法的基本思想是添加恰當(dāng)?shù)臄U(kuò)散項(xiàng)和銳利化項(xiàng)，從而保證界面的數(shù)值解是收斂的(寬度不隨計(jì)算時(shí)間增加而變寬)。目前這幾類方法推廣應(yīng)用到可壓縮多介質(zhì)問題時(shí)，均會遇到一個(gè)共同且依舊沒有徹底解決的問題——當(dāng)描述界面演化的方程(如體積分?jǐn)?shù)方程)采用上述界面壓縮方法后，其他方程(質(zhì)量-動量-能量方程)如何進(jìn)行相容的壓縮。前三類方法為了實(shí)現(xiàn)相容壓縮，目前已有的方案僅有一階精度，我們的數(shù)值實(shí)驗(yàn)證實(shí)過低的精度可能會抑制界面的大變形運(yùn)動。第四類方法的相容壓縮方法已經(jīng)從理論上得到解決(針對五方程模型[43])，但是目前這類方法添加的擴(kuò)散項(xiàng)和銳利化項(xiàng)均來源于level set方法[55]或相場方法[56]，如何設(shè)計(jì)和離散這些項(xiàng)(穩(wěn)定且保界)目前還沒有得到徹底解決[57]。

為了易于實(shí)現(xiàn)高精度計(jì)算且避免如何設(shè)計(jì)擴(kuò)散項(xiàng)和銳利化項(xiàng)的問題，在文獻(xiàn)[37]中，針對五方程模型[40-41]，我們提出了一種基于界面壓縮方案(第三類和第四類方法的雜交方案)。

考慮如下兩種物質(zhì)的五方程模型[40]：

(13)

其中，ρ、u、E、p分別為混合物密度、速度、總能和壓力，ρk、αk分為是第k種物質(zhì)的相密度和體積分?jǐn)?shù)。其中體積分?jǐn)?shù)滿足α1+α2=1。上述方程可以寫成(以一維情況為例)：

其中V=(ρ1α1,ρ2α2,ρu,ρE)T。

采用標(biāo)準(zhǔn)的有限體積方案，上述方程可以離散為：

提出的界面壓縮方案主要包括兩部分[37]：

1) 界面壓縮物理量的選擇。針對五方程模型，采用了如下的原始變量進(jìn)行(相密度、速度、壓力、體積分?jǐn)?shù))重構(gòu)[37](可進(jìn)一步投影至特征空間進(jìn)行重構(gòu))，來獲得單元邊界的左右重構(gòu)狀態(tài)：

(16)

對于這些物理量，可采用上述的各種高精度數(shù)值格式進(jìn)行計(jì)算。具體內(nèi)容詳見文獻(xiàn)[37]。

上述方案不具備界面壓縮效應(yīng)。為了引入界面的壓縮效應(yīng)(即采用THINC格式進(jìn)行計(jì)算)，需要對上述的標(biāo)準(zhǔn)方案進(jìn)行改進(jìn)。對何種物理量采用界面壓縮是需要考慮的首要問題。在文獻(xiàn)[37]中提出了如下方案：只對體積分?jǐn)?shù)采用THINC格式進(jìn)行計(jì)算，從而獲得的單元邊界的左右重構(gòu)狀態(tài)。

(17)

可以看到，在這里獲得的單元邊界的左右重構(gòu)狀態(tài)時(shí)，除了體積分?jǐn)?shù)采用了其他類型界面壓縮方法來計(jì)算，其他物理量均采用原高精度數(shù)值格式方案來進(jìn)行計(jì)算。

提出上述方案的基本思想主要基于兩點(diǎn)考慮：① 我們知道，物質(zhì)界面/接觸間斷是一類非常特殊的流動結(jié)構(gòu)。各種物理量過界面的表現(xiàn)性質(zhì)不一樣。我們可以把這些物理量分為快變量和慢變量。其中慢變量定義為過界面連續(xù)的物理量。與之相反，快變量定義為過界面存在間斷的變化劇烈的物理量。從這個(gè)分類可以看出，上述五方程模型中，相密度、速度、壓力均是慢變量，而只有體積分?jǐn)?shù)是快變量。② 僅對快變量做界面壓縮方法。這個(gè)觀點(diǎn)是基于如下事實(shí)：我們發(fā)現(xiàn)，界面壓縮方法，尤其是THINC類型的具有反耗散性質(zhì)的方法，存在一個(gè)問題：這類方法可以使光滑的物理結(jié)構(gòu)被反物理地壓縮成鋸齒狀間斷結(jié)構(gòu)[42]。如圖4所示，當(dāng)采用THINC格式對Burgers方程進(jìn)行計(jì)算時(shí)，反耗散方法對間斷結(jié)構(gòu)的分辨率非常高。但如果在光滑區(qū)也采用THINC方法，光滑的物理結(jié)構(gòu)亦被非物理地壓縮為間斷結(jié)構(gòu)。 關(guān)于這個(gè)問題的討論，詳見最新的文獻(xiàn)[42]?；谶@樣的實(shí)際情況，可以得出如下結(jié)論：界面壓縮只能應(yīng)用于快變量。正是基于上述觀點(diǎn)，我們提出了僅對體積分?jǐn)?shù)進(jìn)行壓縮處理的界面壓縮方案。

圖4 采用THINC格式計(jì)算Burgers方程[42]

2) 相容性。大量數(shù)值實(shí)驗(yàn)表明：上述方案在工程計(jì)算中可以取得較好的計(jì)算結(jié)果，這種方案也可以應(yīng)用于其他模型。我們所需要做的僅僅是：① 選擇重構(gòu)的原始變量；② 對這些原始變量進(jìn)行區(qū)分：快變量還是慢變量；③ 對快變量進(jìn)行界面壓縮處理。表1給出了四方程模型和五方程模型的一些典型快變量和慢變量。

(19)

有了上述的體積分?jǐn)?shù)方程的源項(xiàng)，就可以直接利用Tiwari的理論結(jié)果：在質(zhì)量方程、動量方程、能量方程上加入能嚴(yán)格保證相容的源項(xiàng)。具體表達(dá)式詳見文獻(xiàn)[37]。圖5給出了激波-氣泡干擾問題(具體設(shè)置參見文獻(xiàn)[37])的模擬結(jié)果(有/無界面壓縮情況)。

從圖5中可以看到，采用了界面壓縮方法后，界面的分辨率得到顯著提高。此外，數(shù)值結(jié)果亦表明界面的寬度基本不隨時(shí)間增加而變寬[37]。

圖5 Air-R22激波-氣泡干擾問題結(jié)果[37]

但是，上述方法也存在一定的缺陷——最終的方案不是嚴(yán)格守恒[43]。此外，文獻(xiàn)[37]中指出，在實(shí)際應(yīng)用中，不同近似黎曼求解器帶來的數(shù)值耗散不一樣。而由這部分因素帶來的數(shù)值耗散不能僅通過修改重構(gòu)方法而徹底回避。因此，如何構(gòu)造守恒且相容的界面壓縮方法仍值得進(jìn)一步深入研究。

1.6 界面-單介質(zhì)分區(qū)計(jì)算方案

開展可壓縮多介質(zhì)流動問題數(shù)值模擬時(shí)，由于我們采用的模型為擴(kuò)散界面模型。其本質(zhì)假設(shè)是每一空間點(diǎn)上同時(shí)存在多種介質(zhì)，介質(zhì)之間滿足一定的關(guān)系(如壓力平衡-溫度非平衡假設(shè))。因此為了避免模型的退化，必須在純流體中添加一個(gè)小的體積分?jǐn)?shù)[37]。

在文獻(xiàn)[44-45]中我們發(fā)現(xiàn)，這種人為處理方法會給計(jì)算帶來非物理的波。這種現(xiàn)象在考慮彈塑性-流體相互作用的問題時(shí)尤為明顯：① 在流體中添加小部分固體，會讓流體產(chǎn)生抵抗剪切的能力；② 固體和流體的聲速差異很大，即使是一小部分的固體，也會影響純流體的聲速。在此發(fā)現(xiàn)基礎(chǔ)上，我們提出了一個(gè)改進(jìn)的壓力平衡-溫度非平衡擴(kuò)散界面模型，該模型在等壓縮假設(shè)的基礎(chǔ)上，導(dǎo)出了和體積分?jǐn)?shù)無關(guān)的相密度方程。采用該模型，即使體積分?jǐn)?shù)為0，也不會導(dǎo)致方程退化，避免了非物理波的產(chǎn)生。圖6給出了新模型和舊模型模擬純剪切問題的對比，新模型可以保持住剪切間斷，不產(chǎn)生額外的剪應(yīng)力。而采用原模型，必須在方程中添加小的體積分?jǐn)?shù)，導(dǎo)致虛假剪切應(yīng)力的產(chǎn)生，以及剪切波向兩側(cè)擴(kuò)散。但上述方案是否能推廣到實(shí)際工程應(yīng)用中，還需進(jìn)一步探索。

圖6 純剪切流數(shù)值模擬結(jié)果[44-45]

另一種解決方案是采用區(qū)域分解這類典型的多尺度問題計(jì)算方案。即將計(jì)算區(qū)域分解為純物質(zhì)區(qū)和界面區(qū)，不同區(qū)域采用不同的方案進(jìn)行計(jì)算。其實(shí)，在上一節(jié)構(gòu)造界面壓縮算法過程中，已經(jīng)人為地區(qū)分了界面區(qū)域和純物質(zhì)區(qū)域。相關(guān)數(shù)值實(shí)驗(yàn)亦證實(shí)了此方案的可行性。但界面區(qū)域與純物質(zhì)區(qū)域的邊界處是否會人為引入一些非物理波，這一問題還需要更深入的研究。

2 應(yīng) 用

通過上述多個(gè)維度的工作，建立了一套適用于可壓縮多介質(zhì)流動問題的低耗散、基本無數(shù)值振蕩的高精度歐拉方法。相關(guān)數(shù)值方法成果已集成至工程問題數(shù)值模擬軟件，為相關(guān)工程任務(wù)提供了重要技術(shù)支撐。在具體應(yīng)用方面，不同的問題往往具有自身的特殊性，需要有針對性地選擇不同的具體離散方案進(jìn)行應(yīng)用。本小節(jié)給出幾個(gè)具體應(yīng)用算例。

第一個(gè)算例是彈塑性介質(zhì)的高速沖擊問題。金屬等固體材料在受到高速沖擊后，會發(fā)生彈性到塑性的轉(zhuǎn)變，甚至產(chǎn)生斷裂、破碎等復(fù)雜的物理現(xiàn)象，對數(shù)值模擬提出了極大的挑戰(zhàn)。我們基于上述五方程模型，采用有限差分方法、高精度(MP)WENO格式、非守恒離散方法以及結(jié)合界面壓縮技術(shù)，并進(jìn)一步對彈塑性效應(yīng)進(jìn)行了建模和數(shù)值離散，成功地實(shí)現(xiàn)了對此類問題的模擬。圖7給出了如下問題的數(shù)值模擬結(jié)果：銅塊以800 m/s的初始速度沖擊一個(gè)靜止的銅靶，銅靶先后經(jīng)歷彈性變形、塑性變形，最終發(fā)生斷裂和射流現(xiàn)象。相關(guān)結(jié)果和已有的文獻(xiàn)結(jié)果一致。

(a) t=0 ms (b) t=1.1 ms

第二個(gè)算例是含顆粒Richtmyer-Meshkov(RM)不穩(wěn)定性問題[9]。激波驅(qū)動的多相流動(如氣固多相流)中經(jīng)常伴隨著多相流RM不穩(wěn)定性的的發(fā)生。在氣固兩相流動體系中，顆粒相在流動體系中顯著影響了界面的動力學(xué)行為。我們基于上述五方程模型，采用有限體積方法、高精度(MP)WENO格式、非守恒離散方法，并進(jìn)一步對跨流態(tài)顆粒運(yùn)動進(jìn)行了數(shù)值建模(CMP-PIC 方法[9])，成功實(shí)現(xiàn)了對此類問題的數(shù)值模擬。圖8具體給出了模擬的含空氣/SF6的多相RM問題(其中顆粒均勻分布于空間)問題的具體情況：空氣以及SF6界面存在單模擾動，入射激波由空氣進(jìn)入到SF6區(qū)域。顆粒均勻布置于計(jì)算域，其分布區(qū)域?yàn)椋?.01 m

圖8 計(jì)算域示意圖

(a) t=0.0301 ms，界面演化過程

圖10 不同顆粒體積分?jǐn)?shù)下混合區(qū)寬度隨時(shí)間演化

表2 不同體積分?jǐn)?shù)下混合區(qū)寬度數(shù)值及理論預(yù)測比較

第三個(gè)算例為含激波二次沖擊的RM湍流混合問題[46]。含激波沖擊的湍流混合問題是工程應(yīng)用及自然現(xiàn)象中的基本物理過程，其典型特點(diǎn)是波系與混合區(qū)的多次作用。針對該問題，開展了直接數(shù)值模擬研究，采用了四方程模型，并且由于考慮物理黏性，所以沒有使用界面壓縮技術(shù)。并且采用的網(wǎng)格規(guī)模達(dá)到11億。初始激波從輕流體向重流體傳播，透射激波在下游壁面反彈后二次作用于混合區(qū)，混合區(qū)快速湍流化(如圖11所示)，隨后混合區(qū)依次被稀疏波和壓縮波加速。稀疏波和壓縮波作用期間，混合區(qū)內(nèi)形成了明顯的平均速度梯度，混合區(qū)被拉伸或者壓縮。采用混合寬度無量綱化的平均組分剖面在二次激波作用后基本保持不變。借助高精度直接數(shù)值模擬結(jié)果，對多介質(zhì)湍流的生成、衰減過程進(jìn)行了收支平衡分析。結(jié)果表明，在稀疏波和壓縮波作用階段，速度梯度生成項(xiàng)和壓力梯度生成項(xiàng)是主導(dǎo)項(xiàng)。壓力梯度生成項(xiàng)和速度梯度生成項(xiàng)的貢獻(xiàn)相反，前者的幅值更大。因此，稀疏波導(dǎo)致湍動能增強(qiáng)而壓縮波導(dǎo)致湍動能減弱。進(jìn)一步統(tǒng)計(jì)了含能結(jié)構(gòu)的尺度變化。結(jié)果表明，混合寬度的增長伴隨著大尺度結(jié)構(gòu)的形成，稀疏波階段，結(jié)構(gòu)被縱向拉伸，而壓縮波作用階段，結(jié)構(gòu)被縱向壓縮。

圖11 三維質(zhì)量組分等值面的時(shí)間演化[46]. (a) t=0.0005 s和(b) t=0.0002 s為二次激波作用前；(c) t=0.003 s 和 (d) t=0.0045 s 為二次激波作用后

最后一個(gè)算例是匯聚幾何不穩(wěn)定及湍流混合問題。在實(shí)際工程應(yīng)用中，RM湍流混合往往發(fā)生在不斷向中心匯聚(如柱面、球面等)的物質(zhì)界面上。相比于平面界面上的RM湍流混合，匯聚界面上的RM湍流混合區(qū)在演化過程中會伴隨有其周向長度的改變，而這種改變將會深刻影響混合區(qū)的后續(xù)發(fā)展。針對該問題，基于集成上述方案所開發(fā)的模擬程序(具體方案同上一個(gè)算例)，對圓柱界面上的RM不穩(wěn)定性開展了數(shù)值模擬研究。初始激波由外部的重流體向內(nèi)部的輕流體傳播，并在隨后往復(fù)運(yùn)動于幾何中心和混合區(qū)之間，使得激波多次作用于混合區(qū)，混合區(qū)快速湍流化(如圖12所示)。對于混合區(qū)寬度，其演化被分解為混合區(qū)整體的拉伸壓縮和兩種流體組分相互穿插兩部分的貢獻(xiàn)。通過與平面算例進(jìn)行對比，我們發(fā)現(xiàn)向內(nèi)匯聚的物質(zhì)界面將整體拉伸混合區(qū)，并在同時(shí)促進(jìn)兩種流體的相互穿插。在激波第三次沖擊混合區(qū)之后，輕重流體穿插效應(yīng)所貢獻(xiàn)的混合寬度W滿足W∝(t-t0)θ，這與平面工況中混合寬度的自相似增長相類似。對于混合區(qū)內(nèi)部的平均組分剖面，經(jīng)過混合區(qū)寬度歸一化的曲線近似重合。通過統(tǒng)計(jì)混合區(qū)前沿位置的分布，本文計(jì)算了氣泡和尖釘結(jié)構(gòu)各自的直徑演化。數(shù)值結(jié)果表明，氣泡和尖釘?shù)闹睆脚c其對應(yīng)的混合區(qū)高度之比在演化后期趨近于常數(shù)。

圖12 三維質(zhì)量組分等值面的時(shí)間演化. (a) t=0.006 s為二次激波沖擊后，三次激波沖擊前和(b) t=0.012 s為三次激波作用后；(c) t=0.020 s為激波強(qiáng)度衰減至可以忽略后

從以上的算例可以發(fā)現(xiàn)，實(shí)際工程問題往往都是耦合多種物理因素和過程的多尺度多物理復(fù)雜流動問題。因此為了更好地模擬實(shí)際工程問題，我們需要進(jìn)一步考慮黏性[58-59]、熱傳導(dǎo)[60-61]等耦合物理效應(yīng)(更多物理效應(yīng)見Saurel的系列工作及其綜述[62])，以及彈塑性[63-64]、輻射擴(kuò)散[65]等多物理耦合過程。即需進(jìn)一步研究上述可壓縮多介質(zhì)數(shù)值方法在這類多物理效應(yīng)/過程情況下的適用性問題。但是到目前為止，這一方面的研究還遠(yuǎn)不夠成熟，我們的研究工作也正在進(jìn)一步開展。

3 結(jié) 論

建立一套高效的可壓縮多介質(zhì)流動問題的高精度數(shù)值模擬方案是一個(gè)系統(tǒng)工程，其中涉及的因素較多，主要有：數(shù)值框架選擇、非守恒方程相容離散、高精度有界格式構(gòu)造、界面壓縮、界面-單介質(zhì)分區(qū)計(jì)算方法等。本文從上述維度出發(fā)，綜述了本研究方向以及課題組的一些工作。主要結(jié)論有：

1) 不同數(shù)值方案可能導(dǎo)致過界面出現(xiàn)非物理振蕩現(xiàn)象。因此，不管選用何種類型的數(shù)值方案，必需對此方案在可壓縮多介質(zhì)流動問題中的適用性先進(jìn)行分析。為此，我們針對有限差分方法在此類問題中的適用性進(jìn)行了分析。

2) 非守恒類型方程必需現(xiàn)相容離散。當(dāng)守恒方程采用確定的雙曲守恒離散方法(如有限差分方法、有限體積方法)后，非守恒方程如何進(jìn)行相容離散。這類方法如果離散不恰當(dāng)，同樣會導(dǎo)致過界面的非物理振蕩。針對此問題，我們提出了一套相對普適的離散方案。

3) 保界/保正是構(gòu)造適用于可壓縮多介質(zhì)流動問題的高精度格式的重要要求。因?yàn)樵趯蓧嚎s多介質(zhì)流動問題進(jìn)行模擬時(shí)，不僅需要考慮熱力學(xué)量的保正問題，還得同時(shí)考慮如質(zhì)量分?jǐn)?shù)這類物理量的保界問題。為此，我們改進(jìn)了經(jīng)典的保單調(diào)限制器，進(jìn)一步提升了此算法的保界性能。

4) 在擴(kuò)散界面模型框架下，界面會一直數(shù)值彌散變寬，并且隨著計(jì)算時(shí)間的增加，這種彌散變寬現(xiàn)象會越來越嚴(yán)重。為了避免或抑制界面的持續(xù)彌散變寬現(xiàn)象，我們提出了一種界面壓縮方法。

5) 界面-單介質(zhì)分區(qū)計(jì)算方案。為了避免模型的退化，必須在純流體中添加一個(gè)小的體積分?jǐn)?shù)。我們發(fā)現(xiàn)，文獻(xiàn)上常用的這種方法會帶來非物理的波。因此，我們建議采用區(qū)域分解的這種多尺度計(jì)算方法，開展界面-單介質(zhì)的分區(qū)計(jì)算。但如何保證邊界處不產(chǎn)生非物理振蕩需要進(jìn)一步研究。

通過上述多個(gè)維度的工作，我們建立了一套適用于可壓縮多介質(zhì)流動的低耗散、基本無數(shù)值振蕩的高精度歐拉方法，相關(guān)數(shù)值方法成果已集成至工程問題數(shù)值模擬軟件，為相關(guān)工程任務(wù)提供了重要技術(shù)支撐。

由于實(shí)際工程問題往往都是耦合其他效應(yīng)和過程的多尺度多物理問題。因此，為了進(jìn)一步提高模擬的置信度，需要進(jìn)一步研究的工作有：

1)多物理效應(yīng)和多物理過程的建模與計(jì)算方案。即需進(jìn)一步研究“在上述多尺度多物理問題情況下的多介質(zhì)數(shù)值方法的適用性”問題。比如當(dāng)存在物理黏性時(shí)，界面壓縮能否與物理黏性自洽識別、可溶與不可溶界面的能否自動區(qū)分等問題。 典型的進(jìn)展如最近Xiao等最近提出了Boundary Variation Diminishing (BVD)算法[66]，在此方面做了較為成功的嘗試。

2)極端條件下的數(shù)值算法研究。相比單介質(zhì)問題，多介質(zhì)問題往往更容易遇到各種極端條件(如復(fù)雜狀態(tài)方程、非牛頓流體、固體的彈塑性效應(yīng)與超低馬赫數(shù)效應(yīng)、物性差異巨大等)，目前的計(jì)算流體數(shù)值方法能否適用于這些極端情況、以及如何優(yōu)化改進(jìn)現(xiàn)有算法值得進(jìn)一步研究。

3)高效的大規(guī)模計(jì)算軟件研發(fā)?？蓧嚎s多介質(zhì)問題中流場每一空間點(diǎn)可能為單介質(zhì)或多介質(zhì)，如何構(gòu)造一種高效表征多介質(zhì)流動問題的數(shù)據(jù)結(jié)構(gòu)，以及大規(guī)模并行計(jì)算情況的下的如何荷載平衡等問題，均對軟件研發(fā)提出了挑戰(zhàn)。此外，多介質(zhì)問題的各類模型并沒有經(jīng)過大量大規(guī)模并行計(jì)算的測試與檢驗(yàn)，進(jìn)一步對研發(fā)高效的大規(guī)模計(jì)算軟件提出了需求。