時(shí)間譜方法中的高效GMRES算法

貢伊明,劉戰(zhàn)合,劉溢浪,張偉偉,*

1.西北工業(yè)大學(xué) 航空學(xué)院,西安 710072

2.鄭州航空工業(yè)管理學(xué)院航空工程學(xué)院,鄭州450046

時(shí)間譜方法中的高效GMRES算法

貢伊明1,劉戰(zhàn)合2,劉溢浪1,張偉偉1,*

1.西北工業(yè)大學(xué) 航空學(xué)院,西安 710072

2.鄭州航空工業(yè)管理學(xué)院航空工程學(xué)院,鄭州450046

研究了時(shí)間譜方法求解周期性非定常流場(chǎng)的計(jì)算效率,并對(duì)時(shí)間譜方法應(yīng)用于周期性非定常流動(dòng)的隱式求解方法進(jìn)行探討。當(dāng)采樣點(diǎn)數(shù)增加或減縮頻率增大時(shí),時(shí)間譜方法對(duì)應(yīng)的雅可比矩陣對(duì)角占優(yōu)性質(zhì)迅速惡化,導(dǎo)致很多傳統(tǒng)的迭代方法失效。為了解決上述問(wèn)題,論文采用帶預(yù)處理的廣義極小殘差(GMRES)算法來(lái)提高雅可比系數(shù)矩陣的計(jì)算收斂性。使用時(shí)間譜方法對(duì)NACA0012翼型強(qiáng)迫振蕩算例進(jìn)行計(jì)算,并與時(shí)域差分方法的計(jì)算效率和精度進(jìn)行對(duì)比。研究表明在保證計(jì)算精度的同時(shí),時(shí)間譜方法普遍可將計(jì)算效率提高一個(gè)量級(jí)左右。對(duì)于跨聲速周期性流動(dòng),廣義極小殘差算法不論是穩(wěn)定性還是收斂性都優(yōu)于對(duì)稱SGS迭代算法。

時(shí)間譜方法;廣義極小值殘差(GMRES)算法;周期性非定常流動(dòng);預(yù)處理;計(jì)算效率

近些年來(lái),隨著計(jì)算機(jī)性能的提高,計(jì)算流體力學(xué)(CFD)在工程中得到廣泛的運(yùn)用。在實(shí)際的工程應(yīng)用中,很多問(wèn)題涉及到的非定常流動(dòng)都具有周期性特征。比如繞直升機(jī)旋翼的非定常流動(dòng)和渦輪機(jī)葉片繞流等。在傳統(tǒng)周期性非定常流動(dòng)的計(jì)算方法中,時(shí)間項(xiàng)均采用時(shí)域差分法(Backward Difference Format,BDF),為保證計(jì)算精確性,一個(gè)周期往往被分為幾十個(gè)時(shí)刻點(diǎn)并且要經(jīng)過(guò)若干周期才能得到收斂的周期性非定常解。因此BDF在解決周期性非定常流動(dòng)問(wèn)題上計(jì)算效率偏低,沒(méi)有利用流動(dòng)的周期性特征。

