宋喜慶 張少雄 胡豐梁

(武漢理工大學(xué)交通學(xué)院1) 武漢 430063) (中國(guó)船級(jí)社2) 北京 100007)

0 引 言

超大型礦砂船(VLOC)與普通散貨船相比在船體結(jié)構(gòu)型式上存在較大差異，近年來礦砂船諸多的傾覆事故及其不斷大型化的發(fā)展勢(shì)頭使其結(jié)構(gòu)安全性備受關(guān)注，一般需采用船舶直接計(jì)算方法進(jìn)行強(qiáng)度校核.載荷的計(jì)算是決定直接計(jì)算準(zhǔn)確與否的關(guān)鍵，礦砂貨物載荷是其重要組成.

干散貨對(duì)于壁面的作用的研究屬于土力學(xué)范疇，擋土墻壓力理論其重要參考.在船舶工程領(lǐng)域， Koichi等[1-2]多次開展實(shí)船航行過程中貨物壓力進(jìn)行監(jiān)測(cè)，然而數(shù)據(jù)顯示礦砂船航行過程中受船體的振動(dòng)、晃蕩等因素的影響，貨物壓力的實(shí)測(cè)結(jié)果與按庫(kù)倫土理論整理的貨物壓力存在較大差異，且無明顯規(guī)律.

目前船舶工程中所應(yīng)用的載荷計(jì)算公式通常以規(guī)范公式為指導(dǎo)，如文獻(xiàn)[3]中對(duì)于干散貨貨物載荷的計(jì)算要求，通過側(cè)壓力系數(shù)計(jì)算非水平板的側(cè)向壓力.這種用于船舶強(qiáng)度直接計(jì)算的內(nèi)部載荷計(jì)算方法不能考慮干散貨顆粒的形狀參數(shù)、艙壁的傾斜程度、滑裂面假設(shè)[4]等，也忽略了干散貨顆粒裝載過程中的隨機(jī)分布、顆粒間的摩擦系數(shù)、顆粒與艙壁的摩擦系數(shù)及顆粒的剪切模量，是一種近似的干散貨內(nèi)部載荷的描述方法.因此有必要對(duì)船舶直接計(jì)算中傳統(tǒng)計(jì)算方法的可靠性進(jìn)行重新考慮.

1 考慮方法

干散貨尤其是礦砂類貨物的的力學(xué)特征體現(xiàn)為“散”“動(dòng)”，使用連續(xù)體力學(xué)分析干散貨壓力問題依賴了太多假設(shè).通過離散元法(DEM)描述干散貨堆積狀態(tài)可以使貨物載荷與實(shí)際更加貼合.隨著計(jì)算機(jī)運(yùn)算性能的快速提升，DEM計(jì)算方法的可行性將逐漸顯現(xiàn)出來.本文基于離散元分析程序EDEM對(duì)這種貨物載荷的分析方法展開探索.

EDEM是一款具有高工程應(yīng)用價(jià)值的離散元分析軟件[5]，它可以對(duì)用戶自定義形狀的顆粒添加力學(xué)參數(shù)并設(shè)置顆粒與顆粒間及顆粒與幾何體間的接觸關(guān)系，同時(shí)可仿真顆粒隨預(yù)定義幾何體的勻速和變速運(yùn)動(dòng)，在降低實(shí)驗(yàn)成本的前提下能夠?yàn)閷?shí)際工程問題提供更好的解決方案.本文采用MSC.Patran建立250 000 DWT VLOC三維有限元模型，提取貨艙的結(jié)構(gòu)有限元網(wǎng)格，通過網(wǎng)格處理程序，為EDEM分析提供網(wǎng)格模型，在EDEM中模擬干散貨的裝載情況，分析貨物載荷產(chǎn)生的堆積壓力場(chǎng).最后通過基于python語言開發(fā)的EDEM與Patran的數(shù)據(jù)接口，實(shí)現(xiàn)離散元壓力場(chǎng)向有限元計(jì)算模型的映射.

2 實(shí)現(xiàn)方法

