ALE方法中一種新的二階保界守恒重映算法

李珍珍，蔚喜軍，賈祖朋，安 娜，黃朝寶

（1.鄭州輕工業(yè)學(xué)院數(shù)學(xué)與信息科學(xué)學(xué)院，河南鄭州 450002；2.北京應(yīng)用物理與計(jì)算數(shù)學(xué)研究所計(jì)算物理實(shí)驗(yàn)室，北京 100088；3.中國工程物理研究院研究生部，北京 10008）

ALE方法中一種新的二階保界守恒重映算法

李珍珍1，蔚喜軍2，＊，賈祖朋2，安 娜3，黃朝寶3

（1.鄭州輕工業(yè)學(xué)院數(shù)學(xué)與信息科學(xué)學(xué)院，河南鄭州 450002；2.北京應(yīng)用物理與計(jì)算數(shù)學(xué)研究所計(jì)算物理實(shí)驗(yàn)室，北京 100088；3.中國工程物理研究院研究生部，北京 10008）

在用拉格朗日格式求解流體力學(xué)問題時(shí)，隨著時(shí)間的推進(jìn)，計(jì)算網(wǎng)格會(huì)扭曲變形，影響格式精度，甚至導(dǎo)致計(jì)算中斷。因此，要在網(wǎng)格變形較大時(shí)進(jìn)行網(wǎng)格重分和物理量重映，以保證網(wǎng)格質(zhì)量和格式精度。針對(duì)間斷有限元方法求解流體力學(xué)問題的二階拉格朗日格式，給出了一種守恒重映算法。該重映算法包括兩步：第一步是用已有重映方法計(jì)算新網(wǎng)格上的單元平均值，并用相應(yīng)修補(bǔ)算法對(duì)單元平均值進(jìn)行調(diào)整，保證單元平均值的保界性；第二步是由已得到的新單元平均值重構(gòu)出新網(wǎng)格上分片一次多項(xiàng)式，再使用Van Leer限制器對(duì)新網(wǎng)格上的梯度進(jìn)行限制，使之不出現(xiàn)新的極值。最后用數(shù)值算例表明了該重映算法的保界性和二階收斂性。

守恒算法；重映；ALE方法

0 引 言

通常ALE（Arbitrary Lagrangian-Eulerian）方法求解流體力學(xué)方程組時(shí)包含三個(gè)步驟：1）拉氏計(jì)算步，即求解Lagrangian形式的流體力學(xué)方程組，將時(shí)間推進(jìn)一步。在此過程中，隨著流體的運(yùn)動(dòng)，網(wǎng)格會(huì)發(fā)生變形，甚至可能出現(xiàn)網(wǎng)格扭曲或翻轉(zhuǎn)的情況；2）網(wǎng)格重分步，即生成一個(gè)質(zhì)量較好的新網(wǎng)格；3）物理量重映步，即將舊網(wǎng)格上的物理量映射到新網(wǎng)格上。Li等人［1－3］研究了拉氏計(jì)算步，主要是針對(duì)可壓縮歐拉方程的半拉格朗日微分形式，推導(dǎo)出它的積分弱形式，并用間斷有限元方法對(duì)其進(jìn)行空間離散，最后得到了中心型拉格朗日格式，避免了求解拉格朗日坐標(biāo)系映射到歐拉坐標(biāo)系下的雅可比矩陣，并且計(jì)算網(wǎng)格隨著流體運(yùn)動(dòng)。當(dāng)網(wǎng)格變形較小時(shí)，可以只進(jìn)行拉氏計(jì)算步，來追蹤物質(zhì)界面；但當(dāng)網(wǎng)格變形較大時(shí)，需要將網(wǎng)格重分，使重分后的網(wǎng)格具有更好的幾何品質(zhì)，并將拉氏計(jì)算步得到的物理量映射到新網(wǎng)格上?？梢娢锢砹恐赜巢绞茿LE［4－5］方法的重要組成部分。

