地震作用下柱承式筒仓动态侧压力计算方法研究

陈佳，丁永刚，2，许启铿，2，刘强，索焕

（1.河南工业大学土木工程学院，河南 郑州 450001；2.河南省粮油仓储建筑与安全重点实验室，河南 郑州 450001）

引言

钢筋混凝土柱承式筒仓是应用于粮食、煤炭、建材等仓储物流行业的通用性构筑物。在地震作用下，筒仓除了受到结构本身的惯性力外，仓壁上还会产生贮料的动态侧压力，引起仓壁的侧向变形和较大的动态超压问题［1］。在地震频发区，已经成为筒仓结构破坏的主要原因之一［2］。国内外相关学者对筒仓结构的动态侧压力问题进行了相应研究，Chowdhury等［3－4］将贮料看作一种集中质量施加到仓壁上，推导了矩形储煤仓结构仓壁动态侧压力计算方法，部分学者通过振动台试验研究证明该方法和欧洲规范规定的计算方法过于保守［5－7］。Silvestri等［8－10］基于贮料散体和仓壁相互作用及颗粒剪切效应提出了新的动态侧压力计算方法。但通过钢筒仓振动台试验［11］研究验证，发现该方法仅适用于浅圆仓（高径比大于0.4且小于1.0），尚不适用于筒仓［12］（高径比大于或等于1.5）。施卫星等［12］采用非线性地震反应分析方法，在筒承式单仓振动台试验的基础上建立了煤仓动态侧压力计算公式。在此基础上，周长东等［13］考虑了贮料高度、地震动强度、贮料内摩擦角和贮料与仓壁之间的摩擦系数对侧压力的影响，通过数值分析确定了贮料侧压力修正系数，进一步推导了筒仓和贮料相关作用的贮料侧压力修正公式。刘海林［14］尝试从离散元角度分析了不同高径比情况下的筒仓卸料压力增加机理，提出了适用于圆形浅仓且考虑仓壁和贮料弹性的侧压力计算方法。

由于缺乏筒仓贮料动态侧压力计算方法，现行规范并未涉及筒仓结构地震作用引起的动态侧压力，而是仅考虑了筒仓卸料为主要因素的超压系数进行结构设计［15］。显然，随着仓储结构抗震设防和防灾减灾能力需求不断提高，急需探索研究地震作用下筒仓结构贮料动态侧压力的计算理论和方法，以进一步提高筒仓结构安全和抗震韧性。

丁永刚等［16－19］进行了仓储结构地震响应研究，通过振动台试验获得了柱承式筒仓、排仓以及群仓结构在地震作用下贮料动态侧压力的分布规律和超压系数。为此，文中在现有研究基础上，以柱承式筒仓结构为研究对象，通过理论分析探究适用于筒仓地震作用下的贮料动态侧压力计算方法。在此基础上，基于柱承式筒仓结构的数值算例与振动台试验数据，验证了该方法的有效性，并分析其在不同加速度峰值下的仓壁动态侧压力分布规律，为完善筒仓结构抗震设计理论和方法提供参考。

1 筒仓仓壁分析模型

1.1 基本假设

基于Airy静压力理论［20］，确定带有失效楔块的柱承式筒仓计算简图，如图1所示。针对高径比较小的浅圆仓图1（a），仓内贮料形成了失效楔块ABC后，仓壁上产生静压力，假设失效楔块形成的角度为θ。根据失效楔块ABC的静力平衡条件可知，作用在仓壁上的水平力P可以表示为：

根据Airy理论［20］，若P取最大值，则应满足，进而可得：

式中：P为仓壁水平侧压力；W为失效楔块的重力；θ为失效楔块形成的角度；μ为仓内贮料颗粒间摩擦系数；μ′为仓内贮料与仓壁间摩擦系数。

为确定基于仓壁壁板振动分析的计算模型，结合上述分析，有下列基本假设成立：

（1）由式（2）可知，筒仓内已经形成楔块的贮料处于极限失效状态，本身没有刚度，在地震作用下，可以认为失效楔块ABC被悬挂在仓壁上，只提供惯性作用。

