李洪濤 楊 勇

(1.海軍工程大學訓練部 武漢 430033)(2.駐上海江南造船(集團)有限公司軍事代表室 上海 201913)

?

面元法中偶極子影響系數(shù)計算分析*

李洪濤1楊 勇2

(1.海軍工程大學訓練部 武漢 430033)(2.駐上海江南造船(集團)有限公司軍事代表室 上海 201913)

分析了不同方法應(yīng)用于偶極子影響系數(shù)的計算,并詳細介紹了Gauss-Bonnet定理在偶極子影響系數(shù)計算中的應(yīng)用。探討了這些計算方法的近似處理及由此造成的誤差,進而以Gauss-Bonnet定理為基礎(chǔ)提出了計算偶極子影響系數(shù)的修改辦法。將此種修改辦法應(yīng)用于面元法的計算中,以×××槳及×××槳為計算對象,分析修改前后的偶極子計算方法對計算結(jié)果精度的影響。從理論上探討了面元法中,有關(guān)偶極子影響系數(shù)計算環(huán)節(jié)的誤差及減小此部分誤差的辦法。

面元法; 偶極子; 影響系數(shù); Gauss-Bonnet定理

Class Number TP206

1 引言

面元法以其準確的數(shù)學和物理模型,較快的計算速度和較高的計算精度,已廣泛應(yīng)用于螺旋槳的水動力性能預(yù)報。目前主要有基于速度和基于速度勢的兩種面元法。以速度勢面元法為例介紹面元法的基本原理,其主要是利用格林公式和拉普拉斯方程基本解的存在,將物面和尾渦面以某種規(guī)律劃分成小面元,在小面元上布置等強度源匯及偶極子,這一過程中以現(xiàn)階段比較普遍的面元劃分方法為例,即劃分成雙曲四邊形面元。以控制點在場點產(chǎn)生的擾動速度勢和流體在物面上法向速度為零的條件建立起離散后針對小面元積分的方程。而這個離散的積分方程每一離散項中,僅與物面以及尾渦面幾何形狀有關(guān)的系數(shù)即被稱之為影響系數(shù),因此在面元法計算分析過程中影響系數(shù)的計算表達是比較重要的環(huán)節(jié)。到目前為止,計算影響系數(shù)主要有兩種方法,一種是Morino[1～2]在1974年、1975年發(fā)展的一套影響系數(shù)計算方法,另一種是Newman[3]在1986年提出來的關(guān)于以立體角概念計算的方法。兩種計算方法在最終的表達形式上沒有太大區(qū)別,同樣兩種方法都有各自的近似成分。在這里主要對偶極子影響系數(shù)的計算表達作詳細的討論分析。

2 不同方法計算影響系數(shù)時的原理分析

2.1 應(yīng)用積分原函數(shù)法計算偶極子影響系數(shù)

偶極子影響系數(shù)可表達為式(1)。

(1)

此類方法是由Morino發(fā)展的一套影響系數(shù)計算方法,基本原理是通過坐標變換將總體坐標中的雙曲四邊形面元對應(yīng)局部坐標的正方形面元,從而在求解的過程中,通過坐標之間的向量關(guān)系,運用積分的基本理論求解得到。在這里,稱它為積分原函數(shù)法,是由于影響系數(shù)表達式在變換后,可改寫為式(2),求解此式時,以類似于原函數(shù)的項表達最終結(jié)果。

(2)

-F(-1,1)+F(-1,-1)

(3)

而其中F(1,1),F(1,-1),F(-1,1),F(-1,-1)的表達式可表示為有關(guān)面元法向、切向以及場點到控制點向徑的關(guān)系式。具體表達見式(4):

(4)

在式(4)中向量R為控制點到場點的向徑,向量t1、t2為面元切向。

在上述對偶極子影響系數(shù)推導過程中,以解析的式子對偶極子影響系數(shù)進行表達。對于式(3)可以從二維泰勒公式展開的角度進行說明(展開到二階即可),因此對于式(3)要滿足的條件即所積分的面元應(yīng)盡量光滑。就螺旋槳表面而言,在面元較小的情況下,這一條件可以近似滿足,但同時要注意,式(3)忽略了高階的泰勒展開項。