重映算法可以分為兩類：積分重映和對(duì)流重映。積分重映對(duì)新舊網(wǎng)格的關(guān)聯(lián)性要求不高，可以出現(xiàn)跨網(wǎng)格的情形。所謂跨網(wǎng)格是指下面的情形之一：一是新舊網(wǎng)格有相同的拓?fù)浣Y(jié)構(gòu)，單元數(shù)目相等，但新網(wǎng)格移動(dòng)的距離較大，出現(xiàn)了跨網(wǎng)格情形；二是新舊網(wǎng)格有相同的拓?fù)浣Y(jié)構(gòu)，但單元數(shù)目不相等，比如說粗網(wǎng)格重分為細(xì)網(wǎng)格；三是新舊網(wǎng)格的拓?fù)浣Y(jié)構(gòu)不一樣，比如說結(jié)構(gòu)四邊形網(wǎng)格重分為任意多邊形網(wǎng)格。積分重映比較自然，但需要計(jì)算新舊網(wǎng)格的交點(diǎn)，計(jì)算量較大。對(duì)流重映的計(jì)算效率比較高，但它要求新舊網(wǎng)格上物理量的交換只在相鄰單元中進(jìn)行，因此要求新舊網(wǎng)格離得比較近，不能出現(xiàn)跨網(wǎng)格的情形。在實(shí)際計(jì)算中，我們還希望重映算法滿足以下幾個(gè)條件：一是守恒性，對(duì)于氣動(dòng)力學(xué)方程來講，要求滿足質(zhì)量、動(dòng)量和能量守恒；二是保界性，是指重映后物理量的值不能超過重映前物理量的最小和最大值；三是相容性，當(dāng)新舊網(wǎng)格重合時(shí)，重映前后物理量的值不變；四是保證二階或二階以上收斂。

Dukowicz［6］和Miller［7］分別給出了二維空間中的二階重映方法；Grandy［8］給出了三維空間中的一階重映方法；Dukowicz［9］和Mosso［10］分別在各自的文章中給出了三維空間中的二階方法；Leslie［11］和Leonard［12］則分別利用一維重映算法給出了多維空間中重映的分離算法。在這些重映方法中，具有代表性的一類局部ALE重映方法是通過適當(dāng)?shù)膶?duì)流方程來構(gòu)造重映過程。研究發(fā)現(xiàn)這類方法可以克服許多實(shí)際困難，所以在處理一些特殊應(yīng)用問題上十分有效。此類方法的文獻(xiàn)有很多［10，13－21］，其中Dukowicz［20］表明當(dāng)時(shí)間步充分小時(shí)，重映過程可以寫成一個(gè)通量形式的對(duì)流算法，這使得重映方法變?yōu)橐粋€(gè)簡單、有效的標(biāo)準(zhǔn)對(duì)流格式，并且這種格式保留了一般重映方法的很多優(yōu)勢。

雖然上述用對(duì)流方程描述的重映算法有很多優(yōu)勢，但是對(duì)流方程和守恒重映之間的聯(lián)系不易理解。Margolin和Shashkov［22］構(gòu)造了與對(duì)流方程無關(guān)的重映算法，首先將新單元表示為新舊單元相交部分的組合，并以此為基礎(chǔ)給出了一種相應(yīng)的守恒積分重映算法。接下來Margolin和Shashkov又將新單元質(zhì)量表示為相應(yīng)舊單元質(zhì)量加上其與周圍單元的質(zhì)量交換，再通過對(duì)密度函數(shù)的分片常數(shù)重構(gòu)，得到了一階格式的對(duì)流重映算法。為了得到二階精度的重映算法，Margolin和Shashkov分析了上述一階對(duì)流重映算法的誤差，通過對(duì)誤差做補(bǔ)償，得到了二階、保正的重映算法。Kucharik等人［23］是首先利用已知的舊網(wǎng)格上單元平均值，重構(gòu)出舊網(wǎng)格上密度的分片一次多項(xiàng)式逼近，再利用該分片一次多項(xiàng)式重映出新網(wǎng)格上密度的單元平均值，最后通過修補(bǔ)算法得到了一種保線性且保界的重映格式。

