一種基于特征關系式的預處理遠場邊界條件

李歡，陳江濤，馬明生，周乃春

中國空氣動力研究與發展中心 計算空氣動力研究所，綿陽 621000

一種基于特征關系式的預處理遠場邊界條件

李歡，陳江濤，馬明生*，周乃春

中國空氣動力研究與發展中心 計算空氣動力研究所，綿陽 621000

高速可壓縮流動的計算方法應用于不可壓低速問題會因為當地速度和當地聲速量級相差較大產生剛性問題并導致數值收斂困難。預處理方法引入預處理矩陣使傳統的可壓縮方法具備了求解不可壓或低速問題的能力。由于預處理方程改變了Navier-Stokes方程的特征值，原Navier-Stokes方程的遠場邊界條件不再適用。目前被廣泛使用的預處理方程的遠場邊界條件為簡化邊界條件。本文以預處理方程為基礎進行分析，推導了一種基于特征關系式的遠場邊界條件。通過典型狀態算例數值實驗，驗證了該邊界條件及預處理方法的有效性。結果顯示，對于定常/非定常低速不可壓流動，合理地設置遠場邊界條件，預處理方法都能夠提高計算的收斂性和精度；對于跨聲速流動，預處理方法和未加預處理方法兩者的計算效率和計算精度相當。本文的預處理遠場邊界條件與簡化預處理遠場邊界條件相比一是能夠進一步提高計算的收斂速度；二是能夠有效降低遠場邊界位置對預處理方法數值計算結果的負面影響。

預處理方法；Navier-Stokes方程；特征關系式；遠場邊界條件；數值收斂性

傳統的可壓縮方法在亞/跨/超聲速流場數值模擬計算中迅速發展并被廣泛運用，但無法有效求解低速不可壓流場，導致長期以來可壓縮流場與不可壓縮流場采用不同的計算體系[1]。但是對于一些不可壓低速流域和可壓縮特性高速流域相互耦合的多種速度尺度并存的復雜流動問題，如直升機懸停狀態[2-3]、大收縮比管道流[4-5]及低速大迎角流動[6-7]等，傳統的不可壓縮方法不能求解。而傳統的可壓縮方法由于穩定性限制和低速時近似因子分裂的離散誤差，導致數值求解收斂困難，計算效率和計算準確性降低，甚至完全錯誤[8]。預處理方法通過改變Navier-Stokes方程的時間導數項使方程的特征速度在低速保持同一數量級，消除數值剛性問題，使可壓縮方法能夠同時求解不可壓流場[9]。但由于預處理方法改變了方程的特征系統，預處理方程的無反射邊界條件復雜繁瑣，目前，國內外研究者[10-18]大都采用簡化遠場邊界條件。簡化的遠場邊界條件由于未考慮流場內信息的傳播方向和物理量的守恒關系，為反射遠場邊界條件。反射遠場邊界條件數值求解收斂過程會小幅震蕩，計算效率較低,只有在遠場邊界位置取得足夠遠才能夠得到精度高的數值解[19]。

無反射特征線(特征關系式)法通過計算和計算域內部節點同樣的物理守恒特征相容關系得到邊界處未知的物理條件，是Euler方程一種準確的邊界條件設定理論[20]。但目前只有少數文獻[21-24]對基于特征關系式確定預處理方程遠場邊界條件的方法開展研究。文獻[21-23]將邊界面上左右兩端的來流值與計算域內的值根據網格邊界面上的特征值進行Roe平均確定邊界值。這種方法相對于簡化邊界條件考慮了流場物理信息傳播的方向，但不能準確描述邊界值、計算域內的值和來流值的守恒關系；文獻[24]以一維預處理歐拉方程為基礎進行分析，將流場信息凍結在最外層網格，推導了一種物理量的守恒關系式，且考慮了流場物理信息傳播的方向。但當預處理方程退化為Navier-Stokes方程時，該預處理邊界條件并不能回歸到原邊界條件，同時作者也未具體分析該邊界條件的收斂性。本文參考文獻[24]的思想，在計算域的邊界上將預處理方程的對流項解耦成一系列對應于相應的物理信息傳播的特征波以求得邊界處的未知物理條件。應用Weiss-Smith[25]預處理矩陣，從可壓縮Navier-Stokes方程出發，推導基于雙時間步法[26]和LU-SGS[27]隱式迭代的定常與非定常預處理方程和一種基于特征關系式的預處理方程的遠場邊界條件，保證遠場邊界條件與方程形式的映射對應關系。同時研究預處理方法在定常/非定常流場中的應用和預處理邊界條件對計算收斂性和計算精度的影響。經過數值試驗驗證，本文基于特征關系式的預處理方程的遠場邊界條件能夠有效降低遠場邊界條件位置對計算精度的影響，進一步提高預處理方法的計算效率和流場適用范圍。

