貢伊明, 張偉偉, 劉溢浪
(西北工業(yè)大學(xué) 航空學(xué)院, 陜西 西安 710072)
?
非定常求解的內(nèi)迭代初值對計算效率的影響研究
貢伊明, 張偉偉, 劉溢浪
(西北工業(yè)大學(xué) 航空學(xué)院, 陜西 西安710072)
摘要:基于非定常流場的雙時間求解方法,提出了一種提高非定常流場求解效率的有效策略。通過對前幾個時刻的流場信息進行外插來確定下一時刻的迭代初值,使之更接近于收斂解,降低內(nèi)迭代初始殘值,進而提高了非定常流場的求解效率。將流場中每個點的守恒量在時間方向上進行泰勒級數(shù)展開,設(shè)計了若干種外插格式。采用繞圓柱非定常流動的求解來驗證本方法的計算效果,并研究了不同初值外插格式、空間離散格式、時間步長和收斂標(biāo)準(zhǔn)對初值外插方法效果的影響。研究表明,在雙時間步法基礎(chǔ)上,采用初值外插策略可普遍提高計算效率,其中交替外插策略可以普遍將求解效率提高1倍左右。相比于迎風(fēng)格式,該方法對中心格式的求解效率提高更顯著,并且對于不同的收斂標(biāo)準(zhǔn)和時間步數(shù)均有非常明顯的效果。
關(guān)鍵詞:初值外插;非定常流場求解;泰勒級數(shù)展開;雙時間步法
近些年來,隨著計算機性能的提高,計算流體力學(xué)(computational fluid dynamics,簡稱CFD)在很多工作中得到廣泛運用。在實際的工程應(yīng)用中,很多問題涉及到非定常流動,并且其計算量遠大于定常流動,比如,氣動彈性動力學(xué)問題、渦輪葉片繞流問題、分離渦脫落問題、動導(dǎo)數(shù)的求解等等。由此而促進了各種加速收斂技術(shù)和時間推進格式的發(fā)展。
為了提高計算效率,前人發(fā)展了很多加速收斂措施,如當(dāng)?shù)貢r間步長、隱式殘值光順方法和多重網(wǎng)格方法[1]。隱式時間推進方法可以在保證穩(wěn)定性的同時大幅提高時間步長,進而提高計算效率,常用的有近似因子分解(AF)[2]、超松弛迭代方法(SSOR)[3]、直接求解方法和Newton-Krylov方法[4]。在近似因子分解方法中應(yīng)用最廣的是LUSGS算法[5]。后來Chen和Wang[6]對LUSGS進行了改進,在有限體積方法中發(fā)展了BLUSGS算法。Newton-Krylov方法中最受歡迎的是Newton-GMRES[7],其中無矩陣存儲的GMRES方法[8]在CFD中應(yīng)用更為廣泛。
非定常流場的計算方法可分為2大類:雙時間步方法(dual-time stepping)和物理時間迭代方法(physical time subiteration)。雙時間步法引入了偽時間項,對于每一個物理時刻看做定常流動來求解。雙時間步法的好處是原來用于定常計算的預(yù)處理、當(dāng)?shù)貢r間步、對角化、多重網(wǎng)格都可以應(yīng)用到計算中;真實時間步長可以取得很大;雙時間法可以降低通量線性化誤差和近似因子分解誤差,放寬了穩(wěn)定性限制[9]。這些方法都能有效的提高收斂速度。
流場求解效率不僅僅與收斂速度有關(guān),還與初值與收斂解間的差距有關(guān)。因此,在基于雙時間方法的非定常流場求解中,為了提高計算效率,通常將上一個時刻的流場作為下一個時刻流場求解的初值。至于從初值角度進一步提高計算效率,相關(guān)研究鮮有報道。本文從該研究角度出發(fā),將前2~3個物理時刻的流場信息外插求解出待求時刻的流場作為偽時間步迭代的初場,這樣就進一步縮小了初值與收斂解的差距,以期減小迭代步數(shù),提高計算效率。
基于以上的思想,本文通過采用內(nèi)迭代初值外插,發(fā)展了一種提高非定常流場求解效率的有效策略。采用繞圓柱非定常流動的求解算例來驗證本方法的計算效果,并研究了不同初值外插格式、空間離
散格式、時間步長和收斂標(biāo)準(zhǔn)對初值外插方法效果的影響。
1數(shù)值方法
1.1雙時間步方法
為了提高非定常流場的時間計算精度,同時又要求具有較高的計算效率,Jameson 提出了一種雙時間步方法,即在凍結(jié)的物理時間點上加入類似 Newton 迭代的虛擬時間迭代過程,通過增加內(nèi)迭代過程提高 LU-SGS 等隱式線性化方法所損失的時間精度。雙時間步法思想簡單,且在定常流計算程序基礎(chǔ)上進行改造的工作量小,得到了廣泛應(yīng)用,其實現(xiàn)形式如下。
非定常流動控制方程離散后的形式為
(1)
(2)
這里的Qn表示第個物理時間步的守恒量,Δt代表物理時間步長。
方程(2)直接求解比較困難,雙時間步法采用虛時間迭代技術(shù)對方程(2)進行求解,引入一個虛擬時間τ,將控制方程(2)改寫為
(3)
(4)
1.2初值外插方法
初值外插方法是通過記錄前幾個時刻的流場的基本信息進行外插得到下一個物理時刻的內(nèi)迭代初值。
一般情況下,采用外插時必須要保證插值對象的連續(xù)性。若其出現(xiàn)比較強的間斷,則外插出現(xiàn)的偏差往往會比較大,這樣就會與我們的初衷背道而馳。因此,合理選擇外插對象至關(guān)重要。相對于原始變量,守恒變量的連續(xù)性更好,因此我們將流場的守恒量進行外插。但對于流場本身存在強間斷情況,比如有激波出現(xiàn),此時用初值外插方法容易造成更大偏差,不宜使用。在時間方向上的泰勒級數(shù)展開如下:
(5)
如果時間方向上取一階精度,時間導(dǎo)數(shù)采用一階向后差分,即
(6)
代入(5)式并忽略二階項可得一階外插格式
(7)
同樣的,如果時間方向上取二階精度,一階時間導(dǎo)數(shù)采用中心差分, 二階導(dǎo)數(shù)采用二階向后差分并忽略三階項可得二階外插格式
(8)
同理可得三階精度的外插格式
(9)
本文還提出了一階與三階交替的外插格式。即在第n個物理時刻,取2m為周期有:
若mod(n,2m) 若mod(n,2m)>m-1,則下一時刻的初值為Qn+1=4Qn-6Qn-1+4Qn-2-Qn-3 根據(jù)多次計算對比結(jié)果,取m=3時效果更好。因此后面的交替外插都采用的是每6個時間點為1個周期,三步一階外插三步三階外插的格式。 2數(shù)值算例驗證 2.1算例介紹 以非定常靜止圓柱繞流為例,對于不同外插格式和空間離散格式下的初值外插方法進行效率對比并分析時間步長和收斂標(biāo)準(zhǔn)對初值外插方法效率的影響。來流馬赫數(shù)為0.1,雷諾數(shù)為100,層流,采用非結(jié)構(gòu)三角形網(wǎng)格,圓柱表面網(wǎng)格節(jié)點數(shù)為100,網(wǎng)格總數(shù)為29 788,計算域網(wǎng)格和圓柱表面局部網(wǎng)格示意圖如圖1所示。 圖1 層流圓柱繞流計算網(wǎng)格示意圖 2.2不同因素影響對比 在初值外插方法介紹中給出了一階、二階、三階和交替外插4種外插格式。本節(jié)針對采用不同的外插格式、空間離散格式對初值外插方法效率的影響進行分析并確定使該方法效率最好的外插格式和空間離散格式,并分析不同的時間步長和收斂標(biāo)準(zhǔn)對初值外插方法效率的影響。 在非定常圓柱繞流算例中,空間離散采用中心格式,取時間步長為1個周期120步,設(shè)1個周期的步數(shù)為Tn,則有Tn=120。設(shè)每個周期平均耗時為Tp。收斂要求為最大殘值小于5.0×10-9。具體狀態(tài)在2.1節(jié)中已經(jīng)交代。對一階、二階、三階和交替外插精度的外插效果進行對比。以下是計算結(jié)果: 表1 4種格式效率對比表 圖2 不同外插精度計算效果圖 由圖2可以看出,無論是初始迭代的最大誤差還是平均誤差,均隨著外插精度的提高而降低。而迭代步數(shù)并不滿足該規(guī)律,初始迭代誤差最小的三階精度外插的迭代步數(shù)反而是最多的,換而言之,采用三階精度外插格式會出現(xiàn)收斂困難的現(xiàn)象。對于出現(xiàn)收斂困難現(xiàn)象的原因,在下一節(jié)中會具體分析。 從表1中不難看出一階精度、二階精度和交替外插格式均可明顯減少每周期迭代步數(shù),但由于采用初值外插方法每一步迭代增加了計算量和存儲因而每一步迭代所用的時間會有所增加。比較后發(fā)現(xiàn)平均每迭代一步所用的時間增加量很小,大約10%左右。因此整體的效率仍然是提高的。在4種格式中采用三步一階三步三階的交替外插格式效率最高。因此,本文后面的計算沒有特殊說明均采用三 圖3 采用初值外插法前后的升力與阻力系數(shù)對比 步一階三步三階的交替外插格式。從圖3中可以看出采用初值外插方法前后計算結(jié)果完全吻合,證明了初值外插方法的準(zhǔn)確性。 采用不同的空間離散格式,分析對初值外插方法的效率是否會有所影響以及影響程度如何。本文對比了中心格式和迎風(fēng)ROE格式2種空間離散格式,時間步長取1個周期120步,收斂標(biāo)準(zhǔn)為5.0×10-9。以下是對比結(jié)果: 表2 不同空間離散格式計算效率表 從表2中可以看出采用不同的空間離散格式對于初值外插方法效率會有一定的影響??臻g離散采用中心格式初值外插方法的效果更好。對比中心格式和迎風(fēng)格式不難發(fā)現(xiàn),中心格式利用了較多的網(wǎng)格系統(tǒng)作為模板,因此對于模板的迭代初值更為敏感。當(dāng)?shù)踔蹈咏諗拷鈺r,采用中心格式比采用迎風(fēng)格式提高的效率就會更多。 時間步長和收斂標(biāo)準(zhǔn)選取的不同對初值外插方法的效率也會有所影響。本文通過改變時間步長和收斂標(biāo)準(zhǔn)觀察初值外插方法的效率變化情況。 不改變其他的參數(shù),僅對時間步長進行改變,分析對該方法效率的影響。以中心格式為例,對于時間步長分別取1個周期60步、120步、200步和400步的情況,收斂標(biāo)準(zhǔn)為最大殘值小于1.0×10-8。圖4和表3顯示了初值外插方法的效率變化情況。 圖4 不同時間步長初值外插法迭代步數(shù)圖 Tn無外插Tp/s外插法Tp/s效率提高/%60630.92340.1685.48120734.42381.1192.71200830.92424.6495.684001014.67537.2888.85 從以上圖表中不難看出,時間步長的改變雖對外插方法的計算效率提升量有一定影響,但影響不大。因此,可以認(rèn)為在很大的時間步數(shù)變化范圍內(nèi)初值外插方法均有非常明顯的效果。 下面研究收斂標(biāo)準(zhǔn)對外插方法效率提升的影響??臻g離散采用中心格式,固定時間步長為1個周期120步,收斂標(biāo)準(zhǔn)分別取最大誤差小于5.0×10-9、1.0×10-8、5.0×10-8和1.0×10-7時觀察初值外插方法的效果。 圖5 不同收斂標(biāo)準(zhǔn)初值外插法迭代步數(shù)圖 收斂標(biāo)準(zhǔn)無外插Tp/s外插法Tp/s效率提高/%5.0×10-9884.21452.2795.501.0×10-8728.50354.17105.695.0×10-8642.26289.72121.681.0×10-7609.29262.40132.20 從以上圖表中可以看出在一定范圍內(nèi)隨著收斂標(biāo)準(zhǔn)降低,初值外插法的效率會提高。但即使對于高精度的數(shù)值求解,初值外插法的效率仍然可以將效率提高將近一倍。因此初值外插法對于不同精度的求解均具有非常明顯的效果。 2.3流場誤差傳播分析 一般情況下,外插精度越高,流場迭代初值與收斂解越接近,初值外插方法計算效率就越高。然而三階精度外插的計算效率卻最低,甚至還沒有不用初值外插方法的效果好。產(chǎn)生這種現(xiàn)象的原因歸結(jié)于誤差的2個分量。 不同的初值與真實解的誤差可以按照誤差各個傅里葉分量的衰減程度不同,分為高頻振蕩誤差和低頻光滑誤差。 在線性空間內(nèi)任何初始誤差均可表示為 (10) 式中,E0為初始誤差,vi表示第i個特征向量,ai表示誤差第i個分量對應(yīng)的系數(shù)。 在迭代過程中 (11) 很明顯,特征向量代表誤差分量,其對應(yīng)的特征值代表誤差分量的放大系數(shù)。誤差每一個分量vi都對應(yīng)一個誤差頻率k,當(dāng)1≤k 高頻振蕩誤差是局部行為,來源于附近幾個網(wǎng)格點之間的相互耦合,與邊界距離較遠的網(wǎng)格點的信息無關(guān);而低頻光滑誤差是全局行為,主要來源于邊界信息。通過局部松弛后誤差呈現(xiàn)光滑性,此時誤差主要來自于邊界。可以設(shè)想二維N×N網(wǎng)格上的點松馳方法, 將邊界信息傳播到所有點至少需O(N)次迭代, 因此收斂速度極慢。低頻光滑誤差量階很小,但在細網(wǎng)格上很難收斂[13]。 在迭代過程中低頻誤差的大小與迭代初值和迭代格式有關(guān)。采用初值外插方法會對迭代初值進行改變,低頻誤差的大小也會隨之發(fā)生變化。在外插過程中,不僅僅是對精確解的插值,由于收斂標(biāo)準(zhǔn)的限制,存在一定的誤差,因此外插的對象還包括殘余的誤差。高精度外插對插值對象的連續(xù)性要求較高,而殘余誤差的變化本身具有隨機性與不連續(xù)性。低頻誤差由于在迭代過程中難以收斂在殘余誤差中占有較高的比例。通過高階精度外插格式得到的初始迭代誤差比低階精度外插格式小,但是高階精度外插格式只可以大幅衰減高頻誤差,低頻誤差反而被放大,低頻誤差收斂速度極慢因此高精度外插會出現(xiàn)收斂困難情況。 一階外插能夠保證誤差變化的單調(diào)性,不會出現(xiàn)低頻誤差的放大。但是精度較低,高頻誤差分量衰減幅度小;高階外插雖然可以大幅度衰減高頻誤差,但是會出現(xiàn)低頻誤差的放大。通過采用三步一階三步三階的交替外插格式,既可以大幅衰減高頻誤差分量,又能保證低頻誤差不會被放大。 3結(jié)論 本文在非定常流場雙時間步法求解的基礎(chǔ)上,提出了將前2~3個物理時刻的流場信息外插求解出待求時刻的流場作為偽時間步迭代初場的初值外插方法,并通過繞圓柱非定常流動算例驗證了本方法的有效性,研究結(jié)果顯示: 1) 相對統(tǒng)一精度的外插格式,交替精度的外插格式具有更高的魯棒性和計算效率,可普遍將計算效率提高一倍左右。 2) 不同的空間離散格式對初值外插方法的效果有一定的影響,相對于迎風(fēng)格式,中心格式采用初值外插方法的效率會更高。 3) 初值外插方法的效率提升量會受時間步長和收斂標(biāo)準(zhǔn)影響,但在較寬的時間步長和收斂標(biāo)準(zhǔn)范圍內(nèi),效率提升仍十分可觀。 參考文獻: [1]燕振國. 高精度混合線性緊致格式的隱式時間推進方法研究[D]. 成都: 中國空氣動力研究與發(fā)展中心, 2013 YanZhenguo.InvestigationofImplicitTimeIntegrationMethodswithHDCSSchemes[D].Chengdu:ChinaAerodynamicsResearchandDevelopmentCenterGraduateSchool, 2013 (inChinese) [2]PulliamTH,ChausseeDS.ADiagonalFormofanImplicitApproximate-FactorizationAlgorithm[J].JournalofComputationalPhysics, 1981, 39: 347-363 [3]YoungDM.AHistoricalReviewofIterativeMethods[C]∥TheACMConferenceonHistoryofScientificandNumericComputation, 1987 [4]LomaxH,PulliamTH,ZinggDW,etal.FundamentalsofComputationalFluidDynamics[M].Springer,Berlin, 2001 [5]YoonS,JamesonA.Lower-UpperSymmetricGauss-SeidelMethodfortheEulerandNavier-StokesEquations[J].AIAAJournal, 1988, 26: 1025-1026 [6]ChenRF,WangZJ.Fast,BlockLower-UpperSymmetricGauss-SeidelSchemeforArbitraryGrids[J].AIAAJournal, 2000, 38: 2238-2245 [7]MichalakC,Ollivier-GoochC.GlobalizedMatrix-ExplicitNewton-GMRESfortheHigh-OrderAccurateSolutionoftheEulerEquations[J].Computers&Fluids, 2010, 39: 1156-1167 [8]KnollDA,KeyesDE.Jacobian-FreeNewton-KrylovMethods:aSurveyofApproachesandApplications[J].JournalofComputationalPhysics, 2004, 193: 357-397 [9]趙慧勇, 樂嘉陵. 雙時間步法的應(yīng)用分析[J], 計算物理, 2008, 25(3): 253-258 ZhaoHuiyong,LeJialing.ApplicationAnalysisonDual-TimeStepping[J].ComputationalPhysics, 2008, 25(3): 253-258 (inChinese) [10]JamesonA.TimeDependentCalculationsUsingMultigridwithApplicationstoUnsteadyFlowsPastAirfoilsandWings[C].AIAAPaper91-1596, 1991 [11] 王剛. 復(fù)雜流動的網(wǎng)格技術(shù)及高效、高精度算法研究[D]. 西安: 西北工業(yè)大學(xué), 2005 WangGang.NewTypeofGridGenerationTechniquetogetherwiththeHighEfficiencyandHighAccuracySchemeResearchesforComplexFlowSimulation[D].Xi′an,NorthwesternPolytechnicalUniversity, 2005 (inChinese) [12] 劉超群. 多重網(wǎng)格法及其在計算流體力學(xué)中的應(yīng)用 [M]. 北京: 清華大學(xué)出版社, 1995 LiuChaoqun.MultigridmethodwithApplicationtoComputationalFluidDynamics[M].Beijing,TsinghuaUniversityPress, 1995 (inChinese) [13] 李曉梅, 莫則堯. 多重網(wǎng)格算法綜述[J]. 中國科學(xué)基金, 1996, 10(1): 4-11 LiXiaomei,MoZeyao.ViewpointsofMultigridAlgorithms[J].BulletinofNationalScienceFoundationofChina, 1996, 10(1): 4-11 (inChinese) Researching how Initial Value of Internal Iteration Impacts on Computational Efficiency in Unsteady Flow Solving Gong Yiming, Zhang Weiwei, Liu Yilang (College of Aeronautics, Northwestern Polytechnical University, Xi′an 710072, China) Abstract:On the basis of the dual time stepping method in unsteady flow, we come up with a strategy to improve the efficiency in unsteady flow solving. By means of extrapolating the flow message of a few moments forward, we get the iterative initial value of the next moment to make it closer to the convergent solution and decrease the initial residual value of internal iteration, thus increasing the computational efficiency of unsteady flow. Taylor expand the convective term of flow field in time direction, then we design some kinds of extrapolation schemes. We use the solution of unsteady flow around circular cylinder to verify the computational efficiency of this method, and research on different initial value extrapolation schemes, spatial discretization schemes, and time steps and convergence criteria impacting on the efficiency of initial value extrapolation method. Researches show that based on the dual time stepping method, using the strategy of initial value extrapolation can improve the computational efficiency generally. Of which the alternate extrapolation strategy can doubled the efficiency generally. Compared with upwind scheme, the efficiency of this method increases more when using center scheme; all have obvious effects when using different time steps and convergence criteria. Keywords:computational efficiency, design, dual time stepping method, errors, extrapolation, flow fields, iterative methods, laminar flow, Mach number, mesh generation, Reynolds number, initial value extrapolation, unsteady flow solving, Taylor series 中圖分類號:V211.3 文獻標(biāo)志碼:A 文章編號:1000-2758(2016)01-0011-07 作者簡介:貢伊明(1993—),西北工業(yè)大學(xué)碩士研究生,主要從事計算流體力學(xué)研究。 收稿日期:2015-09-18