上述重映算法［22－23］可總結(jié)為：已知舊網(wǎng)格上單元平均值，重映后得到新網(wǎng)格上單元平均值。當(dāng)ALE方法的拉氏步采用有限體積法離散時(shí)（例如文獻(xiàn)［24］），只需知道重分網(wǎng)格上物理量的單元平均值，就可以更新tn＋1時(shí)刻物理量的值，因此可用文獻(xiàn)［22］和文獻(xiàn)［23］的重映算法。但當(dāng)ALE方法的拉氏步采用間斷有限元法離散時(shí)（例如文獻(xiàn)［1-3］），必須知道重分網(wǎng)格上物理量的分片多項(xiàng)式，才能更新出tn＋1時(shí)刻物理量的值，此時(shí)上述重映算法不再適用。不過我們可以以上述重映算法為基礎(chǔ)，先重映出物理量的單元平均值，再通過重構(gòu)或者其他方法得到重分網(wǎng)格上的分片多項(xiàng)式。

因此當(dāng)ALE方法的拉氏步采用間斷有限元法離散時(shí)，以文獻(xiàn)［1-3］中的二階拉氏格式為例，本文給出了相應(yīng)重映步的對(duì)流守恒重映算法。首先在已知舊網(wǎng)格上分片一次多項(xiàng)式的基礎(chǔ)上，采用文獻(xiàn)［22］和文獻(xiàn)［23］的重映和修補(bǔ)算法，得到重分網(wǎng)格上物理量的單元平均值，再通過極小化問題和Van Leer限制器，得到重分網(wǎng)格上物理量的梯度，由此得到重分網(wǎng)格上的分片一次多項(xiàng)式。

1 問題描述

考慮二維計(jì)算區(qū)域Ω，假設(shè)它是一個(gè)任意的多邊形區(qū)域，給定Ω上一個(gè)多邊形剖分｛Ci，i＝1，...，imax｝，其中單元Ci可以是非凸的。與單元Ci共點(diǎn)或共邊的單元稱為它的鄰居，單元Ci所有鄰居的集合稱為它的鄰域，記為C（Ci）。我們稱｛Ci，i＝1，...，imax｝為舊網(wǎng)格，舊網(wǎng)格通過節(jié)點(diǎn)的小位移移動(dòng)得到新網(wǎng)格｛，i＝1，...，imax｝。在接下來的討論中，假設(shè)新單元包含在舊單元Ci和它的鄰域中。下面以密度函數(shù)為例，敘述重映問題。

假設(shè)有一個(gè)定義在區(qū)域Ω上的正函數(shù)ρ（x，y）＞0，稱為密度，但是ρ（x，y）的具體表達(dá)式未知，僅知道ρ（x，y）在舊網(wǎng)格上的分片一次多項(xiàng)式逼近：

這表明近似函數(shù)ρi（x，y）保證ρ（x，y）單元平均值不變，是單元密度平均值。區(qū)域Ω上總質(zhì)量為：

本文最終目的是在保證總質(zhì)量M 不變的情況下，找到新網(wǎng)格上密度的分片一次多項(xiàng)式逼近：

2 重映新網(wǎng)格上單元平均值

2.1 重映單元平均值

本節(jié)采用文獻(xiàn)［22］中對(duì)流重映算法，計(jì)算新單元上密度的平均值。下面我們簡要回顧一下該重映算法，具體細(xì)節(jié)參見文獻(xiàn)［22］。

首先如文獻(xiàn)［22］給定一組結(jié)構(gòu)四邊形網(wǎng)格，因此需適當(dāng)調(diào)整相應(yīng)記號(hào)。網(wǎng)格頂點(diǎn)記為Pi，j＝（xij，yij），i＝1，...，m；j＝1，...，n。單元標(biāo)號(hào)記為（i＋，單元的四個(gè)頂點(diǎn)是：｛Pi，j，Pi＋1，j，Pi＋1，j＋1，Pi，j＋1｝。那么舊單元上密度函數(shù)的一次多項(xiàng)式逼近相應(yīng)表示為：