1 控制方程

1.1 Navier-Stokes方程

本文使用的流場解算器為自主研發的MFlow軟件[28]。守恒形式的非定常可壓縮Navier-Stokes方程的積分形式為

(1)

或者原始變量的非定常可壓縮Navier-Stokes方程的積分形式為

(2)

式中：Ω表示控制體的體積；?Ω表示控制體封閉面的面積；S為積分面積；Fc為無黏通量；Fv為黏性通量, 通量定義參考文獻[24]；W為守恒變量；Q為原始變量；M為W與Q的轉換矩陣。W、Q和M的形式為

(3)

W=[ρρuρvρwρE]T

Q=[puvwT]T

(4)

式中：u、v和w為速度在x、y和z方向的分量；cp、ρ、T、p、H和E分別為定壓比熱、密度、溫度、壓力、單位質量氣體的總焓和總能量；ρp、ρT分別為密度對壓力和溫度的偏導數。

1.2 預處理方程

Weiss-Smith預處理矩陣用參數Θ替換式(2)中轉換矩陣M中的參數ρp，得到預處理方程的形式為

(5)

式中：

(6)

(7)

式中：Ur、U∞、k和c分別為參考速度、遠場來流速度、截斷參數和當地聲速。

由于預處理方法改變了原有控制方程組的瞬時特性，對于非定常計算，引入雙時間步法進行時間推進求解。非定常預處理方程形式為

(8)

式中：t和τ分別為真實時間和虛擬時間。

2 方程離散與時間推進

非定常預處理方程左端第一項隱式離散為

(9)

式中：ΔW=Wn+1,m+1-Wn+1,m;Vol為控制單元體積;上角標m為偽時間層，上角標n為物理時間層,φ取不同的值代表不同的時間精度,φ取0為時間一階精度，φ取1/2為時間二階精度。

偽時間導數項顯式離散為

(10)

對流通量采用隱式離散：

(11)

黏性通量采用顯示離散：

(12)

式中：下角標i為單元體面序號。將式(9)～式(12)代入式(8)，控制方程經過數值離散后改寫為

(13)

式中：Res為殘差向量；A為雅克比矩陣，具體形式為

(14)

(15)

LU-SGS方法將隱式算子式(13)分解為

(D+L)D-1(D+U)ΔWn=Res

(16)

式中：因子L由嚴格下三角矩陣組成；U由嚴格上三角矩陣組成；D為對角矩陣。式(16)可以分兩步求逆——前掃描和后掃描，即

(17)

(18)

(19)

(20)

由表1可知，預處理方程和Navier-Stokes方程的數值離散格式僅在矩陣J和方程特征系統的特征值λk處存在差異，而當ε=1或ε′=1時，各方程特征系統的特征值λk相同，矩陣J的形式相似，預處理方程轉化為Navier-Stokes方程。因此定常/非定常預處理方程和Navier-Stokes方程在此種離散推進格式的數據結構上很相似，預處理方法易于在原Navier-Stokes方程程序上編寫調試。

表1 數值離散格式對比Table 1 Comparison of numerical discretization schemes

3 遠場邊界條件

Weiss-Smith預處理矩陣的引入改變了原Navier-Stokes方程的的特征系統，因此原Navier-Stokes方程基于特征關系式的黎曼不變量遠場邊界條件需要進行預處理修正。大多數文獻采用簡化邊界條件(Precodtioned Equation and Simplified farfield Boundary condition，PESB)。PESB采用Turkel等[29]建議的簡化邊界條件，亞聲速入流邊界速度和溫度采用自由來流賦值，壓強采用流場計算域內插值；亞聲速出流邊界速度和密度采用流場計算域內插值，壓強采用流場外賦值；超聲速邊界所有的原始變量則根據流動方向從自由來流或者流場中賦值。