2.2 應(yīng)用立體角概念計算偶極子影響系數(shù)

立體角的概念應(yīng)用于影響系數(shù)的計算是由Newman[3]提出的。主要原理下面將進行分析,對于偶極子影響系數(shù)表達式(1)經(jīng)過簡單的數(shù)學變換可以表示為式(5)的形式。

(5)

在對上式處理時,以場點為球心作一個單位球面,那么變形后的表達式可以認為是離散的面元dS在此單位球面上的投影,但要注意到在這里對于投影面元上每一位置處帶有cos(n,R)的符號。在符號全相同的情況下可以認為影響系數(shù)近似表達為立體角。若是這樣的話,對于離散的面元為雙曲四邊形面元時,根據(jù)立體角的計算表達式,影響系數(shù)可表達為式(6):

(6)

式(6)中Ii為面元的各個角點的外角,且Ii是帶有相同符號的,關(guān)于Ii的表達式與運用積分原函數(shù)解法中的表達式類似,即式(4)。在實際應(yīng)用過程中也有用反余弦表示的。但要注意到,滿足式(6)的立體角計算中,面元的邊應(yīng)是過大圓的(大圓即是過球心的圓)圓弧,在離散面元較小的情況下,可近似認為這是可以滿足的。對于各個角點符號不同時的情況,如上述分析,此時不能直接應(yīng)用立體角計算公式,這時的處理,在文獻[4]中有較詳細的討論。具體方法是,在有一個面元角點符號與其他不同時采用式(7)。

·sgn(I2)+π·sgn(I4)

(7)

在有兩個角點符號與其他兩個角點符號不同時也采用式(7)。

按照立體角的計算方法,在離散面元的邊界不是過大圓弧和離散面元的四個角點符號不全相同,這兩種情況下,就只能采取了近似的處理方法,這要求在面元的劃分過程中要把握好劃分的函數(shù)形式,確保較少的面元出現(xiàn)上述情況,并且面元的邊界在一定范圍內(nèi),可近似認為是過大圓弧的。

2.3 應(yīng)用Gauss-Bonnet定理計算偶極子影響系數(shù)

Newman采用立體角計算影響系數(shù)時,曾在文獻[3]中提到過Gauss-Bonnet定理,但未作詳細說明,在這里首先詳細說明Gauss-Bonnet定理。Gauss-Bonnet定理是微分幾何中關(guān)于測地曲率、高斯曲率及曲面內(nèi)角關(guān)系的定理,其數(shù)學表達形式為式(8):

(8)

式中Ω是曲面上的單連通區(qū)域;n為曲面面元邊界線條數(shù);L為組成Ω區(qū)域的分段光滑的閉曲線;αi為面元的第i個頂點的內(nèi)角;K為曲面的高斯曲率;kg為邊界曲線L的測地曲率。

當測地曲率kg=0(此時組成Ω區(qū)域的分段光滑的閉曲線L均為測地線),高斯曲率K=1時,上式變?yōu)槭?9):

(9)

式(9)的意義為:在測地曲率恒為0,高斯曲率恒為常數(shù)1的情況下,曲面上由式(9)各參數(shù)構(gòu)成的單連通區(qū)域Ω的面積等于曲面面元的各頂點內(nèi)角和減去(n-2)π。

對影響系數(shù)式(5)進行分析,前面的處理與立體角法相似,即以場點為球心作一單位球面,則影響系數(shù)表達式為面積dS在以場點為球心的單位球面上的投影面積。由于dS是投影到以場點為球心的單位球面上,一方面由幾何關(guān)系可知,在近似認為面元較小的情況下,投影到單位球面上面元的各條邊均是過單位球球心的大圓弧。由微分幾何學的相關(guān)知識得到,球面上的所有測地線就是過球心的大圓弧的全體(可用數(shù)學方法證明),因而投影到單位球面的面元的各條邊測地曲率kg恒為0。另一方面,由于面元是投影到單位球面上,故曲面面元的高斯曲率K恒為常數(shù)1。綜合這兩方面的分析可知影響系數(shù)的表達式很好地符合了Gauss-Bonnet定理簡化后的表達式(9)的計算要求,因此可利用Gauss-Bonnet定理簡化式(9)來計算偶極子的影響系數(shù)。但同時需要注意到影響系數(shù)表達式中由于cos(n,R)并不恒為正,在計算時同樣需考慮到其符號變化。

