賀陽映 ,馬 勇 ,張 松 ,藺世杰
帆翼是帆船帆板前進的動力源,其性能直接影響運動帆船帆板航行的速度(Parolini et al.,2005;Persson et al.,2017;Viola,2011)。帆板運動由沖浪運動演變而來(馬勇 等,2013a;Ma et al.,2016),其比賽的基本航線可以分為梯形外繞航線和梯形內繞航線,整個航線涵蓋了迎風、橫風、順風等多個航段(藺世杰等,2017;馬勇等,2013c),運動員需要進行調帆、換舷等操作來完成超越、繞標等動作(馬勇等,2013a,2013c,2016)。在比賽中,不同航段中的每一次操作都會改變帆翼的攻角,而帆翼攻角的改變會影響帆翼的升力與阻力,因此,需要對不同攻角下帆板帆翼的空氣動力性能進行研究(馬勇等,2016;Ma et al.,2016)。對于奧運會級別的 Neil Pryde RS:X帆板,比賽中帆翼的攻角會主動調整或者受風作用改變,帆翼會發生形變,而帆翼外形的變化又影響帆翼空氣動力性能,也就是帆翼與其周圍流場會發生流固耦合(Fluid-Structure Interaction,FSI)作用(雷曉珊 等,2019;Augier,2012;Bak et al.,2013;Durand et al.,2014),因此,對于帆翼空氣動力性能研究需要考慮FSI。
越來越多的學者利用計算流體動力學(Computation‐al Fluid Dynamic,CFD)的方法研究運動帆翼的氣動特性(Parolini et al.,2005;Viola 2011;Lee et al.,2016;Persson et al.,2017)。在試驗研究和求解Navier-Stokes方程數值模擬方面,Giacobone(2017)對“美洲杯”帆船在不同湍流模型(Spalart-Allmaras、k-ε、SSTk-ω)和不同網格下的升阻力系數、力矩系數進行了計算,并將計算結果與風洞試驗結果進行對比。Nava等(2017)基于雷諾時均方程(Reynolds Averaged Navier-Stokes,RANS)和大渦模擬(Large Eddy Simulation,LES),對迎風條件下速度為7.4m/s、攻角為18.7°的美洲杯帆船前帆和主帆的空氣動力性能進行了數值模擬,所得結果與已發表的試驗數據吻合良好。在Neil Pryde RS:X級別帆板帆翼空氣動力性能研究方面,在不考慮流固耦合情況下,何海峰(2012)和羅曉川(2012)基于RANS方法進行了不同攻角和不同風速下的數值模擬,發現帆翼的失速角為20°~30°之間,隨著攻角的增加,帆面壓力中心向前移動至帆翼前方約1/3處。
隨著計算機水平的發展與提升,考慮流固耦合效應的帆翼空氣動力性能研究越來越多。Lee等(2016)對二維的三面水翼帆0°~360°(間隔為15°)攻角范圍空氣動力性能進行數值模擬,發現 45°、90°和 135°3 個常見視風角下三面水翼帆的最佳攻角范圍。Sacher等(2017)考慮到FSI的IMOCA主帆空氣動力性能數值模擬與試驗較吻合,認為通過高斯方程可以將帆翼基本參數和性能參數進行快速優化。Deparday等(2018)對風角為50°~140°下的三角帆帆翼進行了流固耦合研究,發現帆翼氣動彈性相對較弱,非穩定性振動發生在前帆。雷曉珊等(2019)應用彈簧光順和局部網格重構的動網格技術,進行了同一攻角時不同風速下Neil Pryde RS:X級別帆板帆翼的單/雙向流固耦合數值模擬,得到了帆翼氣動力、表面壓力、變形及應變特征。
從國內外研究進展來看,目前對于運動帆翼的研究主要是不考慮流固耦合時帆翼氣動特性的數值模擬與試驗,而對Neil Pryde RS:X級別帆板帆翼流固耦合的問題也僅討論了風速對于帆翼氣動特性的影響,忽略了帆翼在不同攻角時風致渦激振動對帆翼結構的影響。因此,本文擬通過求解RANSE,基于單向流固耦合,對風速為6m/s時迎風、橫風、順風航段常見攻角下Neil Pryde RS:X帆板帆翼的流場及結構場進行數值模擬,分析攻角變化對帆翼周圍流場及帆翼結構的影響。
研究對象是奧運會級別Neil Pryde RS:X帆板帆翼。首先,利用Pro/ENGINEER Wildfire對測量得到的帆翼外型數據進行逆向簡化建模,忽略索具、加強筋等部分,只考慮桅桿、帆面。帆翼模型如圖1和圖2所示,NP帆板帆翼的基本參數如表1所示。Neil Pryde RS:X帆板帆翼的失速角處于20°~30°之間(何海峰,2012;羅曉川,2012),但是由于所選攻角間隔較大,計算得到的帆板帆翼失速角范圍也過于寬泛,因此,本文對該帆翼在風速為6 m/s、攻角為 20°~50°及80°~100°、間隔5°時的工況進行數值仿真。