對(duì)于周期性非定常流動(dòng),Hall和Crawley[1]提出了一種線化頻域快速計(jì)算方法。在該方法中,時(shí)域流場(chǎng)變量分解為時(shí)均量與非定常一階諧波項(xiàng)。其計(jì)算效率很高,但只局限于非定常性較弱和非線性較低的非定常計(jì)算中。1998年Ning和He[2]提出了一種非線性計(jì)算方法,不同之處在于其通過(guò)額外的應(yīng)力建立時(shí)均部分與非定常項(xiàng)之間的關(guān)系,通過(guò)這種方法能夠求解某些比較復(fù)雜和變化劇烈的非定常流動(dòng),但該方法仍然無(wú)法求解非簡(jiǎn)諧周期性流動(dòng)。隨后Hall等[3]提出了基于傅里葉展開(kāi)的諧波平衡法。該方法不直接求解流場(chǎng)變量的諧波系數(shù),其在時(shí)域內(nèi)通過(guò)譜算子將方程進(jìn)行時(shí)間方向的離散。Hall等將該方法成功應(yīng)用于模擬壓氣機(jī)葉柵繞流。Mc Mullen等[4-5]發(fā)展了一種主要在頻域里求解流動(dòng)控制方程的非線性頻域方法。該方法與諧波平衡法的思想一致只是構(gòu)造方法有所不同且成功應(yīng)用到模擬圓柱渦脫落運(yùn)動(dòng)和模擬翼型強(qiáng)迫俯仰振蕩等問(wèn)題。但該方法需要在每步迭代中計(jì)算守恒變量和殘值的傅里葉變換與逆變換。類似的還有Dai等[6]提出的時(shí)間配置法。Gopinath等[7-8]推導(dǎo)了離散時(shí)間導(dǎo)數(shù)譜算子的顯式公式并進(jìn)一步提出時(shí)間譜方法(Time Spectral Method,TSM)。該方法利用流動(dòng)周期性特征,通過(guò)傅里葉變換與逆變換將時(shí)間導(dǎo)數(shù)項(xiàng)通過(guò)其他所有時(shí)刻點(diǎn)的插值進(jìn)行離散,從而將非定常問(wèn)題中各時(shí)刻的順序排列求解轉(zhuǎn)化為多個(gè)時(shí)刻耦合的定常求解問(wèn)題。TSM在工程中得到廣泛應(yīng)用,并取得了很好的實(shí)際效果,如直升機(jī)旋翼問(wèn)題[9]、渦輪機(jī)葉片繞流[10]以及渦脫落[11]問(wèn)題。國(guó)內(nèi)對(duì)TSM研究相對(duì)較少,比較典型的是楊小權(quán)等[12]應(yīng)用基于隱式LUSGS(Lower-Upper Symmetric Gauss-Seidel)迭代的TSM模擬俯仰翼型和機(jī)翼的強(qiáng)迫運(yùn)動(dòng)問(wèn)題,謝立軍等[13]應(yīng)用基于SGS(Symmetric Gauss-Seidel)迭代的TSM計(jì)算動(dòng)導(dǎo)數(shù)問(wèn)題以及詹磊和劉鋒[14]對(duì)基于龍格庫(kù)塔算法的TSM計(jì)算精度進(jìn)行研究。

2009年Jameson和Shankaran[15]研究發(fā)現(xiàn)基于龍格庫(kù)塔算法的TSM計(jì)算撲翼時(shí)會(huì)出現(xiàn)不收斂的情況,計(jì)算效率不如BDF。另有研究表明,對(duì)于非線性強(qiáng)或非簡(jiǎn)諧的周期非定常流動(dòng)問(wèn)題,TSM則需要大量的采樣點(diǎn)才能擬合,TSM計(jì)算效率也大幅下降[16]。上述2種情況的出現(xiàn)是因?yàn)殡S著減縮頻率增大和采樣點(diǎn)數(shù)增加,雅可比矩陣對(duì)角占優(yōu)性質(zhì)惡化,使得計(jì)算效率下降甚至數(shù)值發(fā)散。

廣義極小殘值(Generalized Minimum RESidual,GMRES)算法[17-18]通過(guò)對(duì)迭代矩陣進(jìn)行預(yù)處理,降低迭代系數(shù)矩陣的條件數(shù),可有效提高高維病態(tài)線性方程組的求解效率和魯棒性。本文針對(duì)TSM出現(xiàn)的雅可比矩陣對(duì)角占優(yōu)性質(zhì)惡化的問(wèn)題,通過(guò)引入GMRES算法來(lái)改善TSM求解的魯棒性和計(jì)算效率。

1 數(shù)值方法

1.1 時(shí)間譜方法

如果流動(dòng)具有周期性,則在每一個(gè)計(jì)算單元內(nèi)的守恒量U隨時(shí)間均呈周期性變化。在一個(gè)周期T內(nèi)設(shè)置N 個(gè)等時(shí)間間距的采樣點(diǎn),若N為奇數(shù),對(duì)守恒量U進(jìn)行正向傅里葉變換得

式中:Δt=T/N為時(shí)間間隔;k為諧波數(shù);ω為角頻率;n代表第n時(shí)刻。對(duì)應(yīng)的傅里葉反變換為

將式(2)對(duì)時(shí)間求導(dǎo)可得

將式(1)代入式(3)中,化簡(jiǎn)可得[7]

式中:

若N為偶數(shù),則類似進(jìn)行計(jì)算可得

于是,可將含有N個(gè)采樣點(diǎn)在第n個(gè)時(shí)刻的Navier-Stokes方程的半離散形式寫(xiě)為

式中:V為單元體的體積;Rn為殘差向量。采用中心格式計(jì)算黏性通量與無(wú)黏通量。

在式(5)中引入偽時(shí)間項(xiàng),則該方程進(jìn)一步寫(xiě)為

式中:τn為偽時(shí)間。當(dāng)顯式離散偽時(shí)間項(xiàng)時(shí),為了保證計(jì)算穩(wěn)定性,偽時(shí)間步長(zhǎng)須限制為

式中:CFL為CFL數(shù);λ為譜半徑;k′為最大諧波數(shù):

當(dāng)采樣點(diǎn)數(shù)增加或減縮頻率增大時(shí),偽時(shí)間步長(zhǎng)會(huì)明顯減小從而嚴(yán)重降低計(jì)算效率。采用隱式格式可以放寬偽時(shí)間步長(zhǎng)限制,因此在計(jì)算中多采用隱式格式對(duì)偽時(shí)間項(xiàng)進(jìn)行離散。對(duì)偽時(shí)間導(dǎo)數(shù)進(jìn)行一階向后差分離散得到式(8)所示的全隱式格式:

式中:其中:Jn為第n個(gè)時(shí)刻的通量雅可比矩陣。當(dāng)矩陣A非對(duì)角項(xiàng)為零時(shí),視為進(jìn)行顯式處理。

1.2 隱式求解穩(wěn)定性

常用時(shí)間推進(jìn)格式如LU-SGS方法,對(duì)稱SGS方法均要求系數(shù)矩陣A對(duì)角占優(yōu)。但是TSM增加了非對(duì)角項(xiàng)元素會(huì)減弱矩陣A的對(duì)角占優(yōu)特性,增加矩陣A條件數(shù)。從而導(dǎo)致常用時(shí)間推進(jìn)格式的計(jì)算效率降低甚至?xí)l(fā)生不收斂情況。非對(duì)角項(xiàng)中包含了線性化的空間偏導(dǎo)數(shù)項(xiàng)和時(shí)間偏導(dǎo)數(shù)項(xiàng)。從物理意義上分析,空間通量近似為

式中:S為網(wǎng)格單元表面積;un為當(dāng)前網(wǎng)格點(diǎn)速度;a為聲速。采用TSM引入非對(duì)角項(xiàng)近似為

Atime～Nω (10)則非對(duì)角時(shí)間導(dǎo)數(shù)離散項(xiàng)與原始空間導(dǎo)數(shù)離散項(xiàng)之比近似為

從式(11)中可以看出,角頻率ω增大或采樣點(diǎn)數(shù)N增加,均導(dǎo)致系數(shù)矩陣主對(duì)角占優(yōu)特性變差。若此時(shí)仍采用LU-SGS或?qū)ΨQSGS等要求系數(shù)矩陣主對(duì)角嚴(yán)格占優(yōu)的方法求解,計(jì)算穩(wěn)定性和收斂性會(huì)嚴(yán)重降低。采用不同減縮頻率f流動(dòng)計(jì)算的收斂歷史對(duì)比如圖1所示。采用不同采樣點(diǎn)數(shù)計(jì)算的收斂歷史對(duì)比如圖2所示。圖1和圖2中縱坐標(biāo)lg(Residmax)為迭代殘差的對(duì)數(shù),從圖中可以明顯看出,當(dāng)減縮頻率f增大或采樣點(diǎn)數(shù)增加時(shí),收斂速度明顯減慢。對(duì)于大減縮頻率或采樣點(diǎn)數(shù)多的情況,采用對(duì)稱SGS迭代方法已經(jīng)不能保證TSM的高效,需要尋找一種在系數(shù)矩陣A相對(duì)病態(tài)的情況下仍能高效求解的算法來(lái)進(jìn)行時(shí)間推進(jìn)。

1.3 GMRES算法

