,2
1.上海交通大學(xué) 船舶海洋與建筑工程學(xué)院 上海 200030 2.上海交通大學(xué) 海洋工程國家重點實驗室 上海 200030
自從Hess和Smith[1]用一階面元法求解三維勢流問題以來,面元法作為一種有效的數(shù)值方法已被廣泛地應(yīng)用于求解空氣動力學(xué)和水動力學(xué)中的三維勢流問題。最初的面元法中存在幾個缺點,如物面用平面四邊形面元離散,相鄰面元間存在間隙;在離散大曲率的曲面時必須要有大量的面元,為此需要大量的計算時間;曲面上的源(匯)強度在面元邊界上不連續(xù)。
為了克服低階面元法的不足,人們引入線性或二次的函數(shù)表達面元,并在面元上分布一階或二階多項式源項。近年來,由于三次樣條、B樣條以及非均勻有理B樣條(non-uniform rational B-splines, NURBS)在工業(yè)產(chǎn)品外形設(shè)計中的廣泛應(yīng)用,人們已經(jīng)把他們應(yīng)用到高階面元法中以表達物面幾何和速度勢[2-5]。對于高階面元法中的Rankine源項,Maniar和Danmeier以及Lee等[2-4]均運用三角形細分和級數(shù)展開的方式來近似求解。該方法在計算中較為復(fù)雜。趙成璧[5]的方法是離散第二類Fredholm積分方程,再直接高斯積分。該方法在場點離面元較遠時可獲得較好的結(jié)果,而當(dāng)場點離面元較近時,數(shù)值積分的精度將受到很大的影響。
在處理基于B樣條的高階面元法中的奇異積分時,根據(jù)場點到面元的距離把曲面積分分為遠場積分和近場積分。遠場積分中,直接應(yīng)用高斯-勒讓德公式;近場積分中,采用Telles等[6-7]介紹的方法,先對積分函數(shù)進行三次多項式變換,然后應(yīng)用高斯-勒讓德公式。
考慮物體在無限靜止流體中以勻速U向前運動,定義一個固定在物體上的右手直角坐標(biāo)系o-xyz(圖1)。假設(shè)流體是不可壓縮理想流體,流動無旋,則存在定常擾動勢φ滿足下列邊值問題:
(1)
圖1 坐標(biāo)
(2)
式中:φ(q)——物面上源點q(ξ,η,ζ)處的源(匯)強度;
S——物面。
把式(2)代入式(1) 中得到確定源(匯)強度的第二類Fredholm積分方程。
(3)
有
(4)
式中:u0、v0——代表NURBS曲面中的兩個參數(shù)方向;
dij——物面控制點網(wǎng)格;
m、n——u0、v0方向上物面控制點的數(shù)目;
Ni、g(u0)、Nj、h(v0)——定義在u0、v0上的B樣條基函數(shù);
wij——權(quán)因子;
g、h——u0、v0上B樣條基函數(shù)的階數(shù)。
(5)
式中:u、v——B樣條曲面中的兩個參數(shù)方向;
m′、n′——u、v方向上控制點的數(shù)目;
Ni,g(u)、Nj,h(v)——B樣條基函數(shù);
g、h——u,v上B樣條基函數(shù)的階數(shù)。
把式(4)、(5)代入式(3)中,可以得到用B樣條參數(shù)表達的積分方程:
(6)
(7)
圖2 場點與面元的關(guān)系
在結(jié)合高斯-勒讓德公式積分時,三次多項式變換在不對面元進行任何細分的情況下,可以使高斯點自動向奇異點靠近,這個特性將有助于提高準奇異積分的精確性。因此,對于近場積分,先把式(6)中的積分函數(shù)進行三次多項式變換,然后再結(jié)合高斯-勒讓德公式積分。
為了驗證三次多項式變換法的效果,考慮以下準奇異的二重積分
(8)
其中a,b為常數(shù)。在此使用8×8的高斯-勒讓德公式,分別采用直接積分和經(jīng)三次多項式變換后的積分。表1給出了計算結(jié)果和解析解的比較。
表1 兩種奇異積分計算結(jié)果對比
從表中可以看出,當(dāng)場點在(1.40,0)處時,由于它離面元相對較遠,直接積分法和多項式變換法的計算誤差都較小,但隨著場點離面元距離的減少,當(dāng)減少到(1.05, 0)處時直接積分法的誤差迅速增大到50.80%,而多項式變換法的誤差僅為8.66%,說明后者有著比前者更好的計算精度。
考慮一單位半徑的球體以單位速度沿x軸正向運動,分別采用直接積分法和文中所介紹的方法對定常擾動勢進行了計算。
表2給出了計算結(jié)果。從表2中可以看出,兩種方法的計算結(jié)果均和解析解比較吻合,但本文方法的計算誤差更小,尤其在靠近x=-1.0處。其中一個重要原因是在球的NURBS表達中,靠近x軸兩端的面元較小,場點離面元相對較近,積分的奇異性表現(xiàn)得非常明顯,直接積分將忽略這種特性,導(dǎo)致計算誤差,而本文的計算方法通過對積分函數(shù)進行三次多項式變換,提高了計算精度。
表2 球面速度計算結(jié)果的對比
采用文中介紹的方法,對一系列以單位速度平移的橢球體在第二模態(tài)下的附加質(zhì)量系數(shù)K22=3λ22(4πρabc)進行了計算,其中λ22是附加質(zhì)量,ρ是流體密度,a、b和c分別是橢球體的三個半軸長;在此a∶b∶c的值分別取7∶1∶1, 5∶1∶1, 3∶1∶1, 2∶1∶1和1∶1∶1。根據(jù)三半軸的比值,布置不同的面元數(shù)。計算結(jié)果與Lamb給出的解析解和文獻[7]中的計算結(jié)果(直接高斯積分法)進行了比較,如圖3所示。可以看出,本文的結(jié)果比文獻[7]的計算值更接近于解析解。
圖3 橢球體附加質(zhì)量系數(shù) K22
本方法為解決高階面元法中的奇異積分問題提供了一種有效的手段,在船舶運動性能預(yù)報相關(guān)的m項計算和水動力系數(shù)計算中具有較大的應(yīng)用價值。
[1] Hess J L, Smith A M. Calculation of nonlifting potential flow about arbitrary three-dimensional smooth bodies[J]. Journal of Ship Research, 1964(7):22-44.
[2] Maniar H D. A three dimensional higher order panel method based on B-splines[D], MIT, Massachusetts, 1995.
[3] Danmeier D G. A Higher-Order Panel Method for Larger-Amplitude Simulation of Bodies in Waves[D]. MIT, Massachusetts, 1999.
[4] Lee C H,Newman J N, Solution of Radiation Problems with Exact Geometry, Proc. of 16th International Workshop on Water Wavers and Floating Bodies[C], Hiroshima, Japan, 2001:93-96.
[5] 趙成璧. 船舶曲面設(shè)計現(xiàn)代方法與軟件系統(tǒng)的研究[D].武漢: 武漢交通科技大學(xué), 1999.
[6] Telles J C F,Oliveira R F. Third degree polynomial transformation for boundary element integrals: Further improvement[J]. Engineering Analysis with Boundary Elements, 1994(13):135-141.
[7] 王化明. 基于NURBS高階面元法的有航速船舶輻射問題數(shù)值計算研究[D].武漢:武漢理工大學(xué), 2005.