圖1 Neil Pryde RS:X帆板帆翼模型正視圖Figure 1.Front View of a Sail Wing Model

圖2 Neil Pryde RS:X帆板帆翼模型俯視圖Figure 2.Top View of a Sail Wing Model

表1 Neil Pryde RS:X帆板帆翼參數Table 1 Parameters of the Sail Wing for Neil Pryde RS:X Class
本文在對帆板帆翼進行單向流固耦合數值模擬時,根據實際需要假定帆翼為:帆面由均質材料組成,帆面材料屬于各項同性材料,帆翼所發生的變形均在其屈服強度之內。帆翼的物理參數:帆面彈性模量為14.4 GPa、泊松比為0.33、密度為100.86 kg/m3,桅桿彈性模量為161.7 GPa、泊松比為0.23、密度為100.86 kg/m3(雷曉珊等,2019)。
在對帆翼進行數值模擬時,為了保證計算域既可以滿足計算精度又不浪費計算資源(馬勇等,2013b),根據帆板帆翼尺寸,計算域和邊界條件為:1)前端(x軸正方向)為12 m,邊界條件為速度入口;2)后端(x軸負方向)為20 m,邊界條件為自由出流;3)側邊界(y軸方向左右)為12 m,邊界條件為對稱平面;4)上邊界(z軸正方向)為12 m,邊界條件為對稱平面;5)下邊界(z軸負方向)為1 m,邊界條件為對稱平面;6)帆翼表面邊界條件為wall壁面,如圖3所示。在數值模擬時為了更加真實地模擬帆翼的運動,將桅桿、桅桿底以及帆翼底部設置為固定約束(圖4)。將流場計算所得帆翼表面壓力加載至結構場的帆翼模型,此處流固耦合面的選擇應與流場中保持一致,施加載荷后帆翼表面如圖5所示。

圖3 計算域及邊界條件Figure 3.Computational Domain and Boundary Conditions

圖4 帆翼約束Figure 4.Constraint of the Sail Wing

圖5 帆翼載荷Figure 5.Load of the Sail Wing
流體域及結構域進行計算前需要對其進行離散,即把原有的計算區域劃分成若干個子計算區域,并且確定區域中的各個節點,進而生成網格。本文在數值模擬過程中需要計算多種攻角,且當風吹過帆翼時,帆翼會產生運動或者發生變形。綜合考慮計算工況及計算過程中模型的變化,本文對計算域及帆翼模型進行非結構化網格劃分。為了更加準確地描述帆翼周圍流場,對帆翼周圍進行適當加密,但由于桅桿的存在導致空氣繞流不順暢,使得整個帆翼背風面產生漩渦(馬勇,2009),因此,在對帆翼進行網格劃分時,對桅桿處的網格進行加密。網格劃分完成后帆翼模型如圖6所示,計算域網格情況如圖7所示。

圖6 帆翼網格分布Figure 6.Mesh of the Sail Wing

圖7 計算域網格Figure 7.Mesh of Computational Domain
帆翼周圍空氣在流動過程中遵循質量守恒、動量守恒以及能量守恒定律,但因其周圍風速遠小于聲速的1/3,認為帆板帆翼周圍的空氣為不可壓縮,且帆翼周圍空氣流動過程中熱交換量較小,所以計算時未考慮能量守恒方程。為了避免計算量大等問題,流體域基于雷諾平均法對質量守恒方程和動量守恒方程中各量進行瞬時和脈動疊加處理,其形式如下(林虹兆,2016):