GMRES以Galerkin原理為基礎(chǔ),先建立Krylov子空間1組標(biāo)準(zhǔn)正交基,然后在Krylov子空間中利用最小二乘法得到解向量。其中,可以采用Arnoldi方法或修正的Gram-Schmidt方法(Modified Gram-Schmidt,MGS)計(jì)算Krylov子空間Km(A,v)中的標(biāo)準(zhǔn)正交基。帶預(yù)處理的GMRES算法具體過(guò)程為[19]

1)給定初值x(0),然后計(jì)算r0=b-Ax(0),對(duì)r0預(yù)處理,即r0=ML-1r0,記θ=(r0, r0)1/2,得到第一個(gè)正交基v1=r0/θ。

2)進(jìn)入GMRES(m)迭代,j=1,2,…,kmax。

① 計(jì)算wj=Avj,并對(duì)wj進(jìn)行預(yù)處理,wj=ML-1wj。

② 計(jì)算hij=(wj,vi),i=1,2,…,j。

j

③計(jì)算

⑤ 可通過(guò)hj+1,j＜δ來(lái)判定是否已經(jīng)收斂。δ為給定內(nèi)迭代收斂標(biāo)準(zhǔn),其中在沒(méi)有滿足此條件下記錄所求得的殘值向量作為初始?xì)堉迪蛄坷^續(xù)返回步驟1)重新迭代,直至滿足hj+1,j＜δ,停止該時(shí)間步的迭代,即為重啟型的GMRES算法。

3)求解方程組的近似解方程組新解可表示為

式中:ym=arg minθe1-Hˉmz2為待定的m維實(shí)z型列矢量,θ為殘差范數(shù),Hˉm為Hessenberg矩陣,e1為單位向量,使 θe1-Hˉmz2達(dá)到最小的z向量即為所求的ym。具體求解方法見(jiàn)文獻(xiàn)[19]。

在本文中需要對(duì)式(8)進(jìn)行迭代求解。在每個(gè)偽時(shí)間步將其簡(jiǎn)化為線性方程組Ax=b形式,其中x對(duì)應(yīng)式(8)中的ΔU,A對(duì)應(yīng)式(8)中的系數(shù)矩陣A,b對(duì)應(yīng)式(8)中等式的右邊項(xiàng)。

對(duì)于需要求解的式(8),當(dāng)左端系數(shù)矩陣條件數(shù)cond(A)=A A-1較大時(shí),方程組是病態(tài)的,此時(shí)傳統(tǒng)GMRES算法不具有高速收斂特性,如果要保證計(jì)算的高效性,就必須對(duì)系數(shù)矩陣進(jìn)行有效的預(yù)處理。

所謂預(yù)處理,就是將方程組轉(zhuǎn)換成另外一個(gè)同解方程組,而后者比前者更容易求解。預(yù)處理的實(shí)質(zhì)就是在不改變方程解的情況下,將矩陣特征值集中在某個(gè)值附近,從而改善矩陣的性質(zhì),加快迭代法的收斂速度。預(yù)處理將線性方程組Ax=b轉(zhuǎn)換為同解的方程組M-1Ax=M-1b(左預(yù)處理)或AM-1Mx=b(右預(yù)處理)[20]。

本文采用的是左預(yù)處理方法并通過(guò)LU-SGS方法對(duì)預(yù)處理矩陣的逆矩陣和列向量的乘積進(jìn)行處理。

2 數(shù)值模擬及分析

本文對(duì)強(qiáng)迫俯仰振蕩的NACA0012翼型進(jìn)行數(shù)值模擬。計(jì)算網(wǎng)格如圖3所示,網(wǎng)格單元總數(shù)為8 568,節(jié)點(diǎn)總數(shù)為6 004,翼型表面節(jié)點(diǎn)為180,附面層第一層高度為1×10-5,附面層總高度為0.02。用有限體積法求解,空間離散采用中心格式,時(shí)間推進(jìn)采用雙時(shí)間步法。

2.1 跨聲速算例

通過(guò)跨聲速標(biāo)模AGARD1對(duì)TSM進(jìn)行驗(yàn)證并研究不同算法TSM的計(jì)算效率。AGARD1算例在迎角較大時(shí)會(huì)出現(xiàn)跨聲速。

翼型俯仰運(yùn)動(dòng)規(guī)律為

式中:αm為平均迎角;α0為俯仰運(yùn)動(dòng)振幅;ω為俯仰振蕩的角速度。角速度ω與無(wú)量綱減縮頻率f的關(guān)系為