假設(shè)定向邊Fi，j＋1＝｛Pi，j＋1，Pi，j｝是單元的公共邊，位移向量場使頂點(diǎn)Pi，j移動(dòng)到新位置i，j，Pi，j＋1移動(dòng)到新位置i，j＋1，因此形成一個(gè)定向四邊形i，j＋1，Pi，j＋1｝稱為邊移動(dòng)掃過的區(qū)域。由文獻(xiàn)［22］可知邊Fi，j＋移動(dòng)掃描區(qū)域上帶符號(hào)的面積｜δFi，j＋｜可通過如下線積分來表示：

式（6）的符號(hào)由定向四邊形δFi，j＋＝｛Pi，ji，j，i，j＋1，Pi，j＋1｝的頂點(diǎn)順序決定。由定向的選擇可知圖1中的情形為｜δFi，j＋｜＜0。

圖1 新舊單元征與掃描區(qū)域Fig.1 Old and new cells and swept regions

同樣可以定義任意多項(xiàng)式在多邊形上的符號(hào)積分。在本節(jié)下面描述中，線性函數(shù)在多邊形上積分都是指在這個(gè)多邊形上的符號(hào)積分。

這里f是邊數(shù)值通量。與文獻(xiàn)［22］不同，本文已知的是舊網(wǎng)格上分片一次多項(xiàng)式逼近，所以以fi，j＋為例，定義為：

其中，

這意味著在邊移動(dòng)掃過的區(qū)域中，面積符號(hào)決定了被積函數(shù)取為邊哪一側(cè)單元上的函數(shù)。

2.2 修補(bǔ)算法

具體計(jì)算時(shí)，一般要求重映算法具有保界性，也就是說：如果一個(gè)新單元包含在相應(yīng)舊單元的鄰域中，則有：

其中I（Ci）是C（Ci）中所有單元標(biāo)號(hào)的集合。文獻(xiàn)［23］給出了一種基于質(zhì)量守恒的、局部單元質(zhì)量重新分配的修補(bǔ)方法——局部調(diào)整不在范圍內(nèi)的值到給定界限內(nèi)。將這種方法用在本文重映過程（重映出單元平均值之后），下面給出該方法的簡潔描述。

首先，對(duì)每個(gè)單元Ci選擇一個(gè)決定它單元平均值范圍的鄰近單元集合R（Ci），R（Ci）包含單元Ci本身，那么密度平均的局部最小、最大值分別為：

其中ri是R（Ci）中所有單元標(biāo)號(hào)的集合。

因?yàn)橐罂傎|(zhì)量守恒，所以這個(gè)需要添加的質(zhì)量必須從其它單元獲得。首先檢查在不超過鄰近單元平均值最下界的情況下，可以從周圍取出的總質(zhì)量值。也就是說對(duì)任意k∈ri，可以安全從這個(gè)單元取走的質(zhì)量為：

那么總的可用質(zhì)量是

如果總的可用質(zhì)量足夠提供單元珟Ci的需求質(zhì)量，即：

那么單元珘Ci的質(zhì)量和密度平均值可以修補(bǔ)到最下界，

鄰近單元質(zhì)量適當(dāng)?shù)臏p少為：

并且單元珟Ci和它鄰近單元的總質(zhì)量不變。對(duì)于密度平均值：

更多修補(bǔ)細(xì)節(jié)參見［23］，并且修補(bǔ)算法中相鄰單元間的質(zhì)量轉(zhuǎn)移，與共有邊上的質(zhì)量通量無關(guān)，只是將臨近單元中可以安全取走的質(zhì)量直接加入當(dāng)前單元。

3 重構(gòu)分片一次多項(xiàng)式

3.1 極小化問題求解梯度

那么將式（20）代入式（21）可得一元二次方程組：

3.2 Van Leer限制器