本文從一維定常預處理方程出發，推導一種基于特征關系式的遠場邊界條件(Precodtioned Equation and Preconditioned characteristic variable Boundary condition，PEPB)，并類推到三維定常預處理方程和三維非定常預處理方程。

式(5)忽略黏性項的微分形式為

(21)

式(21)可進一步轉化為

(22)

(23)

式中：Λ為矩陣Γ-1AM的特征值組成的對角矩陣；P為矩陣Γ-1AM的左特征向量組成的矩陣。式(23)為定常預處理方程的特征相容關系式，式(23)展開形式為

(24)

(25a)

(25b)

(26a)

(26b)

式中：

(27a)

A2=ρ(λ4-εu),B2=ρ(λ5-εu)

(27b)

A2A1=1,B2B1=1

(27c)

由式(27c)可知：式(25a)和式(25b)、式(26a)和(26b)本質上具有相同的形式。文獻[24]直接將A2、B2和C1凍結在計算域最外層網格單元I，得到預處理方程的3個不變量。

本文將特征相容關系式(24)進一步化簡為D(p/ργ)/Dt=0,即沿特征線dx/dt=λ1，p/ργ=const1。式(25a)和式 (26a)形式復雜，直接化簡不易，這里將A1和B1進一步分解為

A1=A3/ρc,B1=B3/ρc

A3=(λ4-u)/εc,B3=(λ5-u)/εc

(28)

將式(28)中的A3和B3凍結在計算域最外層網格單元I，特征相容關系式(25a)和式(26a)可以進一步化簡為

(29)

(30)

式中：

A3,I=(λ4,I-uI)/(εIcI)

B3,I=(λ5,I-uI)/(εIcI)

(31)

式中：λ4,I、λ5,I、εI、cI和uI均為網格單元I的體心值。

由特征相容關系式(29)和式(30)，可以得到預處理方程的2個不變量，本文的3個不變量形式為

(32)

預處理方程特征關系式得到了與Euler方程的Riemann不變量(Riemann Invariants，RI)類似的3個不變量，定義為類黎曼不變量(Similar Riemann Invariants，SRI)。SRI、RI和文獻[24]的不變量對比見表2。

表2 不變量對比Table 2 Comparison of invariants

表2中，A0、B0和C0的定義為

(33)

在ε=1時,Γ=M，預處理方程回歸為Euler方程。A3,I=1，B3,I=1，SRI也與RI完全相同，本文的預處理遠場邊界條件回歸到Euler方程的遠場邊界條件，保證了遠場邊界條件與方程形式的映射對應關系。而文獻[24]中的不變量卻與RI不相同，文獻[24]中的預處理遠場邊界條件不能回歸到Euler方程的遠場邊界條件。對于大收縮比噴管流動，大部分區域為不可壓縮流動，但在入流或出流存在超聲速流動，文獻[24]的邊界條件是不合適的；在ε≠1時，式(32)雖然是在一維定常預處理方程基礎上推導得出，但很容易推廣到三維定常/非定常預處理方程。由于預處理方程未改變原方程特征值的符號，因此預處理方程基于SRI推導PEPB邊界條件與Euler方程采用基于Riemann不變量推導無反射邊界條件(No precodtioned Equation and characteristic variable boundary condition，NE)的過程完全類似，具體參考文獻[30]。

4 考核算例

4.1 無黏圓柱繞流

此算例考察預處理方法兩種邊界條件和未加預處理方法計算效率的差別。運用預處理方程和歐拉方程對典型算例無黏圓柱繞流進行計算，網格分布為197×61×2(周向×法向×展向)的同心圓網格。預處理方程采用兩種遠場邊界條件：PESB和PEPB；歐拉方程采用黎曼無反射邊界條件：NE。分別計算來流馬赫數Ma=0.001,0.01,0.1,0.3這4個工作狀態。