由上述分析可知,Gauss-Bonnet定理的計算方法與立體角的計算方法在本質(zhì)上是相同的,不同的是Gauss-Bonnet定理從更基礎(chǔ)的數(shù)學角度,推導出了投影面元面積的計算公式,從而解釋了為什么只有投影的面元邊線過大圓弧時,才能保證計算公式的準確性。

3 面元角點符號不同時的處理

由前面的討論可知,在求解偶極子影響系數(shù)時,若面元角點符號均相同時,運用立體角方法或Gauss-Bonnet定理是沒有太大問題的,而實際計算中會遇到面元角點符號不同的情況,這在文獻[4]中討論過其計算處理方法。但在這里,根據(jù)上述理論的推導,相對于文獻[4]中的處理,主要針對面元角點有一個符號與其他符號不同時作一些修改。處理方法如下。

圖1 面元修改圖

投影面元四個角點的符號有一個與其他角點符號不同時,若沿用式(9)則會不合理。但注意到,由于離散面元在較小的局部范圍內(nèi),若有一個角點符號與其他三個不同,則三個相同符號所組成的三角面元的角點符號仍可認為是相同的。此時處理方法如圖1所示。

當角點1(或角點3)與其他角點的符號不同時,如圖1中左圖所示,此時將四邊形面元分為兩個三角形面元,在這兩個三角形面元中,由2、3、4三個角點所構(gòu)成三角形面元的角點符號可近似認為是相同的,因此可直接應(yīng)用Gauss-Bonnet公式,其表達式為式(10)所示(式中αi不帶符號):

(10)

而另一個三角形面元的表達式則按照近似處理辦法,在文獻[4]有較詳細的說明,其表達式為(式中αi帶有符號):

(11)

同理對于角點2(或角點4)與其他面元角點符號不同時,運用與式(10)、式(11)相同的處理辦法。相對于對四邊形面元作整個的處理,上述的處理方法,保證了面元一部分的偶極子影響系數(shù)的計算結(jié)果在理論上期的準確度有所提高。

對于有兩個角點的符號與其他角點符號不同時,運用了與文獻[4]中計算偶極子影響系數(shù)相同的處理辦法,即與式(7)相同的處理辦法。

以上的計算處理,與文獻[4]中主要的不同在于面元有一個角點符號與其他角點符號不同時,更精細地處理了偶極子影響系數(shù)的計算。從理論上,這樣處理后能夠提高偶極子影響系數(shù)理論上的計算精度。

下面通過將此種修改之后的處理方法應(yīng)用于×××槳、×××槳敞水性能計算,與修改前的計算結(jié)果進行對比,檢驗此種方法的計算精度。具體對比結(jié)果如圖2、圖3所示。在所有的曲線圖中,1為應(yīng)用修改后方法的計算值,2為原計算值,3為試驗值。

圖2 ×××敞水性能比較曲線圖

圖3 ×××敞水性能比較曲線圖

以上對×××槳和×××槳應(yīng)用修改后的方法求解的結(jié)果在某些方面更接近試驗值,×××槳運用修改后的計算方法所得到的結(jié)果與修改前的計算結(jié)果差別不是很大,但在設(shè)計點附近,修改后的計算結(jié)果要更接近試驗值。,×××槳運用修改后的計算方法所得到的計算結(jié)果,在Kq的計算上明顯要比修改前的計算結(jié)果更接近實驗值,同樣在設(shè)計點附近,所得到的計算結(jié)果也是優(yōu)于修改前的計算結(jié)果的。但注意到,對于Kt的計算值,修改后的計算結(jié)果并不比修改前的計算結(jié)果好,而敞水效率的計算結(jié)果也沒有太大變化。初步分析,由于面元法是在離散的基礎(chǔ)上求解迭代的,理論上,網(wǎng)格劃分得越密,從幾何上看就越接近真實,但同時注意到,離散得越密,以公式計算及向量描述的累積誤差則越大。因此合理的劃分網(wǎng)格,同時盡量保證面元投影到單位球面上的邊界線與過球心的大圓弧重合,對于運用Gauss-Bonnet公式求解時的精度控制有著較大的意義,出現(xiàn)上述計算結(jié)果的誤差也是兩者的綜合作用。因此如何選擇合適的網(wǎng)格劃分,對于此類計算影響系數(shù)的方法的應(yīng)用具有很重要的意義,這也是需要進一步分析的問題。