式中:c為翼型根部弦長(zhǎng);U∞為自由來(lái)流速度。記xca為質(zhì)心位置,表1列出了主要的計(jì)算條件。

物理時(shí)間離散分別采用TSM和BDF,收斂標(biāo)準(zhǔn)均為5.0×10-7。采用BDF計(jì)算時(shí),每個(gè)時(shí)間周期分為120個(gè)時(shí)間等份,每一步子迭代最大為1 000步,共進(jìn)行4個(gè)周期的計(jì)算。在后面計(jì)算中收斂標(biāo)準(zhǔn)以及BDF的子迭代條件保持不變。對(duì)于AGARD1算例,圖4給出了TSM與BDF計(jì)算得到的升力系數(shù)CL和力矩系數(shù)Cm并與實(shí)驗(yàn)值[21]對(duì)比。從圖中可以看出對(duì)于該算例的升力系數(shù)CL,TSM僅需要5個(gè)采樣點(diǎn)就可以達(dá)到與BDF同樣精度。而對(duì)于力矩系數(shù),TSM則需要7個(gè)采樣點(diǎn)才可以達(dá)到與BDF同樣的精度。

表1 跨聲速計(jì)算條件Table 1 Transonic computational conditions

圖5給出了NACA0012翼型作俯仰振蕩時(shí)用TSM計(jì)算的收斂歷史。通過(guò)圖5可以看出,采用GMRES算法收斂速度明顯快于對(duì)稱SGS迭代。且隨著采樣點(diǎn)數(shù)增多,GMRES算法的收斂速度并沒(méi)有太大的變化。取采樣點(diǎn)數(shù)為4～9,表2對(duì)比了采用對(duì)稱SGS以及GMRES算法推進(jìn)的計(jì)算效率。從表2中不難看出,隨著采樣點(diǎn)數(shù)增加,系數(shù)矩陣的條件數(shù)增加,對(duì)角占優(yōu)特性減弱,GMRES算法的優(yōu)勢(shì)就體現(xiàn)的更為明顯,甚至采用GMRES算法的計(jì)算效率可以達(dá)到對(duì)稱SGS迭代的6倍。當(dāng)N=9時(shí),采用對(duì)稱SGS迭代計(jì)算會(huì)發(fā)散,對(duì)稱SGS迭代已經(jīng)失效。

表2 不同時(shí)間推進(jìn)格式效率對(duì)比Table 2 Comparison of efficiency of different time marching methods

通過(guò)對(duì)比可以看出,采用基于GMRES算法的TSM相比于BDF方法計(jì)算效率提高了9倍左右,因此該方法有很好的應(yīng)用前景。后面的算例沒(méi)有特殊說(shuō)明,TSM均采用GMRES算法進(jìn)行偽時(shí)間推進(jìn)。

2.2 二維亞聲速俯仰振蕩

翼型的俯仰運(yùn)動(dòng)規(guī)律以及計(jì)算網(wǎng)格均與跨聲速算例相同,表3列出了亞聲速算例的計(jì)算條件。物理時(shí)間離散分別采用TSM和BDF。

圖6給出了NACA0012翼型作俯仰振蕩時(shí)用TSM計(jì)算的收斂歷史,通過(guò)圖6可以看出該算例在采樣點(diǎn)數(shù)N≤7時(shí)其收斂速度并非隨著采樣點(diǎn)數(shù)的增加而單調(diào)下降,這是因?yàn)樵诓蓸狱c(diǎn)數(shù)較少時(shí),雖然Atime有所增加,但Aspatial的影響仍占主導(dǎo),其收斂速度主要與采樣點(diǎn)定常流場(chǎng)求解特性有關(guān)。當(dāng)N＞7時(shí)采用TSM引入的非對(duì)角時(shí)間導(dǎo)數(shù)項(xiàng)Atime占主導(dǎo)作用,具體體現(xiàn)在隨著采樣點(diǎn)數(shù)增加,收斂速度呈現(xiàn)明顯的單調(diào)下降趨勢(shì)。圖7給出了不同采樣點(diǎn)數(shù)下的TSM計(jì)算效率。從圖中可以明顯看出,對(duì)于此算例而言,本文的TSM計(jì)算效率遠(yuǎn)遠(yuǎn)高于BDF。當(dāng)N在4～11范圍內(nèi)時(shí),隨著采樣點(diǎn)數(shù)增加,TSM計(jì)算時(shí)間近似呈線性變化。當(dāng)N=4時(shí)TSM的計(jì)算效率約是BDF的10倍,當(dāng)N=12時(shí)TSM的計(jì)算效率約是BDF的2倍。