圖1為殘差收斂曲線和升力系數收斂曲線歷程。其中升力系數收斂曲線歷程圖中NE方法的橫坐標為圖中上下橫坐標區間，圖上方橫坐標為圖下方橫坐標迭代步數的延續；PEPB和PESB的橫坐標為下橫坐標。圖1顯示，對于Ma=0.1,0.01,0.001這3種工況，不管采用何種遠場邊界條件，預處理方法都能保證數值解的精度并顯著改善數值解的收斂性，升力系數在15 000步已收斂至0，殘值在15 000步下降了8個量級。PEPB的殘值收斂歷程基本無震蕩，PESB的殘值收斂歷程震蕩劇烈；不加預處理的可壓縮方法(NE)不僅收斂慢，且升力系數計算結果不正確。對于工況Ma=0.1，升力系數在150 000步才能收斂至0，計算時間約為預處理方法的10倍。而對于Ma=0.01和Ma=0.001這兩種工況，升力系數在150 000步和300 000步均未收斂至0，計算結果不正確。對于工況Ma=0.3，殘值收斂歷程預處理方法的PEPB的計算收斂速度最快，殘值在15 000步依然下降了8個量級，而預處理方法的PESB的殘值收斂性反而不如未加預處理方法的NE的殘值收斂性。升力系數最終都能收斂到0，但NE在25 000步收斂至0，計算時間約為預處理方法的1.6倍。無黏圓柱繞流算例既考核了預處理方法高效率高精度求解無黏低速不可壓流場的能力,又驗證了預處理方法的PEPB能夠進一步提高預處理方法的計算效率。

圖1 殘差收斂曲線和升力系數收斂曲線歷程
Fig.1 Convergence history of residual and lift coefficients curves

4.2 Re=40的圓柱繞流

此算例考核預處理方法求解黏性流場的能力和遠場邊界位置對計算結果的影響。來流馬赫數為0.02，Re=40的圓柱繞流尾部上下兩側會形成穩定的渦，此狀態為定常狀態。預處理方法仍然采用兩種遠場邊界條件，Navier-Stokes方程為基于黎曼不變量的遠場邊界條件。計算網格采用同心圓網格，網格分布為113×169×2(周向×法向×展向)；壁面最小距離為0.002D(D為圓柱直徑)，計算域的遠場邊界位置取30D(Farfield_30D)和10D(Farfield_10D)。

圖2為采用PEPB計算的圓柱繞流的流線圖，遠場邊界位置取30D(Farfield_30D)。預處理方法能夠準確計算出圓柱繞流尾部上下兩側形成的穩定的渦結構。

圖3為計算的壁面壓力系數和實驗數據[31]的對比圖。表3為該流動狀態的阻力系數CD和前駐點壓力系數Cpf與文獻計算結果的對比。由圖3和表3可以清晰地看到PEPB計算的CD和Cpf與參考文獻符合很好，不受遠場邊界位置影響；PESB計算結果受遠場邊界位置影響較大，在遠場邊界位置Farfield_30D時與參考文獻符合較好；在遠場邊界位置Farfield_10D，計算結果與參考文獻偏差較大，CD偏差約為6%，Cpf偏差約為11.5%。而不加預處理的Navier-Stokes方程計算結果與參考文獻相比都存在較大偏差。

來流馬赫數為0.02、Re=40的圓柱繞流算例考核了預處理方法模擬低速黏性流場的能力，只要遠場邊界位置選取合適，相比未加預處理方法，預處理方法計算結果更準確。如果遠場邊界位置取得不是足夠遠，采用預處理方法的PESB會對計算結果造成一定程度的偏差，但采用PEPB能夠有效降低遠場邊界位置對計算結果的影響。

圖2 圓柱繞流的流線圖
Fig.2 Streamline diagram of flow around cylinder

圖3 壁面壓力系數對比
Fig.3 Comparison of wall pressure coefficients

=
=

計算方法CDCpfPESB＿Farfield＿30D1．5321．165PESB＿Farfield＿10D1．6271．275PEPB＿Farfield＿30D1．5111．143PEPB＿Farfield＿10D1．5361．190NE＿Farfield＿30D1．5901．180NE＿Farfield＿10D1．6711．275文獻[24]1．5171．123文獻[32]1．5491．147

4.3 Re=200的圓柱繞流

此算例考核預處理方法的非定常計算能力。來流馬赫數為0.02、Re=200的圓柱繞流尾部會形成周期性渦脫落，此狀態為非定常狀態。預處理方法仍然采用兩種遠場邊界條件，Navier-Stokes方程為基于黎曼不變量的遠場邊界條件。計算網格采用同心圓網格，網格分布為113×179×2(周向×法向×展向)；壁面最小距離為0.002D，計算域的遠場邊界位置取60D(非定常低速流動遠場邊界位置對結果影響較大)。真實時間步長Δt=0.002 5 s。

