城市抗震弹塑性分析及其工程应用
陆新征，许镇，程庆乐，熊琛，曾翔，田源，顾栋炼
第一届全国动力多灾害工程结构防护学术研讨会
《动力多灾害作用下工程结构防护研究新进展》
中国建筑工业出版社，2018，198-215
1 研究背景
随着城市化的迅速发展，城市人口、建筑和基础设施的数量和密度迅速提高，因而地震对城市的威胁也在不断增加。城市建筑震害模拟可以揭示地震对城市造成的破坏，服务震前防灾减灾规划和震后快速救援，对减轻城市地震灾害风险具有非常重要的意义。城市抗震弹塑性分析则为城市建筑震害的模拟提供了重要方法。
城市抗震弹塑性分析通过将完整的地震动时程记录输入城市建筑群，逐个建筑进行动力弹塑性时程分析，从而可以充分反映不同建筑的抗震特性差别及不同地震动的时域和频域特征。因此，从理论上说城市抗震弹塑性分析方法与已有方法相比有着明显的优势。但是，为了实现城市抗震弹塑性分析，需要解决海量建筑建模、高性能计算、高真实感可视化、次生灾害预测和考虑“场地—城市效应”等一系列关键科学问题。本研究相应的提出了以下解决办法：
(1) 基于物理驱动模型的建筑群多尺度模型；
(2) 基于CPU/GPU异构并行的高性能计算方法；
(3) 基于3D城市模型和物理引擎的震害高真实感展示方法；
(4) 基于精细化模拟和新一代性能化设计的震损预测和次生灾害模拟方法；
(5) 考虑“场地—城市效应”的区域建筑非线性时程分析方法。
2 城市抗震弹塑性分析方法
本研究提出的城市抗震弹塑性分析方法的整体技术框架如图1所示。
图1 城市建筑抗震弹塑性分析技术框架
2.1 基于物理驱动模型的城市建筑群多尺度模型
城市中建筑数量和种类繁多，本研究将城市中的建筑划分成了常规多层建筑、常规高层建筑和特殊建筑三类，并针对这三类建筑提出了相应的基于弹塑性时程分析的震害预测方法。该方法与传统的基于数据驱动（data driven）的易损性矩阵方法有着本质的不同，是一种基于力学/物理驱动（physics driven）的震害预测模型。
城市区域中的常规多层建筑通常表现出较为明显的剪切变形模式，可以将每栋建筑结构简化成图2a所示的多自由度（multiple-degree-of-freedom, MDOF）剪切层模型。该模型假设结构每一层的质量都集中在楼面上，因此可以将每一层简化成一个质点。不同楼层之间的质点通过剪切弹簧连接在一起。楼层之间剪切弹簧的骨架线采用HAZUS报告中推荐的三线性骨架线（图2c），层间滞回模型采用图2d所示的单参数滞回模型。
与多层建筑不同，高层建筑通常表现出较为明显的弯剪耦合变形形态。因此本研究针对常规高层建筑，采用图2b所示的MDOF弯剪耦合模型。该模型每一层分别用一根弯曲弹簧和剪切弹簧来模拟。每层之间用刚性的链杆连接。弯剪耦合模型的弯曲弹簧和剪切弹簧同样采用图2c和图2d所示的骨架线和滞回模型。
图2 建筑计算模型 (a) 常规多层建筑的MDOF剪切层模型 (b) 常规高层建筑的MDOF弯剪耦合模型; (c) 三线性骨架线; (d) 单参数滞回模型.
城市区域中常规多层建筑和高层建筑的数量巨大，每栋建筑可获取的信息较为有限。因此，本研究提出了图2中常规建筑模型计算参数的自动确定方法。其基本原理是：基于容易获取的宏观GIS数据（主要包括每栋建筑的结构高度、结构类型、建设年代、面积、层数、功能等信息），首先根据统计规律确定结构的振动特性；而后根据规范设计方法确定结构的设计抗震性能；最后根据大量试验和计算结果统计确定建筑实际抗震性能和设计抗震性能的比例关系。这样就可以非常高效地建立数量庞大的城市常规建筑计算模型。需要说明的是，对于可以获取更详细设计信息的建筑，还可以根据设计信息直接确定图2中模型的计算参数，从而获得更好的计算精度。详细模型说明和参数标定方法可以参考论文：
Parameter determination and damage assessment for THA-based regional seismic damage prediction of multi-story buildings, Journal of Earthquake Engineering, 2017.
A nonlinear computational model for regional seismic simulation of tall buildings, Bulletin of Earthquake Engineering, 2016.
MDOF模型能较好地模拟城市区域中常规建筑的地震响应。但是除了量大面广的常规建筑，城市区域中同样存在一些大跨空间结构、超高层结构与异型结构等特殊建筑，这些建筑的动力特性更为复杂，MDOF模型无法满足这些建筑的分析需求。因此，对于这些建筑可以采用基于纤维梁和分层壳模型的精细有限元建模方法加以模拟。
2.2 基于CPU/GPU异构并行的高性能计算
城市抗震弹塑性分析带来巨大的计算工作量，需要寻找高效率低成本的解决方案。近年来，图形处理单元（Graphics Processing Unit, GPU）技术飞速发展，相同价格的GPU相对CPU具有更高的计算性能，因而在不同领域得到了大量成功应用。本研究采用CPU/GPU异构并行计算，加速区域建筑震害模拟过程，提出了图3所示的CPU/GPU异构并行计算流程，该流程主要包括三部分。首先是CPU控制区域计算任务分配，将每栋建筑的计算任务分配给各个GPU核心；之后GPU开始区域海量建筑弹塑性时程分析并行计算；最后CPU将GPU的计算结果输出给后续可视化展示。算例表明，采用GPU/CPU异构并行计算，可以在相近的成本下，将计算效率提高39倍以上，满足了城市区域建筑震害模拟低成本-高效率的要求。相关论文参阅：
A coarse-grained parallel approach for seismic damage simulations of urban areas based on refined models and GPU/CPU cooperative computing. Advances in Engineering Software, 2014.
图3 CPU/GPU异构并行计算流程
2.3 基于城市3D模型和物理引擎的震害高真实感展示
由于城市地震模拟的用户，如政府官员、消防队员等，很多并没有专业的地震工程知识，因此高真实感的震害展示对这些非专业用户理解震害模拟结果有着重要价值。本研究采用城市3D模型对区域建筑的地震模拟结果进行动态展示。
基于城市3D模型的展示方法实现流程如图4-6所示。主要包括建筑对象识别、楼层平面多边形生成与位移插值三个部分。(1) 建筑对象识别将每栋建筑的外表面多边形从城市3D多边形模型中提取出来，并与2D-GIS数据中每栋建筑的描述性信息对应生成3D-GIS数据；(2) 通过建筑对象识别，得到了每栋建筑的外表面多边形，如图4d所示。但是生成结构分析的计算模型往往需要建筑各楼层的平面多边形数据(图5d)。为此提出了楼层平面多边形生成方法，对建筑的外表面多边形进行切片，得到每一楼层的平面多边形(图5)；(3) 建筑时程计算通常生成几个离散高程处（如楼层位置）的结构响应结果，为了保证各层之间描述建筑细节的节点跟随各层发生位移，将采用图6所示的线性插值方法，计算位于两层之间所有节点处的建筑响应结果。其最终效果如图7所示。详细可以参考论文：
Building seismic response and visualization using 3D urban polygonal modeling. Automation in Construction, 2015.
图4 建筑对象识别
图5 楼层平面多边形生成
图6 位移插值
图7 北京CBD3D可视化
采用有限元方法实现倒塌模拟计算成本较高。物理引擎是近些年计算机图形学发展的新技术，专门用于计算场景中刚体碰撞等复杂物理行为。本研究提出可以将物理引擎用于城市建筑群倒塌可视化模拟。在MDOF模型中，采用倒塌层间位移角限值判定结构的倒塌状态（图8a）。MDOF模型可以给出不同楼层倒塌发生时的位移和速度，这些位移和速度作为初始状态传给物理引擎（如PhysX）进行后续倒塌模拟（图8b）。物理引擎模拟楼层刚体在重力作用下的运动，直到楼层刚体间相互碰撞或者接触到地面（图8c）。可视化效果如图9所示。详细可以参考论文：
Seismic damage simulation in urban areas based on a high-fidelity structural model and a physics engine. Natural Hazards, 2014
图8 物理引擎倒塌模拟的过程
图9 城市建筑群倒塌可视化模拟
2.4 基于精细化模拟和新一代性能化设计的震损预测和次生灾害模拟
地震可能对受灾区域带来严重的经济冲击。合理的地震经济损失预测可以为决策者提供重要的参考信息，从而有针对性地制定防震减灾规划、地震保险规划等对策。
本研究基于城市抗震弹塑性分析结果，结合FEMA P-58新一代性能化抗震设计方法，开展建筑地震经济损失预测。其基本原理如图10所示：首先通过城市抗震弹塑性分析得到不同建筑、不同楼层的层间位移角和楼面加速度，然后通过FEMA P-58提供的建筑性能模型和构件易损性数据库，确定不同构件的修复费用、修复时间等震损指标。本研究采用该方法，预测了清华大学校园的地震经济损失（图11）。与既有震损预测方法相比，该方法具有以下优势：(1) 基于精细化的结构模拟结果，可以得到不同楼层不同构件的损伤情况；(2) 既可以考虑楼层位移引起的损失，也可以考虑楼层加速度引起的损失，还可以考虑残余变形引起的损失。详细工作可以参考论文：
Application of the FEMA-P58 methodology for regional earthquake loss prediction, Natural Hazards. 2016.
图10 基于新一代性能化设计方法的建筑地震经济损失预测
图11 清华大学校园地震经济损失预测
精细化震损模拟结果可以进一步用于次生灾害的预测。例如，可以用于坠物分布的模拟，以及次生火灾的模拟等。随着建筑抗震安全性的提高，建筑倒塌造成的伤亡在不断下降，但是建筑外围护结构因地震脱落引起的坠物震害造成人员伤亡以及阻碍人员疏散问题日益突出，而现阶段还缺少合适的坠物次生灾害计算方法。本研究基于城市抗震弹塑性分析，可以得到不同建筑不同楼层的层间位移以及楼面速度，由层间位移可以预测外围护结构是否发生破坏，由楼面速度可以预测坠物的影响范围（图12），进而为避难场所规划和疏散道路设计提供参考（图13）。详细工作可以参考论文：
Simulation of earthquake-induced hazards of falling exterior non-structural components and its application to emergency shelter design, Natural Hazards. 2016.
图12 坠物次生灾害模拟示意
图13 坠物危害影响区域与紧急避难场所选址
次生火灾是地震后导致人员伤亡的另一重要次生灾害。本研究基于城市抗震弹塑性分析得到的精细化震害结果，结合起火概率统计模型和火灾蔓延物理模型，可以预测城市次生火灾风险，并可以通过计算流体力学（CFD）模型得到次生火灾场景下的烟雾蔓延情况（图14）。为城市消防规划和灾后应急预案编制提供参考。详细工作可以参考论文：
Physics-based simulation and high-fidelity visualization of fire following earthquake considering building seismic damage, Journal of Earthquake Engineering, 2017
(a) 次生火灾蔓延模拟
(b) 城市次生火灾情境模拟
图14 城市次生火灾模拟
2.5 考虑“场地—城市效应”的区域建筑非线性时程分析
在建筑密集的城市区域，大量多层、高层建筑在空间上紧密分布，由此而产生的土体-结构相互作用（SSI）和结构-土体-结构相互作用（SSSI）将显著的改变场地的特征以及建筑的地震动输入。从宏观的尺度而言，这种由城市与场地之间相互作用带来的影响被称作“场地-城市相互作用”（即SCI效应）。本研究结合开源谱元分析软件SPEED及2.2节所介绍的非线性MDOF建筑模型，提出了考虑“场地－城市效应”的区域建筑震害模拟方法。
图15为考虑SCI效应的区域建筑震害耦合模拟方法示意图。该方法包含两个主要部分：第一部分是在SPEED软件中模拟地震波在土体中的传播；第二部分是采用非线性MDOF模型进行每栋建筑的时程分析。为了将这两个部分耦合，在每个计算时间步，需要提取建筑的基底反力并将其应用于土体分析；同时需要将土体计算得到的建筑所在位置地面运动加速度作为建筑的基底输入用于建筑分析，循环计算直至完成整个地震作用计算。该方法可以准确的把握SCI效应带来的影响，可以为相关工作提供重要的研究工具。详细工作可以参考：
发了那么多SCI，有兴趣了解一下“SCI效应”么？ | 新论文：考虑“场地-城市效应”的区域建筑震害数值模拟方法
A numerical coupling scheme for nonlinear time-history analysis of buildings on a regional scale considering site-city interaction effects, Earthquake Engineering & Structural Dynamics, 2018
图15 SCI效应的耦合计算方法
3 城市抗震弹塑性分析的应用
3.1 地震后灾损近实时预测
在地震发生后及时预测真实地震损失对震后应急救援具有重要价值。本研究基于实测地面运动记录和城市抗震弹塑性分析，提出了一套近实时的地震后灾损评价方法并开发了相应的系统。该方法：(1) 通过地震台站获取发震地区实测地面运动记录；(2) 建立发震地区典型的区域建筑数据库；(3) 运用城市抗震弹塑性分析方法，将实测地面运动记录输入到目标区域建筑分析模型中，根据区域分析结果评价本次地震对该地区建筑的破坏情况。该方法能较好地解决地震输入的不确定性问题；可以充分考虑地震动的幅值、频谱和持时特征以及不同建筑物的刚度、强度和变形特征；可以评价地震对典型建筑物和目标区域建筑群的破坏能力；在地震发生后短时间内给出地震破坏力评估结果，为科学制定抗震救灾决策和普及公众防灾减灾知识提供了有力手段。
表1地震后灾损近实时预测应用案例
3.2 服务震前防震减灾规划
预测城市未来遭受地震灾害时的损失情况，对制定城市防震减灾决策具有重要价值。在2016年唐山地震40周年时，本研究通过和唐山市以及清华同衡规划院合作，采用城市抗震弹塑性分析方法，对今天的唐山市区23万栋建筑物再度遭遇1976年唐山地震可能导致的破坏进行了分析（图16），为城市防震减灾规划提供量化依据。详细工作可以参考
如果40年前的那次地震再次发生，今天的唐山将会怎样? |唐山大地震40周年纪念
图16 唐山市震害情境模拟结果
3.3 城市中心区高层建筑群多尺度震害模拟
2016年，我们与中国地震局等单位合作完成的1679年三河－平谷地震作用于今天CBD地区的震害场景模拟（图17、18）。得到的基本结论为：根据设计资料（不考虑施工质量等非技术问题）分析表明，在1679年三河－平谷地震作用下，CBD地区的建筑破坏程度基本在轻度到中度之间。其原因主要是CBD地区距离震中已经有差不多50km，地震动强度已经发生了相当大的衰减，地面运动强度大致在8度设防地震到8度罕遇地震之间；与此同时，作为8度设防区，CBD的建筑本身也具有较强的抗震能力。特别是一些重点建筑（比如最高的那个），其抗震性能目标更是高于一般建筑物。因此，造成的建筑物破坏也就基本在轻度到中度之间。详细工作可以参考：
陈肇元院士8年前提出的问题：如果重遇1679年康熙年间8级三河-平谷特大地震会成为什么样子？
图17 三河-平谷地震下北京CBD建筑位移响应可视化
图18 三河-平谷地震下北京CBD建筑损伤演化可视化
3.4 城市建筑震害及次生灾害全过程模拟
城市震害全过程模拟包括断层的破裂、地震波的传播和场地放大、建筑地震响应、经济损失和修复时间预测，以及次生火灾等次生灾害模拟。完成这样的全过程模拟对抗震减灾工作意义重大，但是需要多学科交叉，难度较大。本研究利用城市抗震弹塑性分析和基于精细化模拟和新一代性能化设计的震损预测方法，通过与美国国家科学基金重大项目“多灾害模拟平台SimCenter”合作，将城市抗震弹塑性分析方法应用到美国多灾害模拟平台SimCenter上，完成了旧金山湾区180多万栋建筑的震害模拟（图19），并将相关代码开源到GitHub代码托管平台上，其中，核心的建筑模型参数标定模块和损失确定模块均基于本文所提出的方法。详细工作可以参考：
美国NSF SimCenter+清华城市弹塑性分析=旧金山184万建筑地震模拟
图19 旧金山建筑地震场景
另一方面，选择案例区域中的旧金山市中心建筑群，采用2.4节地震次生火灾模拟方法进行了地震次生火灾模拟。根据当地的年均气象统计数据，设置最低气温为Tlow = 10.5 °C, 最高气温为 Thigh = 17.7 °C，风向为西风，风速为 v = 4.8 m/s。总燃烧建筑占地面积如图20所示，对于该案例而言，火灾的蔓延速度基本恒定，在第2 h蔓延速度稍有增加。9 h后，火灾完全熄灭。最终的火灾损毁情况如图21所示。本案例的地震次生火灾并不太严重，其主要原因是所选案例建筑的间距较大，降低了火灾蔓延风险。在次生火灾场景中添加烟气效果后如图22所示，不仅可以提高场景的真实感，还能更明显地标示出燃烧建筑的位置。
图20 总燃烧面积-时间曲线
图21 火灾蔓延情况 (t = 10 h)
图22 旧金山次生火灾模拟结果高真实感显示
4 结论和展望
科学、准确、直观地模拟城市地震情境并预测地震损失对城市防震减灾工作具有重要价值。随着强震台网的建设、数据传输网络的完善以及计算机分析速度的提高，基于逐个建筑动力弹塑性时程分析的城市抗震弹塑性分析方法在提升震害预测和评估的准确性、高效性和真实感方面具有显著的优势和巨大的发展前景。本文介绍了城市抗震弹塑性分析方法的技术框架和典型应用，初步展示了该方法的可行性和优势。由于城市抗震弹塑性分析还是一个新生事物，现有的技术和方法还有诸多不完善之处，未来有必要在以下方面进一步开展深入研究：
(1) 进一步完善建筑模型，特别是考虑建筑年代影响、不同地域特色的城市建筑模型。
(2) 基于新一代性能化设计方法，进一步完善经济损失预测及人员伤亡、次生灾害预测方法。
(3) 进一步考虑建筑以外其他基础设施（如桥梁、生命线等）的地震破坏。
致谢
感谢国家自然科学基金(No. 51578320，71741024)对本项目研究的支持。感谢清华大学叶列平、任爱珠教授，澳大利亚Griffith University H. Guan教授，日本University of Tokyo M. Hori教授，美国UC Berkeley S. A. Mahin教授、F. McKenna博士，美国Stanford University K. H. Law教授，中国地震局地球物理研究所高孟谭教授，中国地震局工程力学研究所林旭川研究员，加拿大University of British Columbia T. Yang教授，以及课题组研究生韩博、杨哲飚、徐永嘉、孙楚津等对本研究的帮助和支持。
参考文献
[1]Wang Z. A preliminary report on the Great Wenchuan Earthquake. Earthquake Engineering and Engineering Vibration, 2008, 7(2): 225-234.
[2]Lu XZ, Ye LP, Ma YH, Tang DY. Lessons from the collapse of typical RC frames in Xuankou School during the Great Wenchuan Earthquake. Advances in Structural Engineering, 2012, 15(1): 139-153.
[3]Stevenson JR, Kachali H, Whitman Z, Seville E, Vargo J, Wilson T. Preliminary observations of the impacts the 22 February Christchurch Earthquake had on organizations and the economy a report from the field (22 February-22 March 2011). Bulletin of the New Zealand Society for Earthquake Engineering, 2011, 44(2): 65-76.
[4]ATC. Earthquake damage evaluation data for California (ATC-13). Redwood, California: Applied Technology Council (ATC), 1985.
[5]FEMA. Multi-hazard loss estimation methodology-earthquake model. HAZUS-MH 2.1 Technical Manual. Washington, DC: Federal Emergency Management Agency (FEMA), 2012.
[6]Lu XZ, Han B, Hori M, Xiong C, Xu Z, A coarse-grained parallel approach for seismic damage simulations of urban areas based on refined models and GPU/CPU cooperative computing. Advances in Engineering Software, 2014, 70: 90-103.
[7]Lu XZ, Guan H, Earthquake disaster simulation of civil infrastructures: from tall buildings to urban areas, Singapore: Springer, 2017.
[8]Xiong C, Lu XZ, Lin XC, Xu Z, Ye LP, Parameter determination and damage assessment for THA-based regional seismic damage prediction of multi-story buildings, Journal of Earthquake Engineering, 2017, 21(3): 461-485.
[9]Xiong C, Lu XZ, Guan H, Xu Z, A nonlinear computational model for regional seismic simulation of tall buildings, Bulletin of Earthquake Engineering, 2016, 14(4): 1047-1069.
[10]Lu X, Lu XZ, Guan H, Ye LP, Collapse simulation of reinforced concrete high-rise building induced by extreme earthquakes, Earthquake Engineering & Structural Dynamics, 2013, 42(5): 705-723.
[11]Förstner W. 3D-city models automatic and semiautomatic acquisition methods. D. Fritsch, R. Spiller (Eds.), Photogrammetric Week 99, Wichmann Verlag, 1999: 291-303
[12]Michihiko S. Virtual 3D models in urban design // Virtual Geographic Environment, Hong Kong, China, 2008.
[13]Batty M, Chapman D, Evans S, Haklay M, Kueppers S, Shiode N, et al.. Visualizing the city communicating urban design to planners and decision-makers// Centre for Advanced Spatial Analysis, University College London, London, UK, 2001.
[14]Shiode N. 3D urban models recent developments in the digital modelling of urban environments in three-dimensions. GeoJournal, 2002, 52(3): 263-269.
[15]Xiong C, Lu XZ, Hori M, Guan H, Xu Z, Building seismic response and visualization using 3D urban polygonal modeling. Automation in Construction, 2015, 55: 25–34.
[16]Xu Z, Lu XZ, Guan H, Han B, Ren AZ, Seismic damage simulation in urban areas based on a high-fidelity structural model and a physics engine. Natural Hazards, 2014, 71(3): 1679-1693.
[17]Zeng X, Lu XZ, Yang T, Xu Z, Application of the FEMA-P58 methodology for regional earthquake loss prediction, Natural Hazards, 2016, 83(1): 177-192.
[18]Xu Z, Lu XZ, Guan H, Tian Y, Ren AZ, Simulation of earthquake-induced hazards of falling exterior non-structural components and its application to emergency shelter design, Natural Hazards, 2016, 80(2), 935-950.
[19]Lu XZ, Zeng X, Xu Z, Guan H, Physics-based simulation and high-fidelity visualization of fire following earthquake considering building seismic damage, Journal of Earthquake Engineering, 2017. DOI: 10.1080/13632469.2017.1351409.
[22]Federal Emergency Management Agency (FEMA). Quantification of building seismic performance factors (FEMA P695). Washington DC: Federal Emergency Management Agency, 2009.
[23]程庆乐, 许镇, 陆新征, 曾翔, 万汉斌, 张孝奎. 采用城市动力弹塑性分析方法预测唐山市区建筑震害. 自然灾害学报, 2018,(01): 071-80.
[24]肖亮. 水平向基岩强地面运动参数衰减关系研究. 北京: 中国地震局地球物理研究所. 2011.
[25]付长华. 北京盆地结构对长周期地震动加速度反应谱的影响. 北京: 中国地震局地球物理研究所. 2012.
[26]Lu XZ, McKenna F, Zeng X, Cheng QL, Mahin S. Seismic damage simulation of 1.8 million buildings in the San Francisco Bay Area using the nonlinear time-history analysis of buildings. International Conference on Continental Earthquakes. 2018.
[27]Rodgers, A J, Pitarka, A, Petersson, N A, Sjögreen, B, Mccallen, D. B. Broadband (0-4Hz) ground motions for a magnitude 7.0 Hayward fault earthquake with three-dimensional structure and topography. Geophysical Research Letters, 2018, 45.
作者简介
陆新征：清华大学土木工程系教授，博士生导师
许 镇：北京科技大学土木与资源工程学院副教授
程庆乐：清华大学土木工程系博士研究生
熊 琛：深圳大学土木工程学院助理教授
曾 翔：清华大学土木工程系博士
田 源：清华大学土木工程系博士研究生
顾栋炼：清华大学土木工程系博士研究生
软件操作视频
为了配合十二五科技支撑计划课题的检查汇报，课题组研究生去年制作了一个“城市抗震弹塑性分析”的软件操作演示视频，供大家了解“城市抗震弹塑性分析”的具体实现过程。需要说明的是，该课题当初立项的名称为《基于并行计算与地理信息技术的城市建筑震害仿真技术》，所以视频也采用该名称。
----------------End----------------------