其中,ui=(u,v,w)是速度在xi=(x,y,z)各方向上的分量,和分別代表時均速度和脈動速度;ρ為流體密度,ν為流體的運動粘度系數,為湍流產生的雷諾應力張量。
增加的雷諾應力項導致雷諾方程組本身不封閉,為了使方程組封閉,本文引入Realizablek-ε模型建立由湍流脈動引起的附加應力與時均應變率之間的聯系,從而使雷諾方程封閉實現求解。Realizablek-ε的湍流動能方程(k方程)和耗散方程(ε方程)的最后形式(Shih et al.,1995):

其中,k為湍流動能,ε為湍流動能耗散率,μt為湍流粘度,μ為湍流粘度系數,Gk為由平均速度梯度而產生的湍流動能項,σk=1.0,σε=1.2,C2=1.9,
升力系數和阻力系數公式如下:

其中,L為升力,D為阻力,ρ為空氣密度,U為風速,S為帆翼面積。
為了驗證數值模擬方法的準確性,在西北工業大學低湍流風洞對縮尺比為1∶15的剛性帆翼進行試驗,帆翼上的力學數據通過ATI nano17系列六分量天平測得。采用Fluent軟件對試驗模型進行數值模擬,計算結果如圖8所示。數值模擬結果與試驗結果對比發現:在風速為5m/s和6m/s、攻角為10°~20°時,試驗得到的帆翼上的升阻力與數值模擬計算結果顯示出良好的一致性,驗證了本文采用的數值方法的可靠性。

圖8 數值模擬與試驗結果對比Figure 8.Comparison of the Numerical Simulation and Experimental Results
升力系數隨攻角的變化分為3個階段,當攻角處于20°~25°范圍內時,升力系數隨著攻角的增加逐漸增大,25°攻角時升力系數達到最大值;當攻角處于25°~90°范圍時,帆翼的升力系數隨著攻角的增大而減小,所以對于女子Neil Pryde RS:X級別帆板帆翼而言,其失速角位于25°附近;當攻角大于90°時,升力系數變為負值,導致升力系數變為反向,此時升力變成使帆板發生橫移的力,隨著攻角的繼續增加,升力系數逐漸增大,其阻礙帆板前進的作用顯著增加(圖9)。

圖9 不同攻角下帆翼升阻力系數Figure 9.Lift and Drag Coefficient at Different Attack Angles
所有航段帆翼背風面均為負壓;迎風及橫風航段攻角從20°增加到50°時,壓力波動較為明顯;攻角為20°時,帆翼背風面帆尾角處以及桅桿上部壓力最大,負壓最低點位于桅桿中下部,形成狹長區域;攻角為25°時,帆尾角處壓力減小,但桅桿上部高壓區域增大,且后帆邊上部壓力也開始增大,桅桿中下部最低負壓區域逐漸減小;攻角為30°~50°時,帆尾角處、桅桿上部以及后帆邊上部壓力減小,但桅桿上部與后帆邊上部所包圍的正壓區域面積逐漸擴大(圖10)。

圖10 不同攻角下帆翼(左)背風面、(右)迎風面壓力云圖Figure 10.Pressure Distribution of the Sail Wing(Left)Leeward Side and(Right)Windward Side at Different Attack Angles
從圖11可以發現,攻角為20°~25°時,流線沿帆翼弦向順著帆翼表面流向后緣,流動順暢且未發生流動分離;經過失速角之后,帆翼背風面流動分離明顯,出現多個回流,并且隨著攻角的增大流線分析區域加大。攻角為40°和100°時,帆翼背風面不僅產生回流,且伴隨明顯的漩渦,流動更加復雜。

