考慮貯料-鋼板筒倉相互作用的靜動態(tài)側(cè)壓力模擬研究*

戴士博 靜行，2 郭旭 歐陽德光

1.河南工業(yè)大學土木工程學院 鄭州450001

2.河南省現(xiàn)代綠色生態(tài)倉儲體系國際聯(lián)合實驗室（河南工業(yè)大學） 鄭州450001

引言

糧食鋼板筒倉作為儲存糧食的一種薄壁圓柱殼體結(jié)構(gòu)，具有輕質(zhì)高強、倉儲容量大、易于機械化操作等優(yōu)點，顯示出良好的發(fā)展前景［1］。但這類倉儲結(jié)構(gòu)由于內(nèi)部貯料的質(zhì)量遠大于結(jié)構(gòu)自重，故在鋼板筒倉的設(shè)計中也存在許多特殊問題，如幾何尺寸和支撐形式、貯料力學特性、貯料與倉體復雜的相互作用及不同類型荷載作用的影響［2］。其中，倉內(nèi)貯料對倉壁產(chǎn)生的側(cè)壓力是筒倉結(jié)構(gòu)受力分析與設(shè)計的關(guān)鍵問題［3］，且目前各國筒倉結(jié)構(gòu)設(shè)計規(guī)范對于側(cè)壓力的計算差異較大［4，5］。因此，有必要對貯料-筒倉結(jié)構(gòu)體系的側(cè)壓力分布規(guī)律進行深入研究。

目前，通過建立精細化的有限元數(shù)值模型進行筒倉結(jié)構(gòu)受力分析已成為該類研究的主要手段。張翀等［6］對地震作用下貯料的側(cè)壓力進行了研究，結(jié)合有限元數(shù)值模型提出了一種便于實際工程采用的筒倉側(cè)壓力計算方法。Mehretehran［7］等利用Drucker-Prager 模型模擬倉內(nèi)貯料，對不同高徑比的筒倉模型側(cè)壓力進行模擬，結(jié)果表明，有限元計算結(jié)果與經(jīng)典Janssen 理論預(yù)測的壓力分布吻合較好。楊鴻等［8］采用Drucker-Prager

模型描述貯料的本構(gòu)關(guān)系，建立了考慮貯料與筒倉相互作用的有限元數(shù)值模型，進行了參數(shù)的拓展分析，討論了貯料泊松比、內(nèi)摩擦角、膨脹角及貯料與倉壁摩擦系數(shù)對側(cè)壓力的影響。Herle［9］等認為顆粒材料的應(yīng)力-應(yīng)變關(guān)系可以用亞塑性模型來描述，并對模型參數(shù)的確定方法進行了研究。周長東［10］等基于亞塑性本構(gòu)理論對貯料的靜態(tài)壓力作用進行有限元模擬，提出貯料側(cè)壓力受初始孔隙比、內(nèi)摩擦角、顆粒硬度和摩擦系數(shù)影響較大。

綜上，為準確計算貯料-鋼板筒倉結(jié)構(gòu)體系的靜、動態(tài)側(cè)壓力，本文采用亞塑性模型描述倉內(nèi)貯料的力學行為，利用ABAQUS軟件建立考慮散-固耦合效應(yīng)的有限元數(shù)值模型，研究淺倉及深倉模型在自重和地震作用下的貯料側(cè)壓力分布，并將分析結(jié)果與中國、歐洲規(guī)范進行對比，以期明確糧食貯料-鋼板筒倉結(jié)構(gòu)體系靜動態(tài)側(cè)壓力的分布規(guī)律。

1 亞塑性本構(gòu)模型

亞塑性本構(gòu)理論是一類新型的散粒型土體本構(gòu)模型理論，其基本思想由德國學者Kolymbas［11］提出。它不是從彈塑性理論發(fā)展而來，沒有屈服面、流動法則、硬化定律以及把應(yīng)變分解為彈性和塑性部分等概念，而是以連續(xù)介質(zhì)力學為基礎(chǔ)，以張量函數(shù)為運算工具，將應(yīng)力速率描述為應(yīng)力、應(yīng)變速率和孔隙比的函數(shù)，主要用于描述砂、碎石等無黏性散粒型土體的非彈性和非線性等主要力學特性。1996 年，Gudehus 和Wu等［12，13］提出了率的相關(guān)方程，隨著研究的深入，von Wolffersdorff［14］對模型提出改進，應(yīng)用于單元試驗的模擬，并將數(shù)值計算與砂土試驗結(jié)果進行對比，結(jié)果表明該模型很好地反映了顆粒材料的非線性特性，但應(yīng)用于預(yù)測小振幅的循環(huán)應(yīng)力或變形時，存在小應(yīng)力循環(huán)下散粒體材料的變形過度累積的棘輪效應(yīng)。為了能夠正確反映土體在小變形范圍內(nèi)的變化，消除棘輪效應(yīng)，Niemunis和Herle［15］基于von Wolffersdorff 提出的本構(gòu)模型，引入了考慮顆粒間應(yīng)變的概念，對模型進行進一步優(yōu)化。

von Wolffersdorff 修正后的亞塑性本構(gòu)模型如下：

式中：φc為材料臨界摩擦角。

特征孔隙比為：

向密性函數(shù)fe、fd的表達式為：

正壓性函數(shù)fb的表達式為：

式中：hs為顆粒硬度；n 為hs與所對應(yīng)的參數(shù)；α、β為確定三軸試驗中內(nèi)摩擦角峰值的相關(guān)參數(shù)；ei、ec和ed分別為顆粒材料的上限孔隙比、臨界孔隙比和下限孔隙比。

由式（1）～（9），von Wolffersdorff 亞塑性本構(gòu)模型共需要8 個參數(shù)：顆粒硬度hs；臨界摩擦角φc；在零壓狀態(tài)下顆粒材料的最大孔隙比ei0、臨界孔隙比ec0和最小孔隙比ed0；參數(shù)n、α、β。

為了消除在小應(yīng)力循環(huán)下散粒體材料變形過度積累的棘輪效應(yīng)，Niemunis 和Herle［15］基于von Wolffersdorff修正的8 參數(shù)本構(gòu)模型，補充5個材料參數(shù)描述材料在荷載逆轉(zhuǎn)后的剛度變化：R、mR、mT、βr、χ。其中，參數(shù)R 定義了顆粒間應(yīng)變的范圍；參數(shù)mR、mT定義了在應(yīng)變空間內(nèi)，應(yīng)變路徑出現(xiàn)90°和180°負載時對應(yīng)于0°的剛度變化；參數(shù)βr、χ 用于對剛度的變化進行平滑處理。

2 數(shù)值模型建立與合理性論證

鋼板筒倉中的貯料是由顆粒組成的帶有空隙的集合體，其動力學行為具有固體和液體雙重特性。本文采用亞塑性模型來表征小麥糧食散體貯料的本構(gòu)模型，綜合考慮了散料的初始孔隙比、顆粒硬度及顆粒間應(yīng)變的影響，具體材料參數(shù)見表1［16］。其中，筒倉倉壁采用4 結(jié)點雙曲線薄殼單元（S4R）建模，貯料采用8 結(jié)點三維實體單元（C3D8）模擬。筒倉倉壁、筒倉底部和筒倉內(nèi)顆粒材料之間的界面由ABAQUS 中提供的“接觸對”算法建模。將筒倉倉壁和倉底表面視為主表面，將與筒倉倉壁和倉底接觸的顆粒材料外表面視為從表面。二者之間的摩擦模型選用庫侖模型，摩擦系數(shù)設(shè)定為0.3［4］。對于接觸約束，考慮了“罰”接觸算法，這類似于在主從表面之間引入剛性彈簧以防止從面穿透主面。為減小誤差，模型中殼單元與實體單元表面的網(wǎng)格劃分為一一對應(yīng)的關(guān)系。

表1 亞塑性模型參數(shù)Tab.1 Parameters of the hypoplastic model

為了驗證以上建模方法的有效性，采用該數(shù)值模擬方法對既有筒倉模型試驗［17，18］進行建模分析。該試驗?zāi)Ｐ蜑槁涞厥綀A形鋼板筒倉，鋼材型號為Q235，筒倉直徑為800mm，高度為900mm，倉壁厚度為1mm，試驗過程中，倉壁始終處于彈性狀態(tài)，故考慮鋼板筒倉為彈性體進行建模計算。

將模擬計算結(jié)果與試驗數(shù)據(jù)進行對比如圖1所示，不同測點位置處的數(shù)值模擬結(jié)果與試驗結(jié)果的最大誤差均在15%以內(nèi)，且二者的分布趨勢基本一致，由此表明亞塑性本構(gòu)模型可以較好地描述小麥貯料的力學行為，從而驗證了利用該建模方法分析貯料筒倉結(jié)構(gòu)體系側(cè)壓力的合理性。

圖1 筒倉倉壁側(cè)壓力試驗值與模擬值對比Fig.1 Comparison of experimental and simulated values of silo wall lateral pressure

3 糧食鋼板筒倉模型建立與結(jié)果分析

3.1 模型建立

中國規(guī)范根據(jù)筒倉高徑比的不同，將筒倉劃分為淺倉和深倉；歐洲規(guī)范劃分為儲存?zhèn)}、粗矮筒倉、中等細長筒倉和細長筒倉，具體劃分標準見表2［4，5］?；谏鲜鼋７椒?，建立淺倉（H／D ＝1.2）和深倉（H／D ＝1.8）兩種有限元模型，淺倉模型高6m，深倉模型高10m，直徑均為5m，厚度0.02m。對模型施加豎直方向上的重力荷載。重力分析完成后，對筒倉底部沿X 方向輸入Elcentro地震波，采用隱式動力分析，地震波總持時為20s，時間間隔為0.02s，參考《建筑抗震設(shè)計規(guī)范》［19］設(shè)計基本地震動峰值加速度，將輸入峰值地面加速度PGA分為0.05g、0.1g、0.2g和0.3g四個等級，其中0.05g對應(yīng)較小地震烈度工況，0.1g和0.2g 對應(yīng)中等地震烈度工況，0.3g對應(yīng)較強地震烈度工況（表3）。圖2 為貯料-鋼板筒倉結(jié)構(gòu)體系模型的網(wǎng)格劃分，其中殼體單元與實體單元的表面網(wǎng)格劃分為一一對應(yīng)關(guān)系。

圖2 有限元模型的網(wǎng)格劃分Fig.2 Meshing of the finite element model

表2 不同規(guī)范的筒倉分類Tab.2 Classification of silos with different specifications

表3 多遇、設(shè)計基本、罕遇地震加速度值Tab.3 Earthquake acceleration values for multiple，basic design and rare intensities

3.2 靜態(tài)側(cè)壓力

對于淺倉而言，中國規(guī)范基于Rankine 理論提出貯料側(cè)壓力的設(shè)計計算公式；而對于深倉側(cè)壓力的計算，兩國規(guī)范都以基于連續(xù)介質(zhì)模型的Janssen理論為基礎(chǔ)，但側(cè)壓力系數(shù)k的取值不同。

中國規(guī)范規(guī)定不同筒倉類型計算深度為s處的貯料靜態(tài)壓力標準值方法如下［4］：

式中：k ＝tan2（45° -φ／2）為側(cè)壓力系數(shù)；γ為貯料的重度；μ為貯料對倉壁的摩擦系數(shù)；ρ ＝dn／4為筒倉凈截面的水力半徑。

歐洲規(guī)范對筒倉類型的劃分不同于中國規(guī)范，本文所建立的有限元數(shù)值模型屬于中等細長筒倉，對于填充結(jié)束時任意深度z處的水平壓力計算方法如下［5］：

式中：K ＝1.1（1-sinφ）為側(cè)壓力系數(shù)；n ＝-（1 +tanφ）（1 -h0／z0）；h0為糧堆錐面重心到錐底的垂直距離；z0為Janssen深度特征值；A為筒倉橫截面面積；U為筒倉橫截面周長；YR為筒倉壓力深度變化函數(shù)。

圖3 為倉內(nèi)貯料的Mises 等效應(yīng)力云圖。可以發(fā)現(xiàn)，沿筒倉深度變化，各筒倉模型貯料的應(yīng)力分布均表現(xiàn)出一定的非均勻性。同一深度面處，靠近結(jié)構(gòu)頂部的貯料應(yīng)力分布較為一致，最大應(yīng)力位于貯料底部位置，且呈現(xiàn)出“中心大，四周小”的趨勢。深倉模型較淺倉模型更為明顯，這表明隨著裝糧高度的增加，受糧食拱效應(yīng)的影響，貯料內(nèi)部應(yīng)力分布的非均勻性更為顯著。

圖3 倉內(nèi)貯料Mises 等效應(yīng)力云圖（單位：Pa）Fig.3 Mises equivalent force cloud of the silo storage material（unit：Pa）

如圖4 所示，將有限元模擬結(jié)果與兩國規(guī)范計算所得不同深度位置處靜態(tài)壓力標準值進行對比分析。可以發(fā)現(xiàn)，淺倉模型與深倉模型側(cè)壓力的數(shù)值模擬結(jié)果均與中國規(guī)范計算所得側(cè)壓力標準值的分布趨勢吻合較好，這也進一步驗證了所建立有限元數(shù)值模型的合理性。

圖4 筒倉靜態(tài)側(cè)壓力分布Fig.4 Static lateral pressure distribution of silos

由圖4 可以看出，對于淺倉模型，根據(jù)中國規(guī)范計算所得設(shè)計側(cè)壓力的分布趨勢更接近于模擬結(jié)果，計算值較有限元結(jié)果最大高出16.5%，具有足夠的安全儲備，表明中國規(guī)范采用Rankine理論對淺倉模型進行側(cè)壓力設(shè)計是合理的，而歐洲規(guī)范則過于保守；對于深倉模型，根據(jù)歐洲規(guī)范計算所得的設(shè)計側(cè)壓力明顯大于中國規(guī)范和數(shù)值模擬的結(jié)果，原因在于中國規(guī)范與歐洲規(guī)范貯料側(cè)壓力計算所采用的基本理論雖然相同，但兩國規(guī)范選取側(cè)壓力系數(shù)不同，使得兩種規(guī)范計算所得的貯料側(cè)壓力存在較大差別。中國規(guī)范計算設(shè)計值與模擬值差距較小，但模型底部的模擬值比規(guī)范值高11.2%，反映出有限元數(shù)值模擬計算時，筒倉底板對倉壁底部約束的增強，使得筒倉倉壁底部的剛度較頂部位置發(fā)生突變而產(chǎn)生應(yīng)力集中現(xiàn)象，從而模型底部位置的側(cè)壓力顯著增大，而規(guī)范中并未考慮筒倉底部約束效應(yīng)的影響。文獻［9］指出伴隨著筒倉高度的進一步增加，此部分差值可能會更大，故對中國規(guī)范側(cè)壓力系數(shù)的取值應(yīng)適當增大。

3.3 動態(tài)側(cè)壓力

除靜態(tài)側(cè)壓力外，糧食散體在地震作用下會對倉壁產(chǎn)生額外的附加法向壓力（相對動態(tài)側(cè)壓力）。歐洲規(guī)范中提出了對地震作用下圓形筒倉倉壁糧食散體產(chǎn)生的動態(tài)側(cè)壓力的估算方法［20］：

式中：Δph，so為基準壓力；θ 為筒倉壁相關(guān)點對應(yīng)的徑線與地震作用水平分量方向之間的角度（0°≤θ ＜360°）。

在離倉底垂直距離為x 處筒倉壁的各點上，基準壓力Δph，so可取值為：

式中：α（z）為離貯料等效面的垂直距離為z處的筒倉反應(yīng)加速度與重力加速度之比；γ 為貯料的重度；＝min（hb；dc／2），hb為筒倉的總高度；dc為與筒倉的內(nèi)徑。

圖5 為淺倉模型和深倉模型在不同PGA作用下沿筒倉高度變化的相對水平動壓力包絡(luò)。由于在對模型進行動力分析前的分析步中考慮了各筒倉模型整體的自重作用，因此，為了得到動態(tài)側(cè)壓力包絡(luò)圖，首先計算各測點的動壓力時程，再用重力分析步末的壓力值減去動壓力時程，最后將各測點動壓力最大值作為動壓力的包絡(luò)值。

圖5 筒倉動態(tài)壓力分布Fig.5 Dynamic pressure distribution of silos

從圖5 可以直觀看出，不同PGA作用下，淺倉模型與深倉模型的相對動態(tài)側(cè)壓力變化趨勢相同，均隨著PGA的增大而逐漸增大，但淺倉模型的相對動側(cè)壓力沿筒倉高度變化較為均勻。在小震與中震區(qū)間內(nèi)（PGA 為0.05g 與0.1g），淺倉和深倉模型最大相對動態(tài)側(cè)壓力均出現(xiàn)在結(jié)構(gòu)的中上部區(qū)域，而頂部較?。浑S著地震動強度的增加（PGA為0.2g與0.3g），深倉模型相對動側(cè)壓力峰值出現(xiàn)的位置較淺倉模型顯著提高，且計算所得動態(tài)側(cè)壓力最大值較淺倉模型分別增大了1.2 倍和2.4 倍。由此可見，PGA 與高徑比的變化對貯料筒倉結(jié)構(gòu)體系動側(cè)壓力分布影響較大；相同PGA作用下，淺倉模型下部的相對動態(tài)側(cè)壓力值要略高于深倉模型。原因在于，隨著高徑比的增大，倉內(nèi)充填貯料隨之增多，此時上部貯料在自重作用下對中下部貯料的壓實作用顯著增大，故在地震作用下，深倉模型中下部的貯料與倉壁更趨于整體運動，貯料對筒倉倉體的碰撞沖擊作用低于淺倉模型；值得注意的是，貯料最大相對動態(tài)側(cè)壓力并不在頂部位置產(chǎn)生，而是在模型的中部及中上部附近位置。這是由于在地震作用下，貯料頂部為自由區(qū)域，可以沿水平及垂直方向運動，而中間位置受到上部及下部貯料的雙重約束，故此時對筒倉側(cè)壁的碰撞及沖擊會較頂部更為劇烈。

圖6 為淺倉及深倉模型最大相對動態(tài)側(cè)壓力平均值（圖5 中壓力分布的平均值）與根據(jù)歐洲規(guī)范計算所得平均值的對比。

圖6 平均相對動態(tài)側(cè)壓力分布與歐洲規(guī)范的對比Fig.6 Comparison of average relative dynamic lateral pressure distribution with eurocode

從圖6 中可以看出，靠近筒倉底部附近，歐洲規(guī)范建議的動壓力值并不保守，這與文獻［21］所得結(jié)果一致。除了底部位置外，淺倉模型的動壓力包絡(luò)均低于歐洲規(guī)范建議的壓力值；而對于深倉模型，靠近結(jié)構(gòu)頂部位置處相對動態(tài)側(cè)壓力分布具有顯著的放大效應(yīng)，最大相對動態(tài)側(cè)壓力超出歐洲規(guī)范計算建議值的47%。這說明，對于淺倉模型，歐洲規(guī)范給出的建議是合理的；在深倉模型中，由于倉內(nèi)貯料對筒倉倉壁的碰撞及沖擊，歐洲規(guī)范假設(shè)α（z）為恒定值（不考慮加速度響應(yīng)沿結(jié)構(gòu)高度的變化）并不合理。

4 結(jié)論

本文采用亞塑性本構(gòu)模型描述倉內(nèi)小麥貯料的力學行為，基于ABAQUS有限元軟件建立了考慮貯料與倉體相互作用的鋼板筒倉有限元數(shù)值模型，研究了淺倉和深倉模型貯料靜動態(tài)側(cè)壓力的分布規(guī)律，并將分析結(jié)果與中國及歐洲規(guī)范進行對比分析，得出以下結(jié)論：

1.鋼板筒倉貯料靜動態(tài)側(cè)壓力有限元模擬結(jié)果與實倉模型振動臺試驗結(jié)果吻合較好，表明利用亞塑性本構(gòu)模型可以較為準確地描述倉內(nèi)小麥貯料的力學行為，同時也驗證了所用建模方法的合理性。

2.自重作用下，貯料內(nèi)部應(yīng)力分布具有明顯的非均勻性，最大應(yīng)力位于筒倉底部，且同一水平面應(yīng)力分布呈現(xiàn)出“中心大，四周小”的趨勢；計算所得淺倉和深倉模型靜側(cè)壓力均與中國規(guī)范吻合較好，Rankine 理論可以較為準確地計算淺倉靜側(cè)壓力；歐洲規(guī)范由于側(cè)壓力系數(shù)偏大導致計算結(jié)果過于保守。

3.同一地震強度下，深倉模型底部附近相對動態(tài)側(cè)壓力較淺倉模型略有降低，隨著地震動強度的增加，淺倉模型的相對動側(cè)壓力沿筒倉高度變化較為均勻，而深倉模型在大震作用下倉頂附近的相對動側(cè)壓力分布具有較大的放大效應(yīng)；歐洲規(guī)范建議的淺倉動壓力分布具有一定的保守性，而深倉動壓力計算中有關(guān)加速度響應(yīng)沿結(jié)構(gòu)高度變化恒定的假設(shè)并不合理。