王 澄，張高峰，袁昌耀，呂曉勇，楊秀萍，王收軍

（1.天津理工大學(xué) 機(jī)械工程學(xué)院天津市先進(jìn)機(jī)電系統(tǒng)設(shè)計(jì)與智能控制重點(diǎn)實(shí)驗(yàn)室，天津 300384；2.南京晨光集團(tuán)有限責(zé)任公司，江蘇 南京 210006）

在全球經(jīng)濟(jì)一體化持續(xù)發(fā)展過程中，各國(guó)之間的貿(mào)易往來需求變得越來越大，這促使船舶運(yùn)輸往大型化、高效化發(fā)展.不管是從船舶能效還是航行安全性考慮，對(duì)影響船舶推進(jìn)效率、穩(wěn)定性、操縱性能的風(fēng)載荷展開研究，都具有實(shí)際的工程意義.

目前，關(guān)于船舶風(fēng)載荷的研究方法有：理論分析、風(fēng)洞實(shí)驗(yàn)、計(jì)算流體力學(xué)（CFD）法.其中，CFD數(shù)值模擬相比于其他方法，具有發(fā)展成熟、費(fèi)用少、周期短、通用性好等優(yōu)點(diǎn)，也常被用于模擬計(jì)算結(jié)構(gòu)的外流場(chǎng)情況.謝華平等[1]采用RNG k-ε湍流模型，研究不同工況下?lián)躏L(fēng)系數(shù)、長(zhǎng)寬比等因素對(duì)角鋼格構(gòu)塔的風(fēng)載阻力系數(shù)影響；朱小海等[2]對(duì)港口起重機(jī)進(jìn)行數(shù)值模擬，將流場(chǎng)各物理量計(jì)算結(jié)果與行業(yè)設(shè)計(jì)規(guī)范對(duì)比，驗(yàn)證數(shù)值模擬方法的可靠性；林一等[3]基于CFD方法對(duì)自升式平臺(tái)所受風(fēng)載荷進(jìn)行數(shù)值模擬計(jì)算，結(jié)果與風(fēng)洞試驗(yàn)相近，相比之下，規(guī)范計(jì)算結(jié)果則較為保守；倪歆韻等[4]采用SST k-ω湍流模型，評(píng)估兩種特殊作業(yè)船所受環(huán)境載荷，模擬結(jié)果與試驗(yàn)結(jié)果比對(duì)吻合較好；賈寶柱等[5]采用Standard kε湍流模型，計(jì)算超大型油輪風(fēng)載荷系數(shù)并與Isherwood公式估算結(jié)果比對(duì)，分析明顯偏差存在的原因.

一般的干散貨船都是單甲板船，所運(yùn)輸?shù)呢浳锲贩N單一，并且不需要被包裝成捆，不怕擠壓，裝卸便宜.本文以干散貨船為例，基于CFD方法對(duì)船舶周圍風(fēng)場(chǎng)進(jìn)行數(shù)值模擬研究，根據(jù)干散貨船型線圖建立三維模型，建立內(nèi)部、外部計(jì)算域并采用混合網(wǎng)格劃分，計(jì)算平均風(fēng)和梯度風(fēng)兩種風(fēng)況下，船舶處于不同風(fēng)舷角下的風(fēng)載荷系數(shù)，并將數(shù)值計(jì)算結(jié)果與Isherwood公式計(jì)算結(jié)果進(jìn)行比對(duì)，分析數(shù)值計(jì)算方法預(yù)估船舶風(fēng)載荷的可行性.

1 風(fēng)載荷系數(shù)與經(jīng)驗(yàn)公式

1.1 船舶風(fēng)載荷系數(shù)

船舶受風(fēng)載荷的影響主要表現(xiàn)在縱向力、橫向力和艏搖力矩，如圖1，為便于在不同類型、尺寸的船舶間比較，常采用無量綱系數(shù)法來處理這三個(gè)載荷，這些無量綱系數(shù)的數(shù)值大小反映風(fēng)場(chǎng)對(duì)船舶的影響程度，表達(dá)式為[6]：

式中：Cx為縱向風(fēng)載荷系數(shù)；Cy為橫向風(fēng)載荷系數(shù)；Cn為艏搖力矩系數(shù)；Fx為風(fēng)載荷縱向分力，F(xiàn)y為風(fēng)載荷橫向分力，Mz為風(fēng)載荷艏搖力矩；ρ為空氣密度；V為計(jì)算風(fēng)速；AT為水線以上船舶正投影面積；AL為水線以上船舶側(cè)投影面積；LOA為船舶總長(zhǎng).