圖4為預處理方法的PEPB計算的某一時刻馬赫數等值線，預處理方法的PEPB能夠捕獲明顯的渦脫落現象。

圖5和圖6為計算的升力系數和阻力系數歷程(圖中時間步為真實時間步)。對于升力系數歷程，預處理方法的PEPB和PESB數值模擬結果基本一致，在2 000步流場達到周期諧振解；未加預處理方法在3 000步流場達到周期諧振解。對于阻力系數歷程，預處理方法的PESB和未加預處理方法的NE的模擬結果的振幅不斷變化，精度不高；而預處理方法的PEPB的模擬結果的周期和振幅保持不變，具有更高的計算精度。升力系數和阻力系數歷程說明PESB和NE能夠定性求解低速非定常流場，但時間精度和空間精度不高。表4為PEPB具體的數值結果與文獻的對比，St表示斯特勞哈爾數(St=nD/U∞，n為諧振頻率),St考核非定常方法的時間精度。預處理方法的PEPB的數值模擬結果和參考文獻吻合得很好，說明了本文發展的預處理方法具有求解非定常流場的能力，時間精度和空間精度均較高。

圖4 瞬時馬赫數等值線云圖
Fig.4 Contours of instantaneous Mach number

圖5 升力系數歷程
Fig.5 History of lift coefficients

圖6 阻力系數歷程
Fig.6 History of drag coefficients

=
=

計算方法CLCDStPEPB±0．681．30±0．0500．188文獻[33]±0．691．31±0．0490．192文獻[34]±0．661．32±0．0400．192

4.4 ONERA-M6機翼跨聲速繞流

此算例考核預處理方法非結構網格求解跨聲速可壓縮流場的能力以及預處理方法相比未加預處理方法對跨聲速可壓縮流場計算結果的影響。計算的湍流模型采用S-A模型[35]，計算狀態為

(34)

式中：雷諾數基于機翼的平均氣動弦長。在該計算狀態下，機翼表面為附著流，機翼上表面有λ型激波[36]。上表面網格分布與空間網格分布如圖7所示，在激波流場特征處未有明顯的加密網格。

圖8為ONERA-M6機翼壁面壓力系數分布對比，x/L為相對弦長，y/b為相對展長。預處理方法和未加預處理方法計算結果保持一致，說明預處理方法相比未加預處理方法對跨聲速流場計算結果無影響；對于機翼下表面，數值計算結果與實驗數據[37]吻合得很好；雖然上表面存在λ型激波，且計算網格分布未在激波處加密，但上表面壓力系數數值計算結果與實驗值仍然吻合得較好，驗證了預處理方法具有求解跨聲速流場的能力。

圖7 ONERA-M6機翼網格
Fig.7 Grid of ONERA-M6 wing

圖8 ONERA-M6機翼壁面壓力系數分布對比
Fig.8 Comparison of wall surface pressure coefficients distribution of ONERA-M6 wing

5 結 論

1) 本文發展的定常/非定常預處理方程和Navier-Stokes方程在離散推進格式的數據結構上相似，預處理方法易于在原Navier-Stokes方程程序上編寫調試。

2) 本文發展的基于特征關系式的預處理遠場邊界條件PEPB的基礎是SRI(類黎曼不變量)，PEPB考慮了特征波的傳播方向和守恒不變量；當預處理方程回歸為原方程，SRI回歸到RI(黎曼不變量)，保證了遠場邊界條件與方程形式的映射對應關系。

3) 對于無黏低速流場，預處理方法的求解效率高于未加預處理方法。計算的馬赫數越低，未加預處理方法的計算時間成本越高。計算的馬赫數越高，本文發展的PEPB的收斂效率越優于PESB。

4) 對于低雷諾數黏性低速流場，遠場位置較近，PEPB比PESB的計算精度高；遠場位置足夠遠，PEPB和PESB的計算精度保持一致。PEPB能有效降低遠場邊界位置對模擬結果的負面影響，節省數值計算的網格成本。

5) 對于低雷諾數黏性非定常低速流場， PEPB高于PESB和未加預處理方法的計算精度與計算效率。對于跨聲速流場，預處理方法與未加預處理方法兩者的求解精度相當，預處理方法對跨聲速流場計算結果無影響。