選取服務(wù)于澳大利亞至中國(guó)的250K DWT VLOC為目標(biāo)船.該船總長(zhǎng)315.29 m、型寬57.00 m，型深25.00 m，設(shè)計(jì)吃水18.00 m，貨艙體積19 544 m3，方形系數(shù)0.863，最大貨物密度2.5 t/m3，貨艙結(jié)構(gòu)模型見圖1.

圖1 貨艙結(jié)構(gòu)模型

2.1 計(jì)算流程設(shè)計(jì)

船舶直接計(jì)算是一項(xiàng)復(fù)雜、高難度的船舶結(jié)構(gòu)安全性評(píng)估方法，涉及多項(xiàng)技術(shù)的綜合應(yīng)用[6]，同時(shí)載荷形式多樣，為保證不對(duì)現(xiàn)有直接計(jì)算流程產(chǎn)生影響，離散元分析得到的靜壓力加載過程設(shè)置在建立船舶有限元計(jì)算模型并完成網(wǎng)格附屬性之后和靜平衡調(diào)整之前完成.

離散元分析程序EDEM和有限元分析工具Patran具有完全不同的結(jié)構(gòu)網(wǎng)格和仿真環(huán)境.其中EDEM僅允許三角形網(wǎng)格，承受載荷并限制顆粒運(yùn)動(dòng)，而Patran建立的有限元網(wǎng)格則以體現(xiàn)結(jié)構(gòu)力學(xué)特性和傳遞外載荷為主要功能，為提高計(jì)算精度，船舶計(jì)算中有限元網(wǎng)格以四邊形網(wǎng)格為主，并帶有少量三角形網(wǎng)格.結(jié)合貨物壓力及程序間結(jié)構(gòu)網(wǎng)格的功能特點(diǎn)，本文設(shè)計(jì)形成的靜載荷加載流程見圖2.

圖2 離散元分析流程圖

2.2 有限元網(wǎng)格及網(wǎng)格處理

通過離散元仿真實(shí)現(xiàn)靜壓力場(chǎng)計(jì)算首先要由目標(biāo)船的有限元模型提取貨艙模型，該模型為由四邊形和三角形網(wǎng)格組成的五面封閉的網(wǎng)格模型，本文直接采用直接計(jì)算中800 mm的網(wǎng)格尺寸即一個(gè)肋距，提取得到的貨艙有限元網(wǎng)格模型見圖3，包含了內(nèi)底板、貨艙底凳前后側(cè)板、槽型艙壁、底邊艙斜板、內(nèi)舷板及艙口圍板.

圖3 有限元分析模型A

在此基礎(chǔ)上，通過網(wǎng)格處理程序?qū)υＰ臀募M(jìn)行轉(zhuǎn)化形成了用于離散元分析的模型數(shù)據(jù)文件.離散元分析中的結(jié)構(gòu)模型保留了模型A的節(jié)點(diǎn)位置信息，建立的適用于EDEM分析的貨艙模型見圖4.該模型具有7 070個(gè)節(jié)點(diǎn)，14 050個(gè)三角形單元.

圖4 離散元分析模型B

2.3 DEM仿真參數(shù)設(shè)計(jì)

DEM仿真參數(shù)包括計(jì)算模型的選擇、材料參數(shù)、接觸關(guān)系模型、裝載方式等.其中在分析模型的選擇上，本文采用的分析模型[7]為Hertz-Mindlin(No Slip)模型，見圖5.

圖5 球形顆粒碰撞接觸模型

當(dāng)量彈性模量滿足：

當(dāng)量半徑滿足：

切向力Ft由切向重疊量δt和切向剛度St決定，滿足

Ft=-Stδt

式中：G*為當(dāng)量剪切模量.離散顆粒法向力和切向力具有阻尼分量，其中阻尼系數(shù)和恢復(fù)系數(shù)相關(guān)，切向摩擦力遵守庫(kù)倫摩擦定律，滾動(dòng)摩擦力通過接觸獨(dú)立定向恒轉(zhuǎn)矩模型實(shí)現(xiàn).