圖1 風(fēng)載荷示意圖Fig.1 Schematic diagram of wind load

1.2 經(jīng)驗(yàn)公式

目前求解風(fēng)載荷系數(shù)的經(jīng)驗(yàn)公式有：Isherwood方法、Hughes方法、Taylor方法等，從中選擇應(yīng)用廣泛、精度較高和適用性好的Isherwood法計(jì)算風(fēng)載荷系數(shù)，作為數(shù)值計(jì)算結(jié)果的參照[7]，表達(dá)式為

式中：A0～A6、B0～B6、C0～C5均為Isherwood常數(shù)；AT、AL、LOA同式（1）、（2）、（3）；B為船舶型寬；S為水線以上部分側(cè)投影面積周長(zhǎng)；C為側(cè)投影面積中心與船首間距；M為側(cè)投影中支柱或桅桿數(shù)目；Ass為上層建筑側(cè)投影面積.

2 數(shù)值模型

2.1 模型建立

干散貨船主要特征參數(shù)尺寸如表1.利用Solid Works繪制船舶水線以上部分結(jié)構(gòu)特征，以1∶100比例建立三維模型，如圖2所示.由于風(fēng)載荷作用在船舶外表面，內(nèi)部結(jié)構(gòu)對(duì)其沒有影響，所以不對(duì)船舶內(nèi)部建模.同時(shí)，忽略船舶甲板上對(duì)風(fēng)載荷計(jì)算影響微小的細(xì)節(jié)特征，比如鏤空欄桿、纜繩鋼索等，以便之后劃分船舶外流場(chǎng)計(jì)算域.

表1 船舶特征參數(shù)尺寸Tab.1 Dimensionsof ship characteristic parameters

2.2 流場(chǎng)多計(jì)算域與網(wǎng)格劃分

通過在ICEM中建立合理的外流場(chǎng)計(jì)算域來模擬船舶周圍風(fēng)場(chǎng)，數(shù)值模擬中的計(jì)算域存在邊界，然而現(xiàn)實(shí)中的風(fēng)場(chǎng)不存在邊界，為了盡可能減小風(fēng)場(chǎng)模型產(chǎn)生的誤差，工程計(jì)算時(shí)數(shù)值風(fēng)洞的阻塞比須低于3%[5]，阻塞比ε計(jì)算公式為.

圖2 船舶水線以上特征模型Fig.2 Characteristic model of shipsabovewaterline

式中：A為船舶正投影于流場(chǎng)入口上的面積；AN為流場(chǎng)入口面積.

根據(jù)船舶模型尺寸，計(jì)算域的長(zhǎng)取10倍船長(zhǎng)，寬取30倍船寬，高取15倍船高，船體距離流場(chǎng)入口為3倍船長(zhǎng)，距離流場(chǎng)出口為6倍船長(zhǎng)，這保證了計(jì)算域中湍流可以完全發(fā)展.外計(jì)算域如圖3所示，同時(shí)，為了便于計(jì)算不同風(fēng)舷角時(shí)的風(fēng)載荷，在船舶周圍建立直徑為2 m的內(nèi)計(jì)算域，內(nèi)計(jì)算域嵌套在外計(jì)算域中，內(nèi)計(jì)算域如圖4所示.在已建立好的流場(chǎng)計(jì)算域中，90°風(fēng)舷角時(shí)，阻塞比最大，為0.7%，滿足工程要求.

圖3 外計(jì)算域示意圖Fig.3 Schematic diagram of external calculation domain

圖4 內(nèi)計(jì)算域示意圖Fig.4 Schematic diagram of the inner computing domain

計(jì)算域網(wǎng)格劃分采用混合網(wǎng)格，外計(jì)算域使用六面體結(jié)構(gòu)網(wǎng)格，如圖5所示，由于船舶特征細(xì)節(jié)多為不規(guī)則表面，因此，內(nèi)計(jì)算域采用四面體非結(jié)構(gòu)網(wǎng)格，如圖6所示，在船舶表面設(shè)置邊界層網(wǎng)格，內(nèi)外計(jì)算域之間設(shè)置interface面來實(shí)現(xiàn)不同網(wǎng)格類型間的數(shù)據(jù)交換，最終生成的計(jì)算網(wǎng)格數(shù)量為557萬，網(wǎng)格整體質(zhì)量?jī)?yōu)良.