本文發展了預處理方法，并在忽略黏性的基礎上發展了基于特征關系式的預處理邊界條件，獲得類黎曼不變量。提高了預處理方法無黏低速流場的計算效率，有效降低了遠場邊界位置對預處理方法計算結果的負面影響。對于跨聲速流動，本文發展的預處理方法具有與未加預處理方法相當的計算效率和計算精度。下一步工作一是開展預處理矩陣研究，提高預處理方法在跨聲速甚至全速域流場的求解效率；二是考慮黏性的影響，進一步開展預處理方法的邊界條件研究。

[1] 肖中云， 江雄， 陳作斌, 等. 低速預處理方法在懸停旋翼模擬中的應用[J]. 航空學報, 2008, 29(2): 321-326.

XIAO Z Y, JIANG X, CHEN Z B, et al. Application of low speed preconditioning to hovering rotor simulation[J]. Acta Aeronautica et Astronautica Sinica, 2008, 29(2): 321-326 (in Chinese).

[2] POTSDAM M, SANKARAN V, PANDYA S. Unsteady low mach preconditioning with application to rotorcraft flows: AIAA-2007-4473[R]. Reston, VA: AIAA, 2007.

[3] XIAO W, SHENG C H, CHEN J P. Numerical study of preconditioned algorithm for rotational flows: AIAA-2005-5126[R]. Reston, VA: AIAA, 2005.

[4] VENKATESWARAN S, WEISS J M, MERKLE C L. Propulsion-related flowfields using the preconditioned Navier-Stokes equations: AIAA-1992-3437[R]. Reston, VA: AIAA, 1992.

[5] TOURRETTE L. Navier-Stokes simulations of air-intakes in crosswind using local preconditioning: AIAA-2002-2739[R]. Reston, VA: AIAA, 2002.

[6] 陳龍, 伍貽兆, 夏健. 基于非定常低速預處理和DES的三角翼數值模擬[J]. 南京航空航天大學學報, 2011, 43(2): 159-164.

CHEN L, WU Y Z, XIA J. Delta wing numerical simulation using unsteady low speed precondition and DES[J]. Journal of Nanjing University of Aeronautics & Astronautics, 2011, 43(2): 159-164 (in Chinese).

[7] CHOI Y H, MRKLE C L. The application of preconditioning in viscous flows[J]. Journal of Computational Physics, 1993, 105(2): 207-223.

[8] VENKATESWARAN S, MERKLE C L. Dual time stepping and preconditioning for unsteady computations: AIAA-1995-0078[R]. Reston, VA: AIAA, 1995.

[9] VENKATAKRISHNAN V. Improved convergence of compressible Navier-Stokes solvers: AIAA-1998-2967[R]. Reston, VA: AIAA, 1998.

[10] 劉中玉， 張明鋒， 鄭冠男, 等. 基于預處理HLLEW格式的全速域數值算法[J]. 計算物理, 2016, 33(3): 273-282.

LIU Z Y, ZHANG M F，ZHENG G N, et al. Preconditioning HLLEW scheme for flows at all Mach numbers[J]. Chinese Journal of Computational Physics, 2016, 33(3): 273-282 (in Chinese).

[11] 白曉征, 郭正, 劉君. 黏性極低速流動的預處理方法研究[J]. 空氣動力學學報, 2009, 27(4): 451-455.

BAI X Z, GUO Z, LIU J. Study of preconditioning methods for very low speed viscous flows[J]. Acta Aerodynamica Sinca, 2009, 27(4): 451-455 (in Chinese).

[12] 劉晨, 王江峰, 伍貽兆. 低馬赫數流動中的預處理Euler方程的收斂特性[J]. 航空學報, 2009, 30(5): 842-848.

LIU C, WANG J F, WU Y Z. Convergence characteristics of preconditioned Euler equations at low Mach numbers[J]. Acta Aeronautica et Astronautica Sinica, 2009, 30(5): 842-848 (in Chinese).

[13] 韓忠華, 宋文萍, 喬志德. 一種隱式預處理方法及其在定常和非定常流動數值模擬中的應用[J]. 計算物理, 2009, 26(5): 679-684.

HAN Z H, SONG W P, QIAO Z D. An implicit preconditioning method for steady and unsteady flows[J]. Chinese Journal of Computational Physics, 2009, 26(5): 679-684 (in Chinese).

