单檐歇山式古建模型动力特性与多遇地震响应
周 乾1，纪金豹2
（1.故宫博物院 北京 100009；2.北京工业大学 工程抗震与结构诊治北京市重点实验室，北京 100124）
摘 要：为了更好的保护古建筑，以某单檐歇山古建筑模型为对象，基于已有的试验成果，采用有限元分析方法，进一步研究了其动力特性及8度多遇地震响应。采用弹簧单元模拟榫卯节点及斗拱构造，并考虑柱础为平摆浮搁，建立了古建筑模型的有限元模型。通过模态分析，获得了模型的基频及主振型；通过对模型输入PGA=0.1g（x向）的水平双向地震波（x:y=0.85:1），获得了典型节点的位移、加速度响应曲线、典型单元的内力响应曲线，评价了模型的抗震性能。结果表明：该单檐歇山模型的基频为1.34Hz，主振型以平动为主；8度多遇地震作用下，模型能保持稳定振动状态，结构的内力和变形均在容许范围内。此外，有限元分析结果与试验结果符合较好。
关键词：古建筑；木结构；动力特性；地震响应；单檐歇山屋顶
中图分类号：TU366.2 文献标识码：A
我国的古建筑以木结构为主，具有重要的历史文化价值，保护意义重大。近年来地震频发，对古建筑的安全构成一定威胁[1-3]。相应的，部分学者开展了古建筑性能评估，提出了科学保护建议[4-7]。然而对于单檐歇山屋顶类古建筑而言，其相关成果较少。笔者曾以故宫某单檐歇山屋顶类古建筑为原型，以红松为主要材料，制作了1：2的缩尺比例模型(图1)，开展了以振动台试验为主的研究，获得了该类古建筑在不同烈度地震作用下的试验结果[8-9]。本研究采用有限元分析方法，对该古建筑模型开展模态分析和时程分析，进一步探讨在8度多遇地震作用下，模型的动力特性及时程响应特点，以丰富试验成果，更好地保护我国古建筑。
图1 模型尺寸图(单位:mm)
1 有限元模型
本研究中的有限元模型的建立，主要考虑以下几个方面的假定：
(1) 采用有限元分析程序ANSYS建立单檐歇山古建模型的有限元模型。
(2) 柱础：约束条件考虑为半刚接。其主要原因在于：柱底平摆浮搁在柱顶石上，中、低烈度地震作用下，地震力超过柱顶与柱顶石之间的静摩擦力后，柱底可绕柱顶石产生相对滑移和转动。采用COMBIN40单元模拟柱顶—柱顶石半刚接特性，主要输入参数包括柱底起滑力F、柱子初始侧向刚度k。考虑柱底与柱顶石之间的静摩擦系数为0.5[10]，计算单柱上部的荷载后，可获得F=2568N。另参照文献[11]的研究成果，取值k=4.5×105N/m。
(3)榫卯节点及斗拱：采用非线性弹簧单元COMBIN39模拟榫卯节点及斗拱。本分析主要考虑中低烈度地震作用，因而主要考虑榫头与卯口在弹性运动阶段的相对转动，以及斗拱构件在弹性运动状态下的相对滑移运动；相应的，COMBIN39单元仅考虑弹性刚度取值。参照文献[12]的试验结果，本研究取值榫卯节点的转动刚度Krotx= Kroty= Krotz=5.755kN∙m/rad；其中，Krotx、 Kroty、Krotz3分别表示节点绕x、y、z轴的转动刚度。参考文献[13]-[14]的研究结果，本研究取值斗拱的平移刚度Kx=Ky=0.3×106N/m，Kz=3.21×106N/m；其中，Kx、Ky分别表示斗拱的水平双向刚度，Kz表示斗拱的竖向刚度。
(4)屋顶：采用SHELL181单元模拟屋面板及上部的灰背、瓦件。根据屋顶分层做法及材料组成，求得SHELL181单元主要输入参数取值为：密度1.8×103kg/m3，等效厚度0.134m。另对于正脊、垂脊和戗脊质量，将其简化为均匀分布的质点单元，附在各条脊的位置，并采用MASS21单元模拟。
(5)梁、柱：采用线性梁单元BEAM189单元模拟，主要输入参数包括弹性模量E=9x109N/m2，泊松比γ=0.3，密度ρ=500kg/m3。
基于上述假定，建立本分析的有限元模型，见图1所示。需要说明的是，本模型为ANSYS程序中，打开形状开关后(/esha,1)的效果，显示为实体方式。后续分析中，为减小模型占用计算内存，采用关闭形状开关(/esha,0)的显示方式，模型显示为线面方式。
图2 有限元模型
2 动力特性
对单檐歇山有限元模型进行模态分析，求得模型的前10阶自振频率，见表1；获得模型的前1、2阶振型见图3。易知：(1) 模型主振型集中在第1、2阶。其中第1阶振型表现为模型纵向（长度方向）振动为主，第2阶振型表现为模型沿横向（宽度）振动振动，且模型在上述两个方向的振动关联很小。(2) 由于单檐歇山屋顶的构造特点，其质量在纵、横两个方向分布不同，导致振动时有轻微的平面摇晃。(3) 由于柱底平摆浮搁在柱顶石上，因而模型在振动时，柱底与柱顶石之间有较为明显的相对滑移运动，但这种运动形式仅见于中柱；但对于边柱（模型两侧的立柱）而言，由于有嵌固墙体的约束作用，因而振动形式表现为柱底绕柱顶石的转动。(4) 墙体底部固定于地面，其振动形式表现为绕地面的纵向、横向摇摆。由于墙体的抗震性能相对木构架要差[15]，因而在地震作用下，通过这种运动方式，先于木构架破坏。
另本研究获得的模型基频为1.34Hz，与现场实测结果偏差小于10%[8]，可近似反映有限元模型是有效的。
表1 模型自振频率
(a) 第1振型正立面
(b) 第2振型侧立面
图3 结构主振型图
3 地震响应
3.1 结构阻尼确定
在进行时程分析时，采用瑞利阻尼作为本研究有限元模型的结构阻尼，计算公式见式(1)：
式(1)中，α、β为比例常数，用于程序输入的阻尼参数，且可按(2)-(3)式确定；[C]、[M]、[K]分别为模型的阻尼矩阵、质量矩阵和刚度矩阵。
式(2)-(3)中，ξi、ξj分别为第i、j振型阻尼比，ωi、ωj分别为第i、j振型圆频率。对式(2)-(3)进行变换，可获得α、β的求解式，见式(4)-(5)：
又
一般取i=1，j=2，ξi=ξj=0.05，将模型前2阶频率f1=1.34Hz及f2=1.64Hz代入(6)式，并联立(4)-(5)式，解得α=0.152，β=0.016。
3.2 地震波的选取
本文主要研究单檐歇山屋顶类古建筑在中、低烈度地震作用下的响应特征。采用与文献[8]相同的地震波作用于模型：适合于故宫场地类型(ІІ类)的1940年El-Centro波，x、y双向，其中在x向的加速度峰值为PGA=0.1g，水平双向加速度峰值比例为x:y=0.85:1，时间间隔为0.01s，持续作用时间为15s。北京为8度抗震设防区，上述地震波加速度峰值略大于2016年版的《建筑结构抗震设计规范》(GB50011-2010)表5.1.2-2中，关于8度多遇地震加速度时程峰值的规定，因而计算结果偏于保守。本分析用的地震波波形见图4。
(a)x向
(b) y向
图4 El-centro波形图
3.3 地震响应曲线
为全面了解地震作用下结构的内力及位移响应状况，结合振动台试验结果，选取代表性的节点或单元进行分析：(1) 明间金柱柱底 (节点编号：125)、明间屋脊正中 (节点编号730) ，以研究结构的位移及加速度响应特征；(2) 明间桃尖梁跨中(单元编号157)，以研究结构的弯矩响应特征；(3) 明间桃尖梁端部(单元编号153)，以研究结构剪力响应特征。上述各节点及单元位置详见图2。
125、730号节点的位移响应曲线见图5，主要特点为：
(1)由于地震波强度较低，且由于模型两侧嵌固墙体的约束作用，浮放柱底与柱顶石之间的摩擦滑移运动不明显，表现以转动为主的运动形式，与文献[8]中的试验现象较符合。125号节点在x向的位移峰值为0.92mm(t=3.06s)，在y向的位移峰值为0.91mm(t=1.76s)，与文献[8]中的试验结果相近。相应的，模型顶部730号节点在x向的位移峰值为7.08mm(t=3.06s)，在y向的位移峰值为10.55mm(t=1.76s)，可反映屋顶位移相对于柱底有较明显的放大。
(2)模型的位移响应曲线表现为以平衡位置为中心，近似往复的均匀的振动，反映了结构在地震作用下处于稳定的运动状态，且与文献[8]中的试验现象较符合。
(a) x向
(b) y向
图5 125、730号节点位移响应曲线
125、730号节点加速度响应曲线见图6，主要特点为：
(1)从加速度响应峰值来看，在x向，125号、730号节点的峰值分别为0.091g（t=1.24s）、0.046g(t=1.07s)；在y向，125号、730号节点的峰值分别为0.101g(t=0.59s)、0.062g(t=1.67s)。
(2)从加速度峰值的变化特点来看，对于同一节点而言，y向的加速度峰值大于x向，其主要原因在于y向的输入地震波峰值要更大。另节点125的加速度峰值要略小于输入地震波的加速度峰值，可反映中低烈度地震波作用下，平摆浮放的柱底可产生一定的摩擦减震作用。节点730的加速度响应峰值要小于节点125，可反映榫卯节点、斗拱等构造可产生耗能作用。
(a) x向
(b) y向
图6 125、730号节点加速度响应曲线
此外，由于节点730位于屋脊部位，其加速度峰值较小，可反映屋顶瓦件受到的地震力较小。而结合试验现象来看，在PGA=0.1g烈度的地震波作用下，屋顶瓦件没有出现松动或掉落现象(图7)。
图7 振动台试验截屏
为探讨不同时刻地震作用下，模型的内力分布特点，以t=1.67s为例进行分析。基于ANSYS程序，获得了该时刻的模型主应力分布图，见图8。易知：模型的大主拉应力峰值（0.325MPa）、主压应力峰值（-0.325MPa）均位于边柱与额枋相交的榫卯节点位置，其周边沿柱身方向应力值较大。此处的榫头、卯口截面均有尺寸削弱，节点在地震力作用下产生的相对运动明显，且该位置为屋顶竖向荷载集中往下传递的交汇处，再加上墙体的挤压，因而相对于其他位置，更容易出现受拉或受压破坏。尽管在本烈度作用下，上述应力峰值低于《木结构设计标准》（GB50005-2017）规定的红松应力容许值（拉应力8MPa，压应力10MPa）[16]，但应力峰值的位置是模型易出现破坏的位置，在高烈度地震作用下容易产生破坏。
(a) 主拉应力
(b) 主压应力
图8 模型主应力图（t=1.67s，单位：MPa）
文献[9]对本模型开展了极端烈度地震作用下的振动台试验，并进行了震后调查。调查照片显示了模型的典型震害部位及震害特征（图9）：模型两侧的立柱与额枋相交处，出现了尺寸较长的粗裂纹（见虚直线表示），且沿柱身方向扩展；另榫头与卯口之间亦出现错动间隙（见虚线圆圈内）。以上可反映上述位置受到了明显的破坏，且与图1所示模型的主应力图特点有着相似之处。
(a) 左侧
(b) 右侧
图9 文献[9]中的模型震害照片
157号单元在x、z向（考虑结构自重）的弯矩响应曲线见图10；弯矩峰值分别为91.8N·m（x向）、25.9N·m（z向），对应截面(240×325mm)最大弯应力值分别为0.03MPa（x向）、0.006MPa（z向）。上述峰值远低于《木结构设计标准》（GB50005-2017）规定的红松弯曲应力容许值（13MPa）[16]。因此，该桃尖梁有充足的的抗弯强度储备。另一方面，研究表明[1,9]：我国古建筑梁架截面尺寸普遍充足，因而很少产生弯曲破坏。
图10 157号单元弯矩响应曲线
153号单元的剪力响应曲线见图11。采取平方和平方根法对上述曲线进行合成，可得单元受剪力峰值为F=324N，相应截面（桃尖梁截面尺寸为240×325mm，考虑上部做成桃尖形、端部刻出卯口，实际截面尺寸取120×160mm）最大剪应力fs=0.017MPa，低于《木结构设计标准》（GB50005-2017）规定的红松剪切应力容许值（1.4MPa）[16]。即该位置不会产生剪切破坏。
图11 153号单元剪力响应曲线
4 结论
本研究基于有限元分析方法，对已开展振动台试验的单檐歇山试古建筑模型，进行了动力特性及8度多遇地震响应分析，得出如下结论：
(1) 模型的计算基频为1.34Hz，主振型以平动为主，屋顶有轻微水平面摆动。
(2) 8度多遇地震作用下，模型柱底与柱顶石之间的相对摩擦滑移不明显，上部结构保持较为稳定振动状态。
(3) 8度多遇地震作用下，模型典型单元的内力峰值均在容许值范围内。
(4) 本研究中的数值模拟结果与振动台试验结果符合较好。
参考文献：
[1]周乾，闫维明，杨小森，等.汶川地震导致的古建筑震害[J].文物保护与考古科学,2010，22（1）：37-45.
[2]宋耀春, 李占国.4.14玉树地震古建筑损害调查及应急性保护对策[J].青海社会科学，2010，（04）：94-98.
[3]安仁兵，游文龙，潘毅，等.泸县6.0级地震文化遗产震害调查与分析[J].土木工程学报，2022，55（12）：13-24.
[4]周乾.基于构造特征的故宫太和殿抗震性能分析[J].建筑技术， 2022，53（3）：281-284.
[5]薛建阳，吴晨伟，张雨森.殿堂式木结构古建筑屋盖梁架体系动力特性及地震响应分析[J].建筑结构学报，2021，42（10）：87-95.
[6]周乾.紫禁城古建筑的传统防震方法[J].工业建筑， 2019，49（6）：210-217.
[7]袁建力.墙体参与工作的木构架古建筑抗震分析方法[J].建筑结构学报， 2018，39（9）：45-52.
[8]周乾，闫维明，纪金豹，等.单檐歇山式木构古建抗震性能振动台试验[J].文物保护与考古科学，2018，30（2）：37-53.
[9]周乾，闫维明，杨娜，等.单檐歇山式木构古建抵抗极端烈度地震试验调查研究[J].文物保护与考古科学,2020,32(2):29-44.
[10]姚侃,赵鸿铁.木构古建筑柱与柱础的摩擦滑移隔震机理研究[J].工程力学，2006，23（8）：127-131.
[11]贺俊筱，王娟，杨庆山.考虑高径比影响的木结构柱抗侧能力试验研究[J].土木工程学报，2018，51（3）：27-35.
[12]周乾，闫维明，周锡元，等.古建筑榫卯节点抗震性能试验[J].振动、测试与诊断，2011，31（6）：679-684.
[13]周乾，杨娜，闫维明,等.故宫太和殿一层斗拱水平抗震性能试验[J].土木工程学报，2016，49（10）：18-31.
[14]周乾，闫维明，慕辰曦，等.故宫太和殿一层斗拱竖向加载试验[J].西南交通大学学报，2015，50（5）：879-885.
[15]周乾，闫维明，纪金豹.含嵌固墙体古建筑木结构震害数值模拟研究[J].建筑结构，2010，40（1）：100-103.
[16]中华人民共和国住房和城乡建设部.木结构设计标准(GB50005-2017)[S].北京：中国建筑工业出版社，2018：21-22.
（《水利与建筑工程学报》2023年第5期刊载）