圖5 外計(jì)算域剖視圖Fig.5 Sectional view of the outer computing domain

圖6 內(nèi)計(jì)算域剖視圖Fig.6 Sectional view of theinner computing domain

2.3 垂直風(fēng)速剖面

大氣邊界層內(nèi)風(fēng)的水平運(yùn)動(dòng)受地表摩擦力影響，風(fēng)速隨高度增大而增大，對(duì)于近地風(fēng)場(chǎng)風(fēng)速的表示有指數(shù)風(fēng)速剖面和對(duì)數(shù)風(fēng)速剖面，我國(guó)風(fēng)工程結(jié)構(gòu)抗風(fēng)設(shè)計(jì)中常用指數(shù)風(fēng)速剖面，故本文采用指數(shù)風(fēng)速剖面，風(fēng)速呈梯度遞增，也稱作梯度風(fēng)，表達(dá)式為

式中：VH為距地面高度H處的平均風(fēng)速；Vref為參考高度Href處的平均風(fēng)速，參考高度Href取10 m，Vref取10 m/s[8]，考慮到模型按1∶100比例建立，Href代入計(jì)算時(shí)取0.1 m；α為地面粗糙度指數(shù)，取0.12.

2.4 湍流模型與邊界條件

湍流模型選擇基于RANS的Realizable k-ε模型，相對(duì)于其它湍流模型，該模型引入旋轉(zhuǎn)和曲率相關(guān)的內(nèi)容，收斂性好并且收斂時(shí)間短.文獻(xiàn)[6]中選取了多種基于RANS的湍流模型，對(duì)船舶進(jìn)行數(shù)值計(jì)算并比對(duì)結(jié)果，得出從模擬精度和計(jì)算時(shí)程兩方面綜合考慮，Realizable k-ε模型最合適.

流場(chǎng)入口設(shè)為速度進(jìn)口邊界，平均風(fēng)速取10 m/s[9]，梯度風(fēng)速為指數(shù)律形式，利用Fluent用戶自定義函數(shù)UDF功能實(shí)現(xiàn).雷諾數(shù)計(jì)算公式：

式中：ρ為空氣密度，取1.225 kg/m3；L為特征長(zhǎng)度，取2 m；μ為空氣動(dòng)力黏度，取1.8×10-5N·s/m2.湍流強(qiáng)度I=0.16（Re）-1/8=2.7%，湍流尺度l=0.007L=0.14 m[10].流場(chǎng)出口設(shè)為完全發(fā)展出流邊界.左右壁面設(shè)為對(duì)稱邊界，這樣可以減少計(jì)算域邊界帶來的計(jì)算誤差.其它面采用默認(rèn)的無滑移壁面.求解器設(shè)置為穩(wěn)態(tài)計(jì)算的壓力基求解，SIMPLEC算法，湍流動(dòng)能和湍流耗散率方程采用二階迎風(fēng)格式.在監(jiān)視器中設(shè)置，每一步迭代時(shí)計(jì)算船舶的縱向、橫向阻力系數(shù)和力矩系數(shù).

2.5 數(shù)值模擬結(jié)果

現(xiàn)實(shí)風(fēng)場(chǎng)中，自然風(fēng)作用于船身時(shí)的風(fēng)向角是隨機(jī)的不確定的，為了分析不同風(fēng)舷角下流場(chǎng)特性，在0°～180°風(fēng)舷角內(nèi)，每隔10°對(duì)船舶進(jìn)行一次數(shù)值計(jì)算，由于船舶整體具有對(duì)稱性，不再模擬計(jì)算180°～360°風(fēng)舷角內(nèi)的船舶風(fēng)載荷.從兩種風(fēng)場(chǎng)中分別選取0°、40°、90°、150°共4個(gè)風(fēng)向角下的計(jì)算結(jié)果，對(duì)比分析不同風(fēng)速下流場(chǎng)內(nèi)的壓強(qiáng)、速度分布情況.

圖7—圖10為風(fēng)舷角不斷增大的過程，縱向?qū)Ρ戎锌梢姶氨筹L(fēng)面的流場(chǎng)分布不再均勻，流線紊亂并且開始相互纏繞，尤其是在船體艉部上層建筑周圍的流場(chǎng)中，漩渦現(xiàn)象更為明顯；橫向?qū)Ρ雀鱾€(gè)風(fēng)舷角下的不同流場(chǎng)，可以發(fā)現(xiàn)相對(duì)于平均風(fēng)風(fēng)場(chǎng)，梯度風(fēng)風(fēng)場(chǎng)的壓強(qiáng)范圍和流速都略大（流速大小由流線顏色深淺判斷）.