在整個(gè)重映算法中，我們還希望新單元上分片一次多項(xiàng)式不出現(xiàn)新的極值。因此使用文獻(xiàn)［20］中Van Leer限制器對(duì)式（4）中梯度）進(jìn)行限制。下面以單元為例，給出限制器的描述。

至此整個(gè)重映過程完成，它把舊網(wǎng)格上分片線性數(shù)值解重映為新網(wǎng)格上分片線性數(shù)值解，并且滿足保界性和守恒性。

4 數(shù)值算例

將上述加修補(bǔ)和限制器的重映方法記為“SRRL”方法，將不加修補(bǔ)和限制器的重映方法記為“SRUN”方法。為了檢驗(yàn)SRRL方法的有效性，對(duì)給定函數(shù)在一系列網(wǎng)格上用SRRL方法做重映，并與SRUN方法的結(jié)果作比較。在下面誤差計(jì)算中L1模和L∞模定義為：

4.1 交替系列網(wǎng)格

初始網(wǎng)格是［0，1］×［0，1］上的一致均勻網(wǎng)格。下一套網(wǎng)格是通過對(duì)該一致均勻網(wǎng)格作一個(gè)小的隨機(jī)位移得到的隨機(jī)網(wǎng)格。交替系列網(wǎng)格就是由一致網(wǎng)格和相應(yīng)隨機(jī)網(wǎng)格交替產(chǎn)生的系列網(wǎng)格，即一致網(wǎng)格——隨機(jī)網(wǎng)格——再恢復(fù)到一致網(wǎng)格，這樣交替進(jìn)行的一系列網(wǎng)格。在數(shù)值算例中作偶數(shù)次網(wǎng)格交替和重映，使最終網(wǎng)格恢復(fù)為初始網(wǎng)格，計(jì)算出整個(gè)重映過程的累積誤差。當(dāng)細(xì)分網(wǎng)格時(shí)，增加重映的步數(shù)，從而使總的位移保持近似相同。

算例1。分別用SRRL、SRUN方法在交替網(wǎng)格上重映函數(shù)（x－0.5）4＋（y－0.5）4。選取三套不同水平的網(wǎng)格：imax＝j(luò)max＝50，100，200和相應(yīng)時(shí)間步數(shù)：nmax＝10，20，40來考查SRRL方法的收斂性。

表1和表2分別給出了SRUN和SRRL方法的數(shù)值結(jié)果。從表中可以看出，在誤差允許范圍內(nèi)，這兩種方法對(duì)于光滑函數(shù)的誤差基本相等。文獻(xiàn)［23］L1模誤差在剖分為100、200時(shí)的收斂階分別是1.862、1.927；L∞模收斂階分別是1.784、1.741。顯然本文SRRL方法收斂更快一些，更趨于2。

表1 SRUN方法的誤差Table 1 Errors of the SRUN method

表2 SRRL方法的誤差Table 2 Errors of the SRRL method

算例2。在50×50剖分的交替系列網(wǎng)格上用SRRL和SRUN方法重映一個(gè)三層塔式函數(shù)，函數(shù)取值范圍是［0，1］。令d（x，y）＝｜x－0.5｜＋｜y－0.5｜，那么這個(gè)函數(shù)表達(dá)式為：

表3給出了在最終網(wǎng)格上用SRRL和SRUN方法重映得到的函數(shù)最大值和最小值。從表3中可以看出，該結(jié)果與文獻(xiàn)［23］中數(shù)值結(jié)論一致：加過修補(bǔ)和Van Leer限制器的SRRL方法是保界的，并且SRRL方法L1模誤差比SRUN方法的誤差小。

表3 SRRL和SRUN方法的數(shù)值結(jié)果Table 3 Results of the SRRL and SRUN methods

4.2 隨機(jī)網(wǎng)格

本節(jié)使用一系列的隨機(jī)網(wǎng)格，其中每一套網(wǎng)格都是通過對(duì)一致均勻網(wǎng)格做一個(gè)隨機(jī)擾動(dòng)得到：