圖11 不同攻角下帆翼XY(Z=1.45m)截面上的速度云圖及流線圖Figure 11.Velocity Distribution and Streamlines on XY(Z=1.45 m)Section of the Sail Wing at Different Attack Angles
從圖12中可以看出,不同攻角下帆翼發生形變的位置基本一致,但最大變形值略有差異。綜合來看,最大變形均發生在帆頂角處,這是由于帆頂角處與流體相對運動速度較高,并且到桅桿以及帆翼底面距離最大,所以最容易發生變形,且變形最大(圖13)。

圖12 不同攻角下帆板帆翼變形云圖Figure 12.Deformation Distribution of the Sail Wing at Different Attack Angles

圖13 帆翼最大變形曲線圖Figure 13.Maximum Deformation of the Sail Wing
等效應力表示各個方向上的材料應力差值,對帆翼表面進行等效應力變化分析,可以快速確定模型中受力最大的區域。根據單向流固耦合計算結果,得到不同攻角下柔性帆翼等效應力云圖(圖14)。從圖14中可以看出,帆翼表面的等效應力分布不均勻,存在著應力集中的現象,桅桿頂端與帆上角連接處等效應力的數值最大,應力分布最為集中。此外,帆翼后帆邊也存在較大的等效應力,而桅桿中下部以及帆翼底部等效應力較小。可見,在本文所涉及的攻角范圍內,帆板帆翼的主要受風區域位于桅桿中上部和后帆邊的中下部,以及兩者之間的區域。

圖14 不同攻角下帆板帆翼等效應力云圖Figure 14.Equivalent Stress Distribution of the Sail Wing at Different Attack Angles

圖15 帆翼最大等效應力曲線圖Figure 15.Maximum Equivalent Stress of the Sail Wing
Neil Pryde RS:X帆板比賽有外繞比賽航線(Ma et al.,2016)(圖16),及與外繞航線對應的內繞航線(馬勇等,2013a,2013c),其不同點在于內繞航線在航標1和航標4之間多進行一次繞標。但無論是內繞航線還是外繞航線,在運動過程中都要完成迎風航段、橫風航段、順風航段以及繞標,航線中涉及到的核心操作就是帆翼攻角調整(馬勇等,2013a,2013c)。其中迎風航段、橫風航段以及繞從1標、繞3標向2標航行、繞4標向1標航行時,帆翼攻角都在20°~50°之間;順風航段帆翼攻角在80°~100°之間;而外繞航線繞2標向3標航行時,帆翼攻角則從20°~50°之間轉向80°~100°之間。在不同航段以及繞航標1、2、3、4時,帆船運動員需要根據外界風況(風速、風向角)變化以及對手所在位置調整帆翼攻角,從而使帆板按照自己的意圖航線。