圖7 0°風(fēng)舷角流線和壓強(qiáng)分布Fig.7 Streamline and pressure distribution at 0°

圖8 40°風(fēng)舷角流線和壓強(qiáng)分布Fig.8 Streamline and pressure distribution at 40°

圖9 90°風(fēng)舷角流線和壓強(qiáng)分布Fig.9 Streamlineand pressure distribution at 90°

圖10 150°風(fēng)舷角流線和壓強(qiáng)分布Fig.10 Streamlineand pressure distribution at 150°

3 風(fēng)載荷系數(shù)對(duì)比分析

圖11—圖13為兩種風(fēng)場(chǎng)數(shù)值模擬所得的風(fēng)載荷系數(shù)與基于Isherwood公式計(jì)算結(jié)果的曲線對(duì)比，在部分風(fēng)舷角下三條曲線出現(xiàn)輕微偏差，但整體上三條曲線走勢(shì)一致，在量值大小上契合度較好；縱向風(fēng)載荷系數(shù)曲線中，極值點(diǎn)出現(xiàn)在40°和140°風(fēng)舷角位置，此時(shí)船舶受縱向風(fēng)載荷影響最大；橫向風(fēng)載荷系數(shù)曲線中，極值點(diǎn)位于90°風(fēng)舷角，此時(shí)兩種風(fēng)場(chǎng)結(jié)果偏差最為明顯，由圖9中船體艉部流線變化可看出，船舶背風(fēng)面出現(xiàn)鈍體繞流現(xiàn)象，兩種風(fēng)場(chǎng)中風(fēng)速的不同，致使形成的氣流漩渦尺度大小不同，能量消耗不同，使得兩種風(fēng)場(chǎng)的計(jì)算結(jié)果出現(xiàn)明顯偏差；艏搖力矩系數(shù)曲線中，40°和150°附件三條曲線相差較大，主要是由于大角度繞流的復(fù)雜性，公式計(jì)算時(shí)無法考慮流場(chǎng)分布的復(fù)雜性.而風(fēng)速不同時(shí)，產(chǎn)生的繞流運(yùn)動(dòng)發(fā)展程度不同，對(duì)船舶的影響也將不同.

圖11 縱向風(fēng)載荷系數(shù)Fig.11 Longitudinal wind load factor

圖12 橫向風(fēng)載荷系數(shù)Fig.12 Lateral wind load factor

圖13 艏搖力矩系數(shù)Fig.13 Coefficient of sway moment

4 結(jié)論

本文對(duì)船舶風(fēng)載荷進(jìn)行數(shù)值模擬計(jì)算，包括不同風(fēng)舷角和不同風(fēng)速下的模擬，將數(shù)值計(jì)算結(jié)果和Isherwood公式計(jì)算結(jié)果對(duì)比分析，所得結(jié)論如下.

1）計(jì)算船舶風(fēng)載荷時(shí)，兩種方法所得結(jié)果雖然在部分風(fēng)舷角下存在一定偏差，但在整體走勢(shì)和量值上吻合較好，表明數(shù)值模擬方法應(yīng)用于船舶風(fēng)載荷計(jì)算的可行性，通過數(shù)值模擬方法計(jì)算船舶風(fēng)載荷，可以直觀了解流場(chǎng)中風(fēng)速、壓強(qiáng)分布，進(jìn)一步分析復(fù)雜風(fēng)況對(duì)船舶的影響，提高預(yù)估船舶風(fēng)載荷的精度，這對(duì)于船舶提高能效和安全性具有實(shí)際的工程意義.

2）在縱向和橫向風(fēng)載荷系數(shù)結(jié)果曲線中，相對(duì)于平均風(fēng)風(fēng)場(chǎng)，梯度風(fēng)風(fēng)場(chǎng)與Isherwood公式的契合度更高，采用垂直風(fēng)速剖面數(shù)值模擬船舶周圍流場(chǎng)，更符合實(shí)際風(fēng)況.

3）40°、140°風(fēng)舷角時(shí)，船舶受橫向風(fēng)載荷和艏搖力矩的影響最大，遭遇惡劣天氣時(shí)，建議盡可能避免船舶風(fēng)舷角處于這兩個(gè)角度左右.