算例3。在此用SRRL方法和SRUN方法在隨機(jī)系列網(wǎng)格上重映激波函數(shù)：

圖2給出了SRUN和SRRL方法計(jì)算得到的激波函數(shù)曲面圖和等值線圖，可以看出SRUN方法的數(shù)值解在間斷處出現(xiàn)了數(shù)值振蕩，并且數(shù)值解超出了給定范圍，出現(xiàn)負(fù)值。而SRRL方法的數(shù)值解在間斷處沒有振蕩，并且數(shù)值解在給定界限內(nèi)。表4給出了在最終網(wǎng)格上兩種重映方法得到的函數(shù)最大和最小值，結(jié)果顯示SRRL方法是保界的，并且SRRL方法L1模誤差比SRUN方法的誤差小。文獻(xiàn)［22］中mc方法在γ＝0.7時(shí)出現(xiàn)負(fù)值，da和mb方法的數(shù)值解仍為正值。通過表4對(duì)比SRRL、da以及mb方法的結(jié)果可以發(fā)現(xiàn)，本文SRRL方法比da和mb方法的結(jié)果更精確。

表4 SRRL、SRUN方法以及文獻(xiàn)［22］中da和mb方法的誤差Table 4 Errors of the SRRL，SRUN methods and the da and mb methods in ref.［22］

圖2 隨機(jī)系列網(wǎng)格上重映激波函數(shù)的結(jié)果Fig.2 Results for the shock function on sequence of random grids

4.3 初始隨機(jī)網(wǎng)格的連續(xù)光滑化

在本次測試中，選取與ALE數(shù)值模擬方法相接近的一系列網(wǎng)格，其中每一組網(wǎng)格都是由前一組網(wǎng)格光滑化得到：

y方向的節(jié)點(diǎn)選取類似。初始網(wǎng)格是通過對(duì)一致均勻網(wǎng)格作隨機(jī)擾動(dòng)得到的：

這里－0.25≤ri，rj≤0.25是隨機(jī)數(shù)，h＝1／（imax－1）。邊界上節(jié)點(diǎn)只有一個(gè)坐標(biāo)被擾動(dòng)。

算例4。用SRRL方法和SRUN方法在連續(xù)光滑化網(wǎng)格上重映激波函數(shù)式（31）。圖3給出了SRUN和SRRL方法計(jì)算的激波函數(shù)側(cè)面圖，可以看出SRUN方法的數(shù)值解在間斷處出現(xiàn)了小振蕩，并且數(shù)值解超出給定范圍，而SRRL方法的數(shù)值解在間斷處沒有振蕩。表5給出了兩種方法的數(shù)值結(jié)果，結(jié)果顯示SRRL方法是保界的，并且SRRL方法L1模誤差比SRUN方法誤差小。比較算例3和算例4的結(jié)果，可以發(fā)現(xiàn)同樣的函數(shù)在連續(xù)光滑化網(wǎng)格上重映結(jié)果更好，這表明本文重映算法適合用于ALE方法。

圖3 連續(xù)光滑化網(wǎng)格上重映激波函數(shù)的側(cè)面圖Fig.3 Lateral views for the shock function on consecutive smoothing of an initially random grid

表5 SRRL和SRUN方法的數(shù)值結(jié)果Table 5 Results of the SRRL and SRUN methods

5 結(jié) 論

