求解定常勢流正命題的全時(shí)-空變分有限元法

第25卷第2期刊Vol 25 No. 22004年6月CHINESE QUARTERLY OF MECHANICSJune 2004求解定常勢流正命題的全時(shí)空變分有限元法陳池,陶毅,劉高聯(lián)(上海大學(xué)力學(xué)所,上海200072)摘要:非定常流動(dòng)變分原理的建立使得用有限元法來(lái)求解多工況點(diǎn)的設計問(wèn)題成為可能。本文在劉高聯(lián)的非定常變分理論的基礎上,對定常變分問(wèn)題進(jìn)行時(shí)間相關(guān)有限元求解。但由于可壓縮非定常位勢流動(dòng)的控制方程是雙曲型的簡(jiǎn)單地把時(shí)間當作同空間一樣的物理維來(lái)求解是不可行的。而現有的時(shí)空有限元法極其復雜增加了計算復雜度,使其很難用于工程設計中。為此,文[2、3]提出了求解一維非定常問(wèn)題的新型時(shí)空有限元法。本文把該方法推廣到二維流動(dòng)用它求解二維彎管內的流動(dòng)和翼型繞流問(wèn)題。計算結果與用定常方法求得的結果幾乎重合,說(shuō)明該方法可以用于多維時(shí)間相關(guān)求解關(guān)鍵詞:非定常;變分;時(shí)空有限元;計算流體力學(xué)中圖分類(lèi)號:0355文獻標識碼:A文章編號:02540053(2004)021576A New Time/Space FEM for Steady Potential FlowCHEN Chi, TAO Yi, LIU Gao-LianShanghai University, Shanghai Inst. App. Math. Mech, Shanghai 200072, China)Abstract: On the basis of Liu Gaolian's unsteady variational theories, the time-marching method is used tosolve steady variational problems. Because the governing equations for compressible unsteady potentialflow is hyperbolic, looking time dimension as space dimension in the same way is never appropriate. Unfortunately, the existing time/ space FEM is quite complicated, and difficult to be put into engineeringuse. A new time/space FEM was presented to calculate the one-dimensional unsteady problems in papers[2,3]. The method is now extended to calculate two-dimensional steady flow in a time-marching way. Asexamples, one pipe flow and one airfoil rounding flow are considered the results agree well with those ofthe steady flows, which demonstrates the usefulness of this method in solving multidimensional problemsKey words: unsteady variational calculus; time/ space FEM; CFD非定常問(wèn)題是當今研究的熱門(mén)問(wèn)題之一,但到目前為止,多數非定常問(wèn)題都是用有限差分或有限體積法求解的。有限元法以其能適應復雜的幾何區域而逐漸在計算流體力學(xué)領(lǐng)域中得到應用。但是,用有限元方法求解非定常問(wèn)題的文獻也主要是基于 Galerkin變分理論的基礎上進(jìn)行的),對時(shí)間偏導數的處理大多數是采用的差分方法。這一方面是因為非定常的變分原理的建立是極其困難的,因為非定常問(wèn)題的控制方程是雙曲線(xiàn)型的,在時(shí)間方向只能給初始條件,而不能給終止條件。劉高聯(lián)根據 Hamilton原理成功地建立了非定常問(wèn)題的變分原理,為非定常問(wèn)題的變分求解建立了完整周密的理論基礎。另一方面是因為差分法相對于時(shí)空有限元法要容易。但我們注意到,有很多與時(shí)間相關(guān)的問(wèn)題是很難用差分法解決的,例如降落傘的降落問(wèn)題但如果用時(shí)空有限元法就容易處理°。變分法以其有強大的變域變分工具,可以為葉輪機進(jìn)行多工況點(diǎn)反設計而逐漸得到重視。遺憾的是,迄今為止,還沒(méi)有發(fā)現非定常變分計收稿日期:2003409-23基金項目:國家自然科學(xué)基金項目(50136030,10372055)作者簡(jiǎn)介:陳池(1973-),男,四川廣安人,博土力學(xué)季刊第25卷算的有效方法。文獻[7對二維振蕩翼型進(jìn)行了求解,但是文中對時(shí)間維的處理方式完全和空間維的處理方式等同。我們知道,非定常絕熱位勢流動(dòng)的控制方程是雙曲型方程,在時(shí)間方向絕對不能采取同空間方向一樣的處理。為此,我們必須采用時(shí)空有限元法,但是現有的時(shí)空有限元法使用很不方便,公式也很復雜,而且也不適合變分計算。為此,我們需要提出新的求解方法。受雙曲型方程特征的啟發(fā),劉高聯(lián)提出在處理單元內的插值函數時(shí)對時(shí)間方向上采用展開(kāi)處理方式的設想。經(jīng)過(guò)一維的計算試驗,證明是有效的,本文將用它來(lái)對二維的絕熱勢流問(wèn)題進(jìn)行時(shí)間相關(guān)求解。1二維非定常勢流的控制方程及其變分理論維非定常勢流其相應的無(wú)量綱(以滯止參數無(wú)量綱化)氣動(dòng)方程為:8t+·()=0(1)中=A(y-1)(atP(3)上述公式中,、A、、P、y分別為無(wú)量綱密度速度、位勢壓力和絕熱指數。有量綱量的定義見(jiàn)文獻[7]計算區域如圖1所示初邊值條件為(1)在初始邊給定勢函數∮和密度值P作為初始條件,穩定解與初始值關(guān)系不大。即φ=f2P=f3。(2)遠場(chǎng)邊界:應用 Riesman不變量來(lái)計算。(3)物型面上:元=0(4)激波面上:以g表示其法向分速,則有 Rankine-Hugniot激波關(guān)系。+P(A -g)[lA=0|p(An-g)=0圖1計算塊[|A,]=0Fig 1 A slab[|H門(mén)]=(y-1)g[|An](5)尾流線(xiàn)上下:an=0(6)尾緣點(diǎn): Kutta無(wú)載荷條件。用變分推導的系統性途徑可以得到以下的變分原理:上述非定常問(wèn)題的解將使a=0,其中(1-(2∵()+2(1ydadydt+‖(pAn),dtp陽(yáng)dady(φ-f2)-2∫3∮dxc第2期陳池,等:求解定常勢流正命題的全時(shí)空變分有限元法2控制泛函的離散對勢函數在有限單元內作如下的離散:1-[中k(,n,5=22+a1(1+51+§1-m+a2(1+5]+251274+a2(1+1-1[+a4(1+5](6)引入二維空間單元的形函數替代上式中含、?的各項,即有(,,3)=N[+a:(1+5](遵守求和約定)其中N=-41=2,M=9+41=2,N2=9+9+2,M,=99+2定義從空間-時(shí)間坐標(x,y,t)到有限單元體積坐標(、7之間的 Jacobian轉換矩陣如下atdx ay atasas as從而有單元內勢函數對空間-時(shí)間的偏導數8沖沖其中J1為 Jacobian矩陣的逆矩陣,J1、J2、J32分別為由 Jacobian矩陣的逆矩陣的第一、二、三列元素組成的行向量。v綽為勢函數對體積坐標的散度,分別為a285[+a,(1+5]a7a71[+a,(1+5]a=Na,(遵守求和約定),1=12,4(下面對六k同)。由此我們可得到單元內泛函對任意a;的偏導數為e(1+5e[+a(1+5](1+5)JAAn2[+a1(1+5)/、、、1+5+J;18N,(1+)[+a/(1+5(1+J 2anD[+a1(1+]·J2aNo(1+5)NN力學(xué)季刊第25卷+‖(aA,)N(1+5)]·E3dd+‖pN(1+5·E1dln(9)上列諸式中J|為 Jacobian矩陣的行列式,E12、E13為面積分的轉換系數為了把(9)式中含a的項線(xiàn)性化處理我們把和中的a提出來(lái),從而可以得到單元內寫(xiě)成矩陣形式的代數方程組:b2b2 ba b2|a2p2(10)b3 b32 b3 b3abb ba b其中的參數為(1+)aN(1+y)(1+ya(1+5aN(1+5+J,iaN(1+)aN(1+)PJ3 aN(1+5)+j,aNJ,1(1+5NN(1+)(1+JIdedrds0N+‖(a,)N(1+5]·E2d5d3+‖pN(1+t·E1dηp式中對k有求和約定合并各單元的上述方程可以得到整個(gè)網(wǎng)格上的關(guān)于a的方程。初始密度假設為來(lái)流的密度,初始勢函數為來(lái)流速度與當地橫坐標的乘積,即:=A2。由于我們的計算方法是隱式的,可以取較大的時(shí)間步長(cháng)。另外,為了消除上述線(xiàn)性方程組的奇異性,我們必須給定某些點(diǎn)上的a值。由此我們以特征線(xiàn)(或面)公式得到周?chē)吔鐁上的速度和密度,從而可以算出當地的值。應用上述變分泛函求解,不需要其它邊界條件,因為在泛函中我們已經(jīng)考慮了所有的邊界條件3數值結果3.1一維問(wèn)題求解為了驗證這種方法我們首先對一維非定常的流動(dòng)進(jìn)行了求解。我們的例子選自文[5],即要計算如圖2中的大容器當閥門(mén)打開(kāi)后管內的流動(dòng)狀況。用本文的方法所得的結果和文獻[5]用特征線(xiàn)所求得的結果見(jiàn)圖3和4二者幾乎完全重合,說(shuō)明我們的方法用于計算一維非定常流動(dòng)是很成功的3.2二維非定常求解用上述方法我們對二維的可壓縮繞流進(jìn)行了非定常計算。計算的處理方式同一維問(wèn)題求解的過(guò)程。第2期陳池,等:求解定常勢流正命題的全時(shí)空變分有限元法161計算管沒(méi)圖2計算管段及相關(guān)參數Fig 2 Calculated pipe and related parameters10000100000圖3壓力比較(用符號衰示的為文[5]的結果圖4速度比較(用符號示的為文[5]的結果Fig 3 Pressure comparison (results denotedFig 4 Velocity comparison(results denotedby character taken from paper [5])by character taken from paper [5])在這里我們做了兩個(gè)算例。其中一個(gè)算例選自NACA的研究報告,文中是用流函數和勢函數計算的定常流動(dòng)。計算區域及網(wǎng)格如圖5所示。進(jìn)口無(wú)量綱速度是0.3993。我們的計算經(jīng)過(guò)12個(gè)周期達到穩定。得到的結果( Time-Marching)和研究報告的結果( Stream-Potential)見(jiàn)圖6和圖7所示?？梢钥闯鏊鼈兪俏呛系孟喈敽玫?。勢流函數法“魯一喜一著(zhù)一著(zhù)一一看圖5二維管道網(wǎng)格圖6下管壁面速度Fig 6 Velocity on the lower wall of pipe另一個(gè)算例是計算的NACA0012翼型的繞流計算區域和網(wǎng)格見(jiàn)圖8。計算經(jīng)過(guò)40多個(gè)時(shí)間步(單步步長(cháng)為10個(gè)無(wú)量綱時(shí)間)計算基本達到穩定(勢函數絕對誤差總和小于105),所得到的結果如圖9所示。為了比較我們用定常的公式計算了同樣的問(wèn)題把它們的結果畫(huà)在了一起?？梢悦黠@看出,本方法所得的結果和定常計算的結果吻合得相當好。實(shí)際的計算顯示它們之間勢函數差值的總和在10以?xún)?。力學(xué)季刊第25卷特流函數法圖7上管壁面速度圖8翼型網(wǎng)格Fig. 7 Velocity on the upper wall of pipeFig8 Airfoil grid勢流函數法圖9壓力系數Fig 9 Pressure coefficient4結論針對非定常絕熱流動(dòng)的控制方程的雙曲型特點(diǎn),本文采用時(shí)空變分有限元來(lái)計算定常的變分問(wèn)題。從所得到的結果可以看出,本方法用于計算二維時(shí)間相關(guān)的定常流動(dòng)是切實(shí)可行的。通過(guò)上面的推導過(guò)程還可以看出,本方法可以推廣應用于求解多維問(wèn)題。參考文獻:[11 Liu G L. A vital innovation in Hamilton principle and its extension to initial-value problems[C]. Proc 4th Intl Conf Nonlinear MechShanghai Univ Press, 2002, 90-97.[2]劉高聯(lián)陶毅,李孝偉陳池新型全時(shí)空有限元法[C]水動(dòng)力學(xué)學(xué)術(shù)會(huì )議論文集深圳,2003,[3]陶毅劉高聯(lián)李孝偉陳池計算非定常流動(dòng)的全時(shí)空有限元法[C]水動(dòng)力學(xué)學(xué)術(shù)會(huì )議論文集深圳,2003[4 Chung T J. Finite element analysis in fluid dynamics[M]. McGraw-Hill Int Book Company, 1978[5 Hirsch C. Numerical computation of internal external flows, vol 2[M). John-wiley, New York, 1990[6] Hughes T J R, Hulbert G M. Space-time finite element methods for elastodynamics, Formulations and error estimates[J]. Comp Meth[7]郭加宏劉高聯(lián)二維振蕩機翼跨聲速非定常繞流的變域變分有限元解空氣動(dòng)力學(xué)學(xué)報[1996,14(4);217-222[8]劉高聯(lián)輛合熱彈性動(dòng)力學(xué)的統一的變分原理族[].力學(xué)學(xué)報,1999,31(2):165172