[14] LIU Z N, RAMAKRISHNAN S V. Preconditioned upwind schemes for simulation of low speed flows on hybrid grids: AIAA-2004-1268[R]. Reston, VA: AIAA, 2004.

[15] PANDYA S A, VENKATESWARAN S. Implementation of preconditioned dual-time procedures in overflow: AIAA-2003-0072[R]. Reston, VA: AIAA, 2003.

[16] TURKEL E. Preconditioning techniques in computational fluid dynamics[J]. Annual Review of Fluid Mechanics, 1999, 31: 385-416.

[17] HAUKE G, HUGHES T J R. A comparative study of different set of variables for solving compressible and incompressible flows[J]. Computer Methods in Applied Mechanics and Engineering, 1998, 153(1): 1-44.

[18] MERKLE C L. Preconditioning methods for viscous flow calculations[J]. Computational Fluid Dynamics Review, 1995, 153: 419-436.

[19] DARMOFAL D L. Eigenmode analysis of boundary conditions for the one-dimensional preconditioned Euler equations[M]. Cambridge, Massachusetts: Academic Press Professional, Inc., 1998.

[20] THOMPSON K W. Time-dependent boundary conditions for hyperbolic systems[M]. Cambridge, Massachusetts: Academic Press Professional, Inc., 1987.

[21] BARTH T J. Aspects of unstructured grids and finite-volume solvers for the Euler and Navier-Stokes equations[R]. Moffett Field, California: Lecture notes for the Von Karman Institute for Fluid Dy-namics Lecture Series,1994.

[22] UNRAU D，ZINGG D W. Viscous airfoil computations using local preconditioning: AIAA-1997-2027[R]. Reston, VA: AIAA, 1997.

[23] 張強， 楊永. 迎風格式的低速預處理及遠場邊界影響研究[J]. 西北工業大學學報, 2012, 30(3): 412-416.

ZHANG Q， YANG Y. A better preconditioning method of low Mach number flow with effects of farfield boundary conditions considered [J]. Journal of Northwestern Polytechnical University, 2012, 30(3): 412-416 (in Chinese).

[24] 肖天航, 昂海松, 余少志, 等. 預處理法求解定常/非定常混合網格的全速流場[J]. 空氣動力學學報, 2007, 25(4): 425-430.

XIAO T H, ANG H S, YU S Z, et al. Using preconditioning method to simulate steady/unsteady flows at all speeds on hybrid grids[J]. Acta Aerodynamica Sinca, 2007, 25(4): 425-430 (in Chinese).

[25] WEISS J M, SMITH W A. Preconditioning applied to variable and constant density flows: AIAA-1995-0033[R]. Reston, VA: AIAA, 1995.

[26] JAMESON A. Time depended calculations using multigrid with applications to unsteady flows past airfoils and wings: AIAA-1991-1596[R]. Reston, VA: AIAA, 1991.

[27] KIM J S, KWON O J. Improvement on block LU-SGS scheme for unstructured mesh Navier-Stokes computa-tions: AIAA-2002-1061[R]. Reston, VA: AIAA, 2002.

[28] 張耀冰, 周乃春, 陳江濤. 小展弦比飛翼標模雷諾數影響數值模擬研究[J]. 空氣動力學學報, 2015, 33(3): 279-288.

ZHANG Y B, ZHOU N C, CHENG J T. Numerical investigation of Reynolds number effects on a low-aspect-ratio flying-wing model[J]. Acta Aerodynamica Sinica, 2015, 33(3): 279-288 (in Chinese).

[29] TURKEL E, VATSA V N, RADESPIEL R. Preconditioning methods for low-speed flows: AIAA-1996-2460 [R]. Reston, VA: AIAA, 1996.

[30] 廖守億, 王承堯. 預處理方法求解低馬赫數流場的邊界條件[J]. 國防科技大學學報, 2001, 23(2): 52-56.

LIAO S Y, WANG C Y. Boundary conditions for solving flowfields at low Mach numbers by preconditioning[J]. Journal of National University of Defense Technology, 2001, 23(2): 52-56 (in Chinese).

[31] GROVE A S, SHAIR F H, PETERSEN E E, et al. An experiment investigation of the steady separated flow past a circle cylinder[J]. Joural of Fluid Mechanics, 2006, 19(1): 60-80.

[32] ROGERS S E, KWAK D. An upwind differencing scheme for the incompressible Navier-Stokes equations: NASA-1988-101051[R]. Washington, D.C.: NASA, 1988.