本文針對(duì)ALE方法中由間斷有限元法求解流體方程的拉氏步，給出了一種新的重映算法。算法由兩步組成：第一步在已有對(duì)流重映算法基礎(chǔ)上，計(jì)算出新網(wǎng)格上物理量的單元平均值，并通過修補(bǔ)算法，使之滿足保界性；第二步是通過重構(gòu)得到新單元上函數(shù)的近似梯度，并通過Van Leer限制器對(duì)梯度進(jìn)行限制，使之不出現(xiàn)新極值。這樣就得到了新網(wǎng)格上的分片一次多項(xiàng)式。最后通過一系列數(shù)值算例表明，加上修補(bǔ)和限制器的重映算法和不加修補(bǔ)和限制器的重映算法，對(duì)于光滑函數(shù)L1和L∞模誤差都是二階收斂，且誤差相差不大。但不加修補(bǔ)和限制器時(shí)，會(huì)出現(xiàn)數(shù)值結(jié)果越界的情形，加上修補(bǔ)和限制器之后，重映算法是保界的，并且誤差更小。本文重映算法可以直接用于Euler方程或NS方程，在重映物理量ρV和ρE時(shí)（其中V是速度，E是單位質(zhì)量的總能量），只是將文中的密度ρ替換為ρV和ρE，并且可以保證動(dòng)量和能量守恒，以及整個(gè)方程組的相容性。

［1］ Li Z Z，Yu X J，Jia Z P.The cell-centered discontinuous Galerkin method for Lagrangian compressible Euler equations in two dimensions［J］.Computers ＆Fluids，2014，96（C）：152-164.

［2］ Li Z Z，Yu X J，Zhu J，et al.A Runge Kutta discontinuous Galerkin method for Lagrangian compressible Euler equations in two-dimensions［J］.Communications in Computational Physics，2014，15：1184-1206.

［3］ Li Z Z，Yu X J，Zhao G Z，et al.A RKDG finite element method for Lagrangian Euler equations in one dimension［J］.Chinese Journal of Computational Physics，2014，31（1）：1-10.（in Chinese）李珍珍，蔚喜軍，趙國忠，等.RKDG有限元法求解一維拉格朗日形式的Euler方程［J］.計(jì)算物理，2014，31（1）：1-10.

［4］ Margolin L G.Introduction to“an arbitrary-Lagrangian-Eulerian computing method for all flow speeds”［J］.Journal of Computational Physics，1997，135：198-202.

［5］ Loubere R，Maire P H，Shashkov M，et al.Re ALE：a reconnection based arbitrary Lagrangian Eulerian method［J］.Journal of Computational Physics，2010，229：4724-4761.

［6］ Dukowicz J K，Kodis J W.Accurate conservative remapping（rezoning）for arbitrary Lagrangian-Eulerian computations［J］.Siam Journal on Scientific ＆Statistical Computing，1987，8：305-321.

［7］ Miller D S，Burton D E，Oliviera J S.Efficient second order remapping on arbitrary two dimensional meshes［R］.UCRL-ID-123530.Livermore，CA：Lawrence Livermore National Laboratory，1996.

［8］ Grandy J.Conservative remapping and region overlays by intersecting arbitrary polyhedra［J］.Journal of Computational Physics，1999，148：433.

［9］ Dukowicz J K，Padial N T.Accurate REMAP3D：a conservative three-dimensional remapping code［R］.LA-12136-MS.Los Alamos，NM：Los Alamos National Laboratory，1991.

［10］Mosso S J，Burton D E，Harrison A K，et al.A second order two and three-timensional temap method［R］.LA-UR-98-5353.Los Alamos，NM：Los Alamos National Laboratory，1998.

［11］Leslie L M，Purser R J.Three-dimensional mass-conserving semi-Lagrangian scheme employing forward trajectories［J］.Monthly Weather Review，1995，123：2551.

［12］Leonard B P，Lock A P，Macvean M K.Conservative explicit unrestricted-time-step multidimensional constancy preserving advection schemes［J］.Monthly Weather Review，1996，124：2588.

［13］Mosso S J，Swartz B.An unsplit，two-dimensional advection algorithm［R］.LA-UR-01-1476.Los Alamos，NM，USA：Los Alamos National Laborstory，2000.

［14］Orourke P J，Sahota M S.Avariable explicit／implicit numerical method for calculating advection on unstructured meshes［J］.Journal of Computational Physics，1998，143：312-345.

［15］Anderson R W，Pember R B，Elliott N S.An arbitrary Lagrangian-Eulerian method with local structured adaptive mesh refinement for modeling shock hydrodynamics［R］.AIAA 2002-0738.LNLL Report UCRL-JC-141625.