材料參數(shù)和接觸關(guān)系模型的設(shè)置中，本文根據(jù)實(shí)船裝載特點(diǎn)及裝載的貨物類型，材料參數(shù)根據(jù)實(shí)際裝載貨物類型進(jìn)行設(shè)置，采用CSR規(guī)定的休止角和材料特性設(shè)置接觸參數(shù)，同時(shí)，為提高仿真效率，鋼質(zhì)圍壁與鐵礦砂的滾動(dòng)系數(shù)設(shè)置為0.由此建立表1～2的材料參數(shù)及接觸關(guān)系模型.

表1 材料參數(shù)

表2 接觸關(guān)系

在裝載方式上[8]，采用在以艙口中心為原點(diǎn)的6 m邊長(zhǎng)的8邊形范圍內(nèi)隨機(jī)生成顆粒并在重力加速度作用下裝載堆積的方式.由此建立的貨艙重載裝載模型見圖6.

圖6 EDEM裝載模型

2.4 壓力場(chǎng)映射

離散元分析中考慮了顆粒的彈塑性及運(yùn)動(dòng)、碰撞，在顆粒裝載完畢后通過總體坐標(biāo)系下三個(gè)方向上的合力隨時(shí)間的變化曲線可取出相對(duì)穩(wěn)定態(tài)下的壓力場(chǎng)分布.裝載量為93.4萬個(gè)顆粒的貨艙裝載模型X方向合力隨時(shí)間的變化曲線見圖7.

圖7 貨艙X方向合力隨時(shí)間的變化曲線

離散元分析中干散貨顆粒與三角形網(wǎng)格在接觸平衡模型可求得三角網(wǎng)格節(jié)點(diǎn)上的集中力.由于輸出格式及節(jié)點(diǎn)重編號(hào)的問題，網(wǎng)格節(jié)點(diǎn)上的集中力難以直接用于船體強(qiáng)度的直接計(jì)算.利用python語言開發(fā)的EDEM與Patran的數(shù)據(jù)接口，可從EDEM數(shù)據(jù)庫(kù)讀取并完成節(jié)點(diǎn)位置及載荷信息的篩選，在將有效信息加工為Patran的命令流字符串后輸出創(chuàng)建PCL腳本文件[9-10]，以python的OS程序模塊建立程序控制關(guān)系，使節(jié)點(diǎn)位置及壓力場(chǎng)在模型間建立映射關(guān)系，實(shí)現(xiàn)離散元分析得到的靜載荷在Patran中的程序化快速施加.其中，由于模型A與模型B節(jié)點(diǎn)編號(hào)的改變，程序?qū)χ亟ü?jié)點(diǎn)的位置關(guān)系進(jìn)行搜索判斷，消除舊節(jié)點(diǎn)，將重建節(jié)點(diǎn)上的載荷映射到原節(jié)點(diǎn)所在位置上，保留了原網(wǎng)格的形狀及質(zhì)量.詳細(xì)實(shí)現(xiàn)流程見圖8.

圖8 壓力映射程序流程設(shè)計(jì)

3 建立對(duì)比條件

為驗(yàn)證通過離散元方法開展貨物載荷計(jì)算的可行性，本文加載至貨艙的貨物載荷按本文第二部分描述的方法完成，并以本節(jié)定義的載荷總量控制模型和載荷輸出形式為標(biāo)準(zhǔn)與CSR法貨物壓力施加方法建立對(duì)比條件.

3.1 載荷總量控制模型

圖9 貨物堆積截面參數(shù)設(shè)置

由此建立堆積上表面的高度曲面方程為

將Hsur代入得到側(cè)壓力方程

根據(jù)確定的貨物上表面，由此建立控制關(guān)系模型為

式中：N為顆粒的裝載數(shù)量；P為球形顆粒的孔隙率；ρ為初始設(shè)定的貨物密度；Mt為貨物裝載總重；Fz為貨艙重載工況下貨物的垂向總壓力.

3.2 載荷形式

DEM法通過節(jié)點(diǎn)力的方式考慮顆粒與結(jié)構(gòu)網(wǎng)格的接觸關(guān)系，輸出結(jié)構(gòu)模型各節(jié)點(diǎn)上壓力、摩擦力的合力在總體坐標(biāo)系三個(gè)坐標(biāo)軸方向上的等效分量.CSR將貨物靜載荷分解為兩部分載荷，即垂直于貨艙圍壁的壓力、傾斜板沿板面向下的剪切力.兩種計(jì)算方法下貨物靜載荷的體現(xiàn)見圖10.

圖10 兩種計(jì)算方法下貨物靜載荷

本文將CSR計(jì)算得到的兩部分載荷按總體坐標(biāo)系三個(gè)方向等效分解到節(jié)點(diǎn)上，與DEM法得到的節(jié)點(diǎn)力進(jìn)行對(duì)比.

4 分析結(jié)果

過程中選取半徑0.95，0.85，0.75，0.65，0.55，0.45，0.35，0.25和0.15 m的顆粒按2 000個(gè)/s的裝填速度分別開展了裝載實(shí)驗(yàn)，提取得到相對(duì)穩(wěn)定態(tài)下作用于整個(gè)貨艙垂向合力隨顆粒半徑的變化關(guān)系見圖11.

圖11 貨艙Z向合力隨顆粒半徑的變化關(guān)系

由圖11可知，當(dāng)半徑達(dá)到0.15 m，網(wǎng)徑比(網(wǎng)格尺寸/半徑)達(dá)到5.3時(shí)，顆粒與單個(gè)網(wǎng)格最大接觸數(shù)量達(dá)到8.4個(gè)，DEM法計(jì)算得到的總體載荷與CSR法計(jì)算得到的總體載荷已基本一致，差距僅為0.36%.在此條件下形成的重載狀態(tài)下CSR法和DEM法貨物靜載荷等效為節(jié)點(diǎn)載荷后的壓力云圖見圖12.

圖12 兩種計(jì)算方法下貨物靜載荷等效節(jié)點(diǎn)力云圖

由圖12可知，貨物顆粒在內(nèi)底板等構(gòu)件上出現(xiàn)了明顯的隨機(jī)壓力分布，與CSR法中貨物載荷理想化的均勻的載荷分布有較大出入，這與實(shí)際礦砂裝載過程中因貨物顆粒間的隨機(jī)性碰撞后的堆積行為相符.

在整體上，DEM計(jì)算方法在作用于貨艙結(jié)構(gòu)網(wǎng)格上的X,Z方向合力與CSR法差距較小，Y向合力雖然與CSR法有量級(jí)上的差距，但與X方向合力同量級(jí)，這與粒徑半徑的選擇相關(guān).貨艙重載工況下三個(gè)方向的合力表現(xiàn)見表3.

在該離散元仿真參數(shù)下，本文提取了半徑在0.15～0.95 m范圍內(nèi)，9種顆粒半徑裝填完成在內(nèi)舷板上產(chǎn)生的側(cè)壓力隨高度的變化曲線，并與CSR法進(jìn)行對(duì)比，離散元分析方法得到的側(cè)向壓力隨裝填顆粒半徑的縮小，在內(nèi)舷板上的側(cè)向壓力逐漸增大，側(cè)壓力變化逐漸平滑,見圖13.

表3 兩種計(jì)算方法下貨物載荷合力對(duì)比

圖13 內(nèi)舷板Y方向沿高度方向壓力變化

5 結(jié) 論

1) 本文開展了離散元與有限元計(jì)算模型結(jié)構(gòu)網(wǎng)格間壓力場(chǎng)映射方法研究，建立了基于Python語言和Patran PCL命令流的壓力場(chǎng)映射方法的可行的流程設(shè)計(jì)，可在不改變?cè)Ｐ凸?jié)點(diǎn)位置及網(wǎng)格質(zhì)量的情況下，可通過場(chǎng)函數(shù)的形式進(jìn)行大量節(jié)點(diǎn)力的高效加載.

2) 從離散元分析結(jié)果看出，使用網(wǎng)徑比達(dá)到5.3的分析模型進(jìn)行散貨裝載仿真還是會(huì)在內(nèi)底板等構(gòu)件上出現(xiàn)壓力場(chǎng)的概率性分布，這與CSR理想化的規(guī)律的載荷分布模型不符.

3) 本文建立了一套基于離散元的貨物壓力計(jì)算流程，并以實(shí)際貨艙模型為例進(jìn)行了裝載測(cè)試，為貨物載荷的研究提供了思路.