[33] LIU C, ZHENG X, LIAO C. Preconditioned multigrid methods for unsteady incompressible flows: AIAA-1997-0445[R]. Reston, VA: AIAA, 1997.

[34] 牟斌, 肖中云, 劉剛, 等. 基于Roe格式的低速預處理方法研究[J]. 空氣動力學學報, 2010, 28(2): 125-131.

MOU B, XIAO Z Y, LIU G, et al. Low-speed preconditioning for Roe scheme[J]. Acta Aerodynamica Sinica, 2010, 28(2): 125-131 (in Chinese).

[35] SPALART P, ALLMARAS S. A one-equation turbulence model for aerodynamic flows : AIAA-1992-0439[R]. Reston, VA: AIAA, 1992.

[36] 崔鵬程, 鄧有奇, 唐靜, 等. 基于伴隨方程的網格自適應及誤差修正[J]. 航空學報, 2016, 37(10): 2992-3002.

CUI P C, DENG Y Q, TANG J, et al. Adjoint-based grid adaptation and error correction[J]. Acta Aeronautica et Astronautica Sinica, 2016, 37(10): 2992-3002 (in Chinese).

[37] SCHMITT V, CHARPIN F. Pressure distributions on the ONERA M6-wing at transonic mach number[R]. Paris: AGARD, 1979.

Afarfieldboundaryconditionforpreconditioningmethodbasedoncharacteristicrelations

LIHuan,CHENJiangtao,MAMingsheng*,ZHOUNaichun

ComputationalAerodynamicsInstitute,ChinaAerodynamicsResearchandDevelopmentCenter,Mianyang621000,China

Duetolargedifferencebetweenlocalflowspeedandlocalsoundspeed,thenumericalmethodoriginallydevisedforcomputationofhigh-speedcompressibleflowswillsufferfromnumericalstiffproblemstoslowdowntheconvergencerate,whenthemethodisextendedtocomputationoflow-speedincompressibleflows.Afterthepreconditioningmatrixisintroduced,theoriginalmethodwillbecapableofhandlinglow-speedincompressibleproblems.However,thewidely-usedfarfieldboundaryconditionofNavier-Stokesequationsisnolongerappropriateforpreconditioningsystemssincetheeigenvaluesoftheequationsarealtered.Currently,thefarfieldboundaryconditionappliedtopreconditioningequationsismuchsimplified.Inthispaper,animprovedfarfieldboundaryconditionisproposed,whichisdevisedbasedoncharacteristicrelations.Thevalidityoftheproposedboundaryconditionisdemonstratedbyseveraltypicalcases.Itisprovedthatthepreconditioningmethodwithappropriatefarfieldboundaryconditionswillimprovetheconvergencerateandaccuracyoflow-speedincompressibleflowcomputations.Whenthepreconditioningmethodisappliedtotransonicflowproblems,theefficiencyandaccuracyofthecomputationsreachthesamelevelasthetraditionalcompressiblemethod.Theadvantagesoftheproposedfarfieldboundaryconditionoverthesimplifiedconditionlieintwoaspectstheconvergencerateisaccelerated,andthenegativeeffectofinsufficientfarfieldextentisreduced.

preconditioningmethod;Navier-Stokesequations;characteristicrelations;farfieldboundarycondition;numericalconvergencerate

2017-04-27；

2017-05-23；

2017-06-26；Publishedonline2017-07-041701

URL：http://hkxb.buaa.edu.cn/CN/html/20171206.html

.E-mailma_mingsheng@sina.cn

http://hkxb.buaa.edu.cnhkxb@buaa.edu.cn

10.7527/S1000-6893.2017.121364

2017-04-27；退修日期2017-05-23；錄用日期2017-06-26；網絡出版時間2017-07-041701

http://hkxb.buaa.edu.cn/CN/html/20171206.html

.E-mailma_mingsheng@sina.cn

李歡，陳江濤，馬明生,等．一種基于特征關系式的預處理遠場邊界條件J． 航空學報,2017,38(12):121364.LIH,CHENJT,MAMS,etal.AfarfieldboundaryconditionforpreconditioningmethodbasedoncharacteristicrelationsJ.ActaAeronauticaetAstronauticaSinica,2017,38(12):121364.

V211.3

A

1000-6893(2017)12-121364-13

李明敏)