4 結(jié)語

總結(jié)了不同方法計算偶極子影響系數(shù),并簡要分析了這些計算方法中近似處理。詳細介紹了Gauss-Bonnet定理,從微分幾何學的角度,更本質(zhì)地揭露了偶極子影響系數(shù)計算的實質(zhì),因此對于分析計算中,由偶極子影響系數(shù)計算過程產(chǎn)生的誤差,能夠更深刻的把握,從而為提高這一部分的計算精度作了有益的探索。通過理論分析,提出了計算影響系數(shù)修改的方法,從理論上說明了其可行性,并通過對×××槳、×××槳對比計算,詳細分析了結(jié)果差異的原因,從而揭示出了面元法的計算過程中偶極子影響系數(shù)環(huán)節(jié)的誤差來源,進一步明確了面元法在應(yīng)用過程中有關(guān)偶極子影響系數(shù)計算方面需要注意的問題。

[1] Morino L, Ching-Ching Kuo. Subsonic Potential Aerodynamics for Complex Configurations: A General Theory[J]. AIAA Journal,1974,12(2):191-197.

[2] Morino L, Chen L, Suciu E O. Steady and oscillatory Subsonic and Supersonic Aerodynamics around complex Configurations[J]. AIAA Journal,1975,13(3):368-374.

[3] Newman J N. Distribution of sources and normal dipoles over a quadrilateral panel[J]. Journal of Engineering Mathematics,1986,120:113-126.

[4] 譚廷壽.非均勻流場中螺旋槳性能預(yù)報和理論設(shè)計研究[D].武漢:武漢理工大學,2003.

[5] 胡健.螺旋槳空泡性能及低噪聲螺旋槳設(shè)計研究[D].哈爾濱:哈爾濱工程大學,2006.

[6] 王國強,董世湯.船舶螺旋槳理論與應(yīng)用[M].哈爾濱:哈爾濱工程大學出版社,2007.

[7] 徐國良.計算幾何中的幾何偏微分方程方法[M].北京:科學出版社,2010.

[8] 蔣志超,王國強,等.導管面元網(wǎng)格劃分和影響系數(shù)計算的研究[J].船舶力學,2000,4(5):1-8.

[9] Hess J L. Calculation of Potential Flow about Arbitrary Tree Dimensional Lifting Bodies[R]. Technical Report MDC J5679-01, McDonnel Douglas, Oct 1973.

[10] Kinnas S A, Hsin C Y, Keenan D A. Potential Based Panel Method for the Unsteady Flow Around Open and Ducted Propellers[C]//The 18thSymposium on Naval Hydrodynamics,1991.

Calculation and Analysis of Dipole Influence Coefficient in Surface Panel Method

LI Hongtao1YANG Yong2

(1. Administrative Office of Training, Naval University of Engineering, Wuhan 430033)(2. Navy Representative Office in Shanghai Jiangnan Shipyard(Group) Co., Ltd, Shanghai 201913)

Different methods used in the calculation of dipole influence coefficient were analyzed and the Gauss-Bonnet theory was introduced in detail. The approximate treatment applied to these methods was discussed and the error result from the treatment was analyzed. Then a modified calculation about dipole influence coefficient based on Gauss-Bonnet theory was proposed. The modified method was used in surface panel method, and the propellers DTMB××× and DTMB××× were analyzed by the surface panel method. The result was compared to the calculation without modification. In the end, the error induced by the calculation of dipole influence coefficient was discussed.

surface panel method, dipole, influence coefficient, Gauss-Bonnet theory

2014年9月14日,

2014年11月3日

李洪濤,男,工程師,研究方向:裝備綜合保障。

TP206

10.3969/j.issn1672-9730.2015.03.011