表3 亞聲速計(jì)算條件Table 3 Subsonic computational conditions

圖8給出了TSM與BDF計(jì)算得到的升力系數(shù)CL和力矩系數(shù)Cm對(duì)比。可以看出對(duì)于該亞聲速俯仰振蕩算例,當(dāng)TSM采樣點(diǎn)數(shù)取N=5時(shí)就可以達(dá)到與BDF相當(dāng)?shù)木?而計(jì)算效率提高了約10倍。因此,對(duì)于本算例而言,TSM的計(jì)算效率比BDF提高了一個(gè)量級(jí)。

2.3 非簡(jiǎn)諧周期運(yùn)動(dòng)算例

對(duì)于翼型強(qiáng)迫非簡(jiǎn)諧俯仰振蕩算例,本文借鑒文獻(xiàn)[16]中計(jì)算條件,通過(guò)兩個(gè)諧波的周期性運(yùn)動(dòng)為例來(lái)研究TSM在非簡(jiǎn)諧周期運(yùn)動(dòng)中的計(jì)算效率。采用Euler方程求解。翼型俯仰運(yùn)動(dòng)規(guī)律為式中:ω1與ω2分別為翼型作強(qiáng)迫俯仰振蕩的兩個(gè)頻率;α1與α2為對(duì)應(yīng)的俯仰運(yùn)動(dòng)振幅。計(jì)算狀態(tài)如表4所示。

圖9給出了BDF與TSM計(jì)算出的升力和阻力隨時(shí)間變化曲線。從圖10中可以看出采用TSM計(jì)算時(shí)需要至少7個(gè)點(diǎn)升力系數(shù)才與BDF的計(jì)算結(jié)果吻合,而要保證阻力系數(shù)與BDF計(jì)算結(jié)果吻合則需要13個(gè)點(diǎn)。表5對(duì)比了采用13個(gè)點(diǎn)的TSM與BDF的計(jì)算效率,從表中可以看出,雖然TSM需要13個(gè)點(diǎn),但采用TSM相比于傳統(tǒng)BDF計(jì)算效率仍能提高6倍。

表4 非簡(jiǎn)諧周期運(yùn)動(dòng)計(jì)算條件Table 4 Non harmonic periodic motion computational conditions

表5 非簡(jiǎn)諧周期運(yùn)動(dòng)計(jì)算時(shí)間對(duì)比Table 5 Comparison of computational time in non harmonic periodic motion

因此,對(duì)于非簡(jiǎn)諧周期性運(yùn)動(dòng),雖然需要較多的采樣點(diǎn)數(shù),但帶預(yù)處理的GMRES方法仍然可以保證TSM有較高的計(jì)算效率。

3 結(jié) 論

本文將GMRES算法應(yīng)用到時(shí)間譜方法中,通過(guò)LU-SGS預(yù)處理降低系數(shù)矩陣條件數(shù),解決了當(dāng)采樣點(diǎn)數(shù)增加或減縮頻率增大所導(dǎo)致的普通迭代法失效的問(wèn)題。本文采用了跨聲速AGARD1翼型俯仰振蕩、亞聲速翼型簡(jiǎn)諧和非簡(jiǎn)諧周期運(yùn)動(dòng)3個(gè)標(biāo)準(zhǔn)算例,計(jì)算結(jié)果顯示:

1)采用GMRES算法進(jìn)行求解比對(duì)稱SGS迭代計(jì)算效率有很大提高,且隨著采樣點(diǎn)數(shù)增加,采用GMRES算法效率優(yōu)勢(shì)更明顯。

2)一定范圍內(nèi)隨著采樣點(diǎn)數(shù)增加,采用GMRES算法求解時(shí)間譜方法的收斂速度有所減慢但變化不大,計(jì)算時(shí)間近似呈線性變化。