［16］Benson D J.An efficient，accurate，simple ALE method for nonlinear finite element programs［J］.Computer Methods in Applied Mechanics and Engineering，1989，72：305-350.

［17］Benson D J.Momentum advection on a staggered mesh［J］.Journal of Computational Physics，1992，100：143-162.

［18］Benson D J.Computational methods in Lagrangian and Eulerian hydrocodes［J］.Computer Methods in Applied Mechanics and Engineering，1992，99：235-394.

［19］Colella P.Multidimensional upwind methods for hyperbolic conservation laws［J］.Journal of Computational Physics，1990，87：171-200.

［20］Dukowicz J K，Baumgardner J.Incremental remapping as a transport／advection algorithm［J］.Journal of Computational Physics，2000，160：318-335.

［21］Margolin L G，Bason C W.Remapping on the staggered mesh［R］.UCRL-99682.Lawrence Livermore National Laboratory，1988.

［22］Margolin L G，Shashkov M.Second-order sign-preserving conservative interpolation（remapping）on general grids［J］.Journal of Computational Physics，2003，184（1）：266-298.

［23］Kucharik M，Shashkov M，Wendroff B.An efficient linearityand-bound-preserving remapping methods［J］.Journal of Computational Physics，2003，188：462-471.

A new second-order bound-preserving conservative remapping algorithm in the ALE method

Li Zhenzhen1，Yu Xijun2，＊，Jia Zupeng2，An Na3，Huang Chaobao3
（1.College of Mathematics and Information Science，Zhengzhou University of Light Industry，Zhengzhou 450002，China；2.National Key Laboratory of Computational Physics，Institute of Applied Physics and Computational Mathematics，Beijing 100088，China；3.Graduate School，China Academy of Engineering Physics，Beijing 100088，China）

When the Euler equations are solved using Lagrangian scheme，the fact that computational cells exactly follow fluid particles may result in severe grid deformation，even more cause inaccuracy and even breakdown of the computation in some cases.So it needs to rezone meshes and remap physical quantities when the deformation of computational grid is severe.Based on the second-order Lagrangian schemes that using the discontinuous Galerkin method to solve the Euler equations，a conservative remapping scheme is proposed.This remapping scheme has two steps：the first step is using the existing remapping method to obtain the approximate average values in the new cells，then using the repair algorithm to ensure the average values in the range of local bounds；the second step is using the average values to reconstruct the linear polynomial in the new cell，and using Van Leer limiter to limit the gradient of this linear polynomial to ensure no new extremum.Results of some numerical tests are presented and demonstrate that this remapping scheme is second-order accurate，conservative and bound-preserving.

conservative algorithm；remapping；ALE methods

V211.3；O242.2

：Adoi：10.7638／kqdlxxb-2014.0078

2014-08-07；

2015-03-20

國家自然科學(xué)基金（11571002）；中國工程物理研究院科學(xué)基金（2013A0202011，2015B0101021）；國防基礎(chǔ)科研計(jì)劃資助（B1520133015）；鄭州輕工業(yè)學(xué)院博士科研基金（2014BSJJ089）

李珍珍（1985-），河南濮陽，女，博士，研究方向：計(jì)算流體力學(xué).E-mail：lizhenzhenpuyang＠163.com

蔚喜軍＊（1959-），內(nèi)蒙古人，男，研究員，研究方向：計(jì)算流體力學(xué).E-mail：yuxj＠iapcm.ac.cn

李珍珍，蔚喜軍，賈祖朋，等.ALE方法中一種新的二階保界守恒重映算法［J］.空氣動(dòng)力學(xué)學(xué)報(bào)，2015，33（6）：765-771.

10.7638／kqdlxxb-2014.0078 Li Z Z，Yu X J，Jia Z P，et al.A new second-order bound-preserving conservative remapping algorithm in the ALE method［J］.Acta Aerodynamica Sinica，2015，33（6）：765-771.

0258-1825（2015）06-0765-08