圖16 奧運會比賽外繞航線圖(Ma et al.,2016)Figure 16.Olympic Outer Circuit in the Sailing Regatta
在迎風、橫風航行階段航行時,帆翼迎流攻角調整范圍一般為20°~50°,結合圖9(a)的計算結果,攻角增加使帆翼流動分離導致帆翼壓差阻力迅速增加;而在順風航段中,從圖9(b)結果可知,阻力系數變化曲線較為平緩,升力最大時攻角為25°。從迎風、橫風航行階段屬于典型的升力型航段分析,為獲取更快的航行速度,需要使帆翼上的升力盡可能大而阻力盡量小,因此,攻角選擇為20°~25°。本文考慮單向流固耦合影響并進行了攻角每隔5°的計算,確定最佳攻角為25°左右,這與已有研究最佳攻角在20°~30°之間(何海峰,2012;羅曉川,2012)范圍符合,并進一步精確。所以,帆板運動員在迎風航段、橫風航段以及繞從1標、繞3標向2標航行、繞4標向1標航行時,將帆翼攻角調整到25°左右可以使帆板航速更快。
隨著攻角的增大,帆翼背風面負壓最低點沿桅桿逐漸下移,最低負壓區域也逐漸減小,帆翼迎風面與來流夾角較小。由于桅桿與來流垂直,所以其受壓最大,而后帆邊受壓最小是桅桿和帆翼相互影響造成的(馬勇,2009)。此外,帆翼迎風面桅桿處速度最大,均大于來流速度(6 m/s),后帆邊處速度最小,背風面后帆邊尾流處速度最大,帆翼拱度最大處速度最小,這與其他級別帆板帆翼的空氣動力性能結果類似(馬勇等,2016;Ma et al.,2016)。由于帆翼拱度特點,順風航段帆翼背風面壓力變化較小,正壓區域位于帆翼中部,隨著攻角的增大,正壓中心逐漸向桅桿方向移動;帆翼迎風面與來流幾乎垂直,導致帆翼整個迎風面均受正壓。
帆船帆板運動員在比賽中要隨時注意風況(風速或者風向角)變化,變化后的風速和風向角度可能使得帆翼產生明顯的渦脫現象,從而引起帆翼的變形和振動或者失速,此時,帆船帆板運動員要及時調整帆翼參數,使帆翼處于較優攻角狀態(雷曉珊等,2019)。
有關奧運會Neil Pryde RS:X級別帆板帆翼與其周圍流場的流固耦合作用的結論與其他帆翼流固耦合作用類似(Augier,2012;Bak et al.,2013;Durand et al.,2014)。從結構域來看,帆翼發生形變的位置基本一致,最大變形基本均發生在帆頂角處,且迎風及橫風航段、順風航段帆翼最大變形與帆翼的升力、阻力變化有關(圖9、圖13)。在迎風航段,當攻角小于失速角(25°附近)時,隨著攻角的增加升力迅速增大,而此時阻力較小,因此,在攻角為20°~25°時,最大變形量隨攻角逐漸增大,且增大幅度較為明顯;當攻角達到失速角之后,升力突然減小,阻力繼續增大,此時最大變形迅速減小;當攻角繼續增大,升力緩慢減小,阻力小幅度增大,此時作用于帆翼面上的升力仍大于阻力,所以帆翼最大變形仍繼續增加,但幅度變小。在順風航段,當攻角<90°時,作用于帆翼上的升力逐漸減小,阻力基本無變化,但此時阻力大于升力,所以帆翼的最大變形量逐漸減小;當攻角>90°時,帆翼拱度發生轉向,阻力成為構成帆板前進的推進力的一部分,升力成為阻礙帆板前進的力,此時盡管阻力遠大于升力,但是升力隨攻角不斷增加,所以最大變形與之前相比有明顯增大,但其趨勢為逐漸減小。
此外,帆翼表面應力分布并不均勻,后帆邊處存在較大的應力集中,桅桿頂端與帆上角連接處應力分布最為集中,但桅桿中下部以及帆翼底部等效應力較小(圖14)。由計算結果對比可知,與順風航段相比,迎風及橫風航段各個攻角下帆翼所受最大等效應力值偏小(圖15)。迎風、橫風航段,20°~30°攻角范圍內,最大等效應力值相差較小;攻角為45°時等效應力最大,約為4.5×105 Pa。順風航段下,85°攻角時最大等效應力最小,約為5.2×105 Pa,90°時最大,約為6.3×105 Pa;95°和100°攻角時,最大等效應力幾乎相等。帆船帆板運動員在比賽中要根據風況(風速和風向角)變化,隨時觀察帆翼應力集中位置與大小,及時調整帆翼攻角,避免應力集中與失速導致渦激振動等疊加作用。
1)迎風及橫風航段,帆板帆翼失速角處于25°附近,其最佳攻角為20°左右;在順風航段,阻力系數均隨攻角的增加而增大,應盡量減小攻角。
2)帆板帆翼最大變形發生在帆頂角處;迎風航段當攻角達到失速角之前,最大變形量隨攻角逐漸增大,且增大幅度較為明顯;當攻角達到失速角之后,最大變形迅速減小。順風航段當攻角大于90°時,最大變形與之前相比有明顯增大。
3)帆翼表面存在應力集中的現象,桅桿頂端與帆上角連接處等效應力的數值最大,應力分布最為集中,帆翼后帆邊也存在較大的等效應力,而桅桿中下部以及帆翼底部等效應力較小。
4)帆船帆板比賽中運動員可將迎風航段和順風航段攻角分別調整為20°~25°之間和90°左右,這樣能夠有效避免應力集中與失速導致的帆翼面渦激振動,從而避免桅桿斷裂或者帆面撕裂的情況出現。