3)對(duì)于線性或非線性較弱的算例,無(wú)論是簡(jiǎn)諧或非簡(jiǎn)諧周期性運(yùn)動(dòng),采用GMRES算法的TSM相比于時(shí)域差分法具有明顯優(yōu)勢(shì),可將計(jì)算效率提高一個(gè)量級(jí)。

[1] HALL K C,CRAWLEY E F.Calculation of unsteady flows in turbomachinery using the linearized Euler equations[J].AIAA Journal,1989,27(6):777-787.

[2] NING W,HE L.Computation of unsteady flows around oscillation blades using linear and nonlinear harmonic Euler methods[J].Journal of Turbomachinery,1998,120(3):508-514.

[3] HALL K C,THOMAS J P,CLARK W S.Computation of unsteady nonlinear flows in cascades using a harmonic balance technique[J].AIAA Journal,2002,40(5):879-886.

[4] MCMULLEN M,JAMESON A,ALONSO J J.Acceleration of convergence to a periodic steady state in turbomachinery flows:AIAA-2001-0152[R].Reston:AIAA,2001.

[5] MCMULLEN M,JAMESON A,ALONSO J J.Application of a nonlinear frequency domain solver to the Euler and Navier-Stokes equations:AIAA-2002-0120[R].Reston:AIAA,2002.

[6] DAI H H,YUE X K,YUAN J P,et al.A time domain collocation method for studying the aeroelasticity of a two dimensional airfoil with a structural nonlinearity[J].Journal of Computational Physics,2014,270:214-237.

[7] GOPINATH A K,JAMESON A.Time spectral method for periodic unsteady computations over two-and threedimensional bodies:AIAA-2005-1220[R].Reston:AIAA,2005.

[8] VAN DER WEIDE E,GOPINATH A K,JAMESON A.Turbomachinery applications with the time spectral method:AIAA-2005-4905[R].Reston:AIAA,2005.

[9] CHOI S,POTSDAM M,LEE K,et al.Helicopter rotor design using a time-spectral and adjoint-based method:AIAA-2008-5810[R].Reston:AIAA,2008.

[10] GOMAR A,BONVY Q,SICOT F,et al.Convergence of Fourier-based time methods for turbomachinery wake passing problems[J].Journal of Computational Physics,2014,278(C):229-256.

[11] GOPINATH A K,JAMESON A.Application of the time spectral method to periodic unsteady vortex shedding:AIAA-2006-0449[R].Reston:AIAA,2006.

[12] 楊小權(quán),程蘇堃,楊愛(ài)明,等.基于時(shí)間譜方法的振蕩翼型和機(jī)翼非定常黏性繞流數(shù)值模擬[J].航空學(xué)報(bào),2013,34(4):787-797.YANG X Q,CHENG S K,YANG A M,et al.Time spectral method for numerical simulation of unsteady viscous flow over oscillating airfoil and wing[J].Acta Aeronautica et Astronautica Sinica,2013,34(4):787-797(in Chinese).

[13] 謝立軍,楊云軍,劉周,等.基于時(shí)間譜方法的飛行器動(dòng)導(dǎo)數(shù)高效計(jì)算技術(shù)[J].航空學(xué)報(bào),2015,36(6):2016-2026.XIE L J,YANG Y J,LIU Z,et al.A high efficient method for computing dynamic derivatives of aircraft based on time spectral method[J].Acta Aeronautica et Astronautica Sinica,2015,36(6):2016-2026(in Chinese).

[14] 詹磊,劉鋒.應(yīng)用傅里葉時(shí)間譜方法求解二維跨音速流動(dòng)問(wèn)題解的精度研究[J].航空工程進(jìn)展,2015,6(4):395-404.ZHAN L,LIU F.Study on accuracy of the solutions using the Fourier time spectral method for two-dimensional transonic flows[J].Advances in Aeronautical Science and Engineering,2015,6(4):395-404(in Chinese).

[15] JAMESON A,SHANKARAN S.An assessment of dualtime stepping,time spectral and artificial compressibility based numerical algorithms for unsteady flow with applications to flapping wings[C]//19th AIAA Computational Fluid Dynamics.Reston:AIAA,2009.

[16] MAVRIPLIS D J,YANG Z.Time spectral method for periodic and quasi-periodic unsteady computations on unstructured meshes[J].Mathematical Modelling of Natural Phenomena,2011,6(3):213-236.

[17] SAAD Y,SCHULTZ M H.GMRES:A generalized minimal residual algorithm for solvingnonsymmetric linear systems[J].SIAM Journal on Scientific and Statistical Computing,1986,7(3):856-869.

[18] JAMESON A,YOON S.Lower-upper implicit schemes with multiple grids for the Euler equations[J].AIAA Journal,1987,25(7):929-935.

[19] 李春娜,葉正寅,王剛.基于二維非結(jié)構(gòu)網(wǎng)格的GMRES隱式算法[J].西北工業(yè)大學(xué)學(xué)報(bào),2007,25(5):630-635.LI C N,YE Z Y,WANG G.GMRES implicit algorithm based on 2D unstructured meshes for solving Euler equations[J].Journal of Northwestern Polytechnical University,2007,25(5):630-635(in Chinese).

[20] 劉巍,張理論,王勇獻(xiàn),等.計(jì)算空氣動(dòng)力學(xué)并行編程基礎(chǔ)[M].北京:國(guó)防工業(yè)出版社,2013.LIU W,ZHANG L L,WANG Y X,et al.Computational aerodynamics parallel programming basis[M].Beijing:National Defence Industry Press,2013(in Chinese).

[21] LANDON R H.NACA0012 oscillatory and transient pitching:AGARD-R-702[R].Paris:AGARD,1982.

Efficient GMRES algorithm in time spectral method

GONG Yiming1,LlU Zhanhe2,LlU Yilang1,ZHANG Weiwei1,*

1.School of Aeronautics,Northwestern Polytechnical University,Xi'an 710072,China
2.School of Aeronautic Engineering,Zhengzhou University of Aeronautics,Zhengzhou 450046,China

ln this paper,the computational efficiency of the time-spectral method for solving the periodic unsteady flow field is studied,and the implicit method of time spectral method for solving the periodic unsteady flow is discussed.When the number of sampling points increases or the reduced frequency magnifies,the diagonal dominant property of the Jacobian matrix corresponding to the time spectral method deteriorates rapidly,resulting in the failure of many traditional iterative methods.ln order to solve the problems above,the generalized minimum residual(GMRES)algorithm with preprocessing is used to improve the computational convergence of the Jacobian matrix.The time spectral method is used to compute the NACA0012 airfoil forced oscillation,and the computational efficiency and accuracy is compared with that of the time-domain difference method.The results show that the time spectral method can generally improve the computational efficiency an order of magnitude with saturated computational accuracy.For the transonic periodic flow,the GMRES algorithm is superior to SGS iterative algorithm both in stability and computational convergence.

time spectral method;generalized minimum residual(GMRES)algorithm;periodic unsteady flow;preprocess;computational efficiency

2016-10-27;Revised:2016-12-04;Accepted:2016-12-21;Published online:2016-12-26 15:17

URL:www.cnki.net/kcms/detail/11.1929.V.20161226.1517.004.html

s:National Natural Science Foundation of China for Excellent Young Scholars(11622220);Programme of lntroducing Talents of Discipline to Universities(B17037)

V211.3

A

1000-6893(2017)07-120894-09

10.7527/S1000-6893.2016.120894

2016-10-27;退修日期:2016-12-04;錄用日期:2016-12-21;網(wǎng)絡(luò)出版時(shí)間:2016-12-26 15:17

www.cnki.net/kcms/detail/11.1929.V.20161226.1517.004.html

國(guó)家自然科學(xué)基金優(yōu)秀青年基金(11622220);高等學(xué)校創(chuàng)新引智計(jì)劃(B17037)

*通訊作者.E-mail:aeroelastic@nwpu.edu.cn

貢伊明,劉戰(zhàn)合,劉溢浪,等.時(shí)間譜方法中的高效GMRES算法[J].航空學(xué)報(bào),2017,38(7):120894.GONG Y M,LlU Z H,LlU Y L,et al.Efficient GMRES algorithm in time spectral method[J].Acta Aeronautica et Astronautica Sinica,2017,38(7):120894.

(責(zé)任編輯:李明敏)

*Corresponding author.E-mail:aeroelastic@nwpu.edu.cn