（2）与失效楔块相比，柱承式筒仓仓壁较薄，假设仓壁只对筒仓结构提供刚度贡献。

（3）考虑到筒仓顶部通常设有仓上设备，将仓壁视为四边为固定约束的壁板。

对于图1（b）所示的筒仓结构，仓壁上产生的静压力是由于仓内贮料形成了失效楔块ABCD后引起的。对于地震作用下筒仓动态侧压力问题可以等效为支撑失效楔块ABCD的四边固定壁板的振动问题。同样，可以得到：

图1 带有失效楔块的柱承式筒仓简图Fig.1 Structure diagram of a column-supported silo with failure wedge mass content

式中：H为仓壁高度；d为筒仓直径；θ为失效楔块形成的角度；μ为仓内贮料颗粒间的摩擦系数；μ′为仓内贮料与仓壁间摩擦系数。

筒仓中失效楔块ABCD由三角形和矩形组成，其惯性质量对筒仓结构在地震作用下的动态侧压力响应有贡献。为方便公式推导，将失效楔块ABCD转化为具有相同高度（H）但直径不同的面积相等的等效三角形。根据面积关系S′Δ=SABCD和tanθ′=H/d′Δ，可得：

由式（5）可以看出，对于筒仓结构仓壁来说，失效楔块形成的夹角θ′受高径比影响。

1.2 计算模型

为方便计算，针对柱承式筒仓仓壁，取简化后的壁板模型进行分析。该简化模型可以实现将仓壁动态侧压力问题简化为支撑三角形失效楔块的四边固定壁板的振动问题。取筒仓仓壁部分弧形壁面进行分析，展开得到高宽为H×b的壁板。地震作用下仓壁壁板振动分析模型如图2所示。

图2 地震作用下支撑失效楔块的深仓壁壁板（四边固定）Fig.2 Slender silo wall with failure wedge mass content under earthquakes（four edges fixed）

由振动的模态叠加原理［21］可知，地震作用下，仓壁的变形按照振型函数展开为：

式中：qi(t)为广义坐标；ϕi(z,y)为仓壁壁板在地震作用下的振型函数。

筒仓结构的动能可以表示为：

将式（6）代入式（7），得到：

式（8）第3项为质量系数表达形式，即：

对于仓壁壁板来说，其弯曲变形能为：

根据拉格朗日方程，可以证明式（10）关于qj(t)对于各向同性的仓壁壁板的刚度系数可以表示为：

式中：D为薄板的弯曲刚度，D=，其中，E为材料弹性模量；δ为仓壁厚度；υ为材料泊松比。ϕi,ϕj为仓壁壁板在在地震作用下z和y方向的振型函数。

选取仓壁壁板y和z方向上的微条单元进行分析，如图3所示。根据假设条件（3），微条被简化为两端固定的梁单元。通过求解两端固定梁弯曲振动微分方程，可得到振型函数，即：

图3 有限元模型及网格划分（单位：mm）Fig.3 FEM model and grid division（Unit：mm）

2 仓壁动态侧压力计算方法

2.1 等效质量和刚度系数

深仓仓壁的质量系数和刚度系数是仓壁动态侧压力计算公式推导的基础物理量。为了方便后期推导，首先确定支撑失效楔块的仓壁（简化计算模型）的质量和刚度系数。通过对微条单元取基振型进行分析，ϕi(z,y)=ϕj(z,y)，将其代入式（9），得到：

将式（12）代入式（13），并令z=Hξ,y=bη即可得到：

同样，将式（12）代入式（11）中可得到，

为了方便计算，将振型函数F(z)和F(y)转化为自变量在（0～1）范围内的函数f(ξ)和f(η)。令z=Hξ,y=bη，则可得到：

将式（16）代入式（15），得刚度系数为

式中：C1∼C5可采用辛普森积分法求解；γ为储料重度；H为振动模型壁板高度；b为振动模型壁板宽度；g为重力加速度；θ为失效楔块形成的角度；D为仓壁壁板的抗弯刚度。

2.2 仓壁动态侧压力

基于弹性力学理论，将地震作用下筒仓仓壁的动态侧压力问题简化为支撑失效质量块的贮料壁板的振动问题。因此，根据薄板弯曲的基本微分方程求解筒仓仓壁在地震作用下的动态侧压力，即：

式中：D为薄板的弯曲刚度；w(z,y,t)为地震作用下仓壁的位移函数

式中：Sa为柱顶加速度反应谱；ω为固有频率；ki为振型参与系数，且由式（14）计算可取，ki=，可采用辛普森积分法求解。

式中：a为仓壁壁板的高径比（H/b）；C1∼C5可采用辛普森积分法求解。

令f(a)=，则f(a)作为已知量代入式（22），结合式（16）和式（21），得到：

将式（23）代入式（20），最终得到仓壁上产生的动态侧压力：

由式（24）可以看出，动态侧压力P受筒仓结构高径比、仓内贮料颗粒间、贮料与仓壁间摩擦系数、柱顶加速度反应谱、仓壁振型函数的影响。其中，柱顶加速度反应谱可以通过振动台试验实测获得。基于四边固定壁板振动模型的动态侧压力随仓壁振型函数选取的不同而变化，是基本形状函数及其及其导数、柱顶加速度反应谱和高径比的函数。该公式反映了高径比大于或等于1.5条件下柱承式筒仓仓壁振动对筒仓结构的动态侧压力贡献大小。该公式为计算柱承式筒仓贮料动态侧压力提供了理论参考。

3 数值算例与对比验证

3.1 数值算例

选取几何尺寸如图3（b）所示的深仓条件下的柱承式筒仓进行有限元ABAQUS数值分析。筒仓高径比为2.56，筒仓壁厚取10 mm。有限元模型材料与丁永刚等［17］开展的振动台试验模型保持一致。结构和贮料均采用实体单元。贮料与仓壁之间考虑接触关系，并设置切向接触为摩擦，法向接触为硬接触，贮料与仓壁之间摩擦系数取0.45。仓壁约束条件结合筒仓实际支撑条件和我国规范确定为：支承柱顶部与仓壁底部接触位置为绑定关系；支承柱底部为固定约束；考虑到仓上设备的存在，设置仓壁顶部为绑定。贮料采用D－P本构模型。三维有限元模型网格划分如图3（a）所示。文中根据丁永刚等［17］开展的振动台试验所确定的贮料材料参数取值如表1所示。

表1 贮料材料参数Table 1 Parameters of storage materials

3.2 对比验证

由式（24）可以看出，地震作用下的动态侧压力受筒仓结构高径比、仓内贮料颗粒间、贮料与仓壁间摩擦系数、柱顶加速度反应谱、仓壁振型函数的影响。对于以上4个方面的影响因素均应采用试验验证其影响规律。限于篇幅，对于某一个确定的筒仓（深仓）结构而言，地震作用对其动态侧压力分布规律起主导作用。因此，文中进行了不同地震加速度峰值下筒仓侧压力理论计算值、有限元值和试验结果对比分析，用以验证文中提出的基于筒仓壁板分析的动态侧压力计算方法的正确性。

依据文中建立的有限元数值分析模型，在筒仓结构的x方向分别输入不同加速度峰值的El Centro波、唐山波和人工波，加速度峰值分别为0.062 5、0.125、0.25 g。计算时取10 s地震波记录。考虑筒仓结构处于线弹性变形范围，得到柱承式筒仓结构和仓壁侧压力计算结果。限于篇幅，图4和图5所示为El Centro波作用下的仓壁侧压力和筒仓结构等效应力云图。

图4 筒仓壁板侧压力Fig.4 Lateral pressure of slender silo wallboard

图5 筒仓结构等效应力Fig.5 Equivalent stress of slender silo

由图4和图5可以直观地观察到，考虑筒仓结构处于线弹性范围，柱承式筒仓仓壁在地震加速度峰值分别为0.062 5、0.125、0.25 g时，其仓壁壁板的应力峰值均出现在壁板底部。这是由于有限元数值分析提取的筒仓壁板侧压力结果是基于静态侧压力的叠加值。图4筒仓壁板侧压力云图趋势表现为地震作用下仓壁受到的动态侧压力与静态侧压力的叠加。鉴于此，建议在柱承式筒仓结构的柱顶和环梁位置采取加强措施，以保证结构的整体安全性。

在团队前期研究成果中，丁永刚等［17］开展了与文中数值模型和相关材料参数一致的单向地震作用下柱承式筒仓振动台试验，获得的动态侧压力分布规律沿筒仓仓壁深度方向呈现“中间大，两端小”。试验所测侧压力为动态侧压力净值，不包括静态侧压力。为了更清晰的呈现对比分析结果，采用文中提出的仓壁动态侧压力计算方法得到的柱承式筒仓仓壁壁板侧压力理论值。结合有限元计算结果，提取了壁板动态侧压力峰值。最终将文中得到的理论值、有限元计算结果和丁永刚等［17］开展的柱承式筒仓模型振动台试验结果进行对比分析，分析结果如表2所示。

由表2数据可以得到，在El Centro波的地震作用下，筒仓结构的壁板应力随着地震加速度峰值的增加而不断增大，且理论值、试验值和有限元值吻合良好，相对误差均小于4.0%，这表明文中提出的柱承式筒仓动态侧压力计算方法合理可行。该方法可以为筒仓结构在地震作用下的动态侧压力计算提供理论参考。

表2 不同地震加速度峰值下筒仓壁板侧压力峰值Table 2 Peak lateral pressure of silo wallboard under different seismic acceleration peaks Pa

为了更直观地观察动态侧压力沿筒仓深度方向的变化规律，提取了柱承式筒仓结构在El Centro波、唐山波和人工波作用下加速度峰值分别为0.062 5、0.125、0.25 g时仓壁不同位置处（0.3、0.6、0.9、1.2 m）侧压力计算结果，绘制侧压力随深度变化图，如图6所示。

图6 3种地震波作用下不同地震加速度峰值的仓壁侧压力沿深度分布Fig.6 Lateral pressure along the depth of slender silo wall under different seismic acceleration peaks with the action of 3 different seismic records

由图6可以得到，考虑筒仓结构处于线弹性范围，柱承式筒仓仓壁在地震加速度峰值分别为0.062 5、0.125、0.25 g时，在一定深度范围内仓壁壁板侧压力变化规律是随深度增加而增加。但超过一定范围后，随深度的增加而减小，但呈非线性减小趋势，表现为动态侧压力响应规律沿筒仓仓壁深度方向呈现“钟形分布”（中间大，两端小）。然而，对于不同地震波作用下柱承式筒仓仓壁壁板侧压力变化规律整体趋势一致，但人工波作用下不同地震加速度峰值的仓壁壁板侧压力响应较为敏感。这是由于贮料特性、贮料与仓壁摩擦特性的影响，使筒仓结构受到的地震反应随不同地震波类型和加速度峰值的变化有所差异。

4 结论

通过将地震作用下筒仓动态侧压力问题简化为支撑失效质量块的贮料壁板振动问题，提出了基于四边固定弹性壁板振动模型的柱承式筒仓动态侧压力计算方法，并建立有限元数值分析模型，进行不同加速度峰值情况下的结果验证与对比分析。主要结论如下：

（1）提出的柱承式筒仓动态侧压力计算方法能够合理地反映地震作用下仓壁侧压力分布规律，计算精度较高，为筒仓结构地震作用下的动态侧压力计算提供了理论依据。

（2）柱承式筒仓结构的动态侧压力随地震加速度峰值的增加而增大。且沿筒仓深度方向呈现“钟形分布”趋势（中间大，两端小），其峰值向仓壁中下部位置集中。动态侧压力是筒仓结构设计不可忽视的一个重要作用力。

（3）针对设防烈度较高地区，建议在柱承式筒仓结构柱顶和环梁部位采取加强措施，以保证筒仓结构的整体安全性。