A highperformance quadrilateral flat shell element for seismic collapse simulation of tall buildings and its implementation in OpenSees Xinzheng Lu^{1,*}, Yuan Tian^{2}, Song Cen^{3}, Hong Guan^{4}, Linlin Xie^{5}, and Lisha Wang^{3} Journal
of Earthquake Engineering, ^{1} Key Laboratory of Civil Engineering Safety and Durability of China Education Ministry, Department of Civil Engineering, Tsinghua University, Beijing, P.R. China, 100084 ^{2} Beijing Engineering Research Center of Steel and Concrete Composite Structures, Tsinghua University, Beijing 100084, China ^{3} Department of Engineering Mechanics, School of Aerospace Engineering, Tsinghua University, Beijing, P.R. China, 100084 ^{4} Griffith School of Engineering, Griffith University Gold Coast Campus, Queensland 4222, Australia ^{5} School of Civil and Transportation Engineering, Beijing University of Civil Engineering and Architecture, Beijing, P.R. China, 100044 Shear walls are important lateralforce resistant components of tall buildings. Hence, a reliable numerical model that can accurately represent the mechanical characteristics and large deformations of shear walls is critical for realistic collapse simulation of tall buildings. Based on the theory of generalized conforming element, a highperformance quadrilateral flat shell element, NLDKGQ, accounting for the large deformation using the updated Lagrangian formulation, is proposed herein and implemented in OpenSees. The reliability of NLDKGQ is validated using classical benchmark problems and reinforced concrete (RC) specimens. In addition, its capability in simulating the collapse of a tall building is also demonstrated. Key words: flat shell element; OpenSees; shear wall; large deformation; collapse simulation. 
1 IntroductionRecently, research on the seismic performance and collapse resistance of tall and supertall buildings has become increasingly popular in earthquake engineering. The importance of these research efforts is evident, given the rapid growth of tall and supertall building construction and an increased frequency of strong earthquakes. Research to date has indicated that numerical simulation using various finite element (FE) software packages is one of the most prevalent methods for the seismic performance analysis and collapse simulation of tall and supertall buildings [Fan et al., 2009; Lu et al., 2011, 2013a, 2013b, 2016]. One such open source FE program, OpenSees, has been widely used to investigate the seismic performance of tall and supertall buildings [Lu et al., 2015b], as well as super long span bridges [Lin et al., 2015]. The advantages of using OpenSees for seismic simulation have also been thoroughly discussed by Lu et al. [2015b]. OpenSees has also been adopted to conduct the collapse simulation of RC frames [Xie et al., 2015]. However, research on the collapse simulation of RC tall buildings using OpenSees has rarely been reported. Beams, columns and shear walls are the primary components of RC tall buildings. To realistically simulate the collapse of RC tall buildings, the numerical model is required to have the capacity to accurately and reliably simulate the collapse behavior (which exhibits significant material and geometric nonlinearities, concurrently) of these components. Collapse simulations of beams and columns have been successfully conducted using OpenSees [Chen et al., 2008; Lignos et al., 2011; Xie et al., 2015]. The material nonlinearities of the shear walls have also been successfully simulated by the authors using the shell elements in OpenSees [Lu et al., 2015b], based on a previously selfdeveloped multilayered shell section model. However, the said shell element, i.e., the MITC4 shell element in OpenSees, has limited capacity for large deformation problems (i.e. geometric nonlinearity). As a result, such a shell element in OpenSees cannot be used to simulate the large deformation behavior of shear walls during the process of seismic collapse. Note that the collapse simulation of simple shear wall structures has been conducted in OpenSees by using the beam¨Ctruss model [Lu et al., 2014] and the modified multiple vertical line element model [Bao and Kunnath, 2010; Fischinger et al., 2012]. However, Lu et al. [2015b] has emphasized that the shell element is advantageous in investigating the seismic performance of complicated shear wall structures. Hence, a new and improved shell element is necessary to be developed to account for large deformations associated with the seismic collapse. In addition to satisfying the demand of large deformations, several other demands are also required for the new shell element: (1) The shear walls in a tall building are usually presented with various types of sections and shapes with complicated openings. The new shell element is expected to be highly adaptive for different types of shear walls. In addition, the element must be insensitive to mesh distortion and free from locking problems. It should be noted that an ideal mesh is not always possible for actual shear walls with irregular shape and openings. (2) The geometric and material nonlinearities should be taken into consideration concurrently for the new shell element. This is critical for simulating the failure and collapse of shear walls. (3) A balance between the modeling workload and computational accuracy must be considered. Specifically, the MITC4 shell element has 5 degrees of freedom (DOFs) per node while the beam/column element has 6 DOFs per node. When a shell element is connected to the beam/column elements directly, an embedded beam element along the boundary of the MITC4 shell element is normally required [Lu et al., 2015b]. Such an embedded element strategy results in additional modeling and computational efforts. To address this problem, the shell element should also be assigned 6 DOFs per node, which can greatly simplify the connection modeling between the shell elements and the beam/column elements. (4) The new shell element should be adaptive to the overall framework of OpenSees, so that it can be implemented into OpenSees with minimal impact on the entire software. To meet the abovementioned demands, it is essential to choose an appropriate shell element type. In this study, the flat shell element is selected, owing to its satisfactory performance in the FE analysis [Long et al., 2009]. In addition, due to the simple and concise equations and geometric descriptions of the flat shell elements [Long et al., 2009], the large deformation algorithm can be conveniently developed and implemented. Several successfully developed flat shell elements [Madenci and Barut, 1994; Zhang et al., 1998; Zhang and Yang 2006] are proven to be free from mesh distortion and locking problems with acceptable modeling workload and satisfactory accuracy. However, several limitations of the currently available flat shell elements still remain: (1) Most of the flat shell elements are triangular in shape. However, quadrilateral ones are normally preferred in earthquake engineering research due to their high performance in computation and reduced workload in modeling. (2) Most of such elements have not taken the coupling effect of material and geometric nonlinearities into consideration. Therefore, a novel highperformance quadrilateral flat shell element, the NLDKGQ element, is proposed and implemented in OpenSees in this work. The NLDKGQ element consists of a planar membrane element GQ12 [Xu and Long, 1993] and a plate bending element DKQ [Batoz and Tahar, 1982]. The GQ12 element is a highperformance quadrilateral membrane element, based on the theory of generalized conforming element by introducing rotational degrees of freedom at each node [Xu and Long, 1993]. The DKQ element, which is immune from locking problems, is a quadrilateral plate element with high accuracy and stability [Batoz and Tahar, 1982]. Given that the GQ12 element has two translational (along x and y) and one rotational (around z) DOFs per node, and the DKQ element has one translational (along z) and two rotational (around x and y) DOFs per node, the combination of these two elements makes the NLDKGQ element having 6 DOFs per node. Consequently, the NLDKGQ elements can be readily connected to the beam elements without introducing extra embedded beam elements, thereby greatly reducing the modeling workload. The geometric nonlinearity is further considered using an updated Lagrangian formulation [Belytschko et al., 2013]. The reliability of the NLDKGQ element for elastic small deformation and large deformation problems is validated using classical benchmark examples. Subsequently, typical RC specimens are simulated using the NLDKGQ element to validate its reliability in representing the complex mechanical behavior of the shear walls. Finally, an incremental dynamic analysis (IDA) is performed for a 42story RC framecore tube building with a height of 141.9 m to demonstrate the capability of the element in the collapse simulation of tall buildings. The research outcome is expected to provide an effective tool for further indepth investigation of the collapse mechanisms of tall and supertall buildings using OpenSees. 
2 Formulation of the NLDKGQ element2.1 Basic formulation under small deformationThe quadrilateral flat shell element NLDKGQ in the local coordinate system is shown in Figure 1.
The definition of the nodal displacement vector q (Figure 1) is: _{ } (1) in which the nodal displacement components of the planar membrane element and the plate bending element are, respectively, represented by q_{i}^{m} and q_{i}^{b} , as follows: _{ } (2) _{ } (3) The displacement field u of GQ12 consists of the bilinear displacement field u_{0} and the additional rotation displacement field u_{¦È}, which can be determined from the interpolation functions N_{i}^{0}, N_{i}^{u}^{q}, and N_{i}^{v}^{q} [Xu and Long, 1993], as given in Equations 4 to 7. _{ } (4) _{ } (5) where _{ } (6) _{ } (7) Then, the strain field e^{m} can be obtained from the derivative of the displacement field. _{ } (8) where _{ } (9) The stiffness matrix of the planar membrane element GQ12 is: _{ } (10) where D_{mm }denotes the material matrix of the planar membrane element. For isotropic linear elastic materials, D_{mm }is given in Equation 11, where h denotes the thickness of the element, and E and u denote the elastic modulus and Poisson¡¯s ratio, respectively. _{ } (11) Equation 3 presents the definition of the rotational DOFs of the DKQ element [Batoz and Tahar, 1982]. The relationship between the rotational strain fields c_{b} and nodal displacements q^{b} is: _{ } (12) where H^{x} and H^{y} are functions of N_{i}, x_{i} and y_{i}. Details about Equation 12 are provided by Batoz and Tahar [1982]. The stiffness matrix of the plate bending element DKQ is: _{ } (13) where D_{bb} denotes the material matrix of the plate bending element. For isotropic linear elastic materials, D_{bb} has the following form. _{ } (14) Under the small deformation, the stiffness matrix of the NLDKGQ element in the local coordinate system can be obtained through the composition of the membrane and the plate stiffness matrices K_{m} and K_{b} following the DOFs sequencing shown in Equation 1. The element stiffness matrix in the global coordinate system can be subsequently obtained after the coordinate transformation. 2.2 Geometric nonlinearityBased on the deformation calculated at the current time step, the stresses and strains are updated in incremental forms using the updated Lagrangian formulation. According to Kirchhoff¡¯s and von Karman¡¯s assumptions [PodioGuidugli, 1989], the strain increment of the shell element (De) can be separated into a linear part (De) and a nonlinear part (Dh). The linear part can be obtained from the composition of the membrane strain increment De_{m }and the plate rotational strain increment Dc_{b}, which are presented in Equations 15 and 16. _{ } (15) _{ } (16) From the current time step t to t+dt, the stress tensor of the shell element is updated according to Equation 17: _{ } (17) where D_{tan} is the tangential constitutive matrix of the material at time t. Using the updated Lagrangian formulation, the system of equations of the shell element in the local coordinate system is: _{ } (18) where _{t}_{+dt}F denotes the external force vector at t+dt, and _{t}R denotes the internal force vector at t. The linear part K_{l} and the nonlinear part K_{nl} of the tangential stiffness matrix of the element are given in Equations 19 and 20. _{ } (19) _{ } (20) where D_{mm}, D_{mb}, D_{bm}, and D_{bb} can be obtained through the integration of D_{tan} over the thickness of the element (Equation 21). The matrix G can be obtained from the interpolation function of the bending plate element according to its definition [Batoz and Tahar, 1982]. The matrix _{ } is composed of variables corresponding to the internal force vector of the membrane element at time t (Equation 22). _{ } (21) _{ } (22) The elemental internal force vector _{t}R given in Equation 18 can be calculated as follows: _{ } (23) where _{ } (24) _{ } (25) 2.3 Implementation in OpenSeesThe formulation of the new shell element is added into OpenSees through the new class named ShellNLDKGQ, which is under the class of shell. The new element perfectly fits the framework of OpenSees without changing any source code beyond the shell element domain. The source code of the element can be downloaded from the official website of OpenSees (http://opensees.berkeley.edu). 
3 Validation through classical benchmark problemsTo validate the performance of the flat shell element NLDKGQ, several classical numerical benchmark problems are analyzed. For convenience of discussion, all parameters of these problems are made nondimensional. 3.1 Pinched cylinder problemThe pinched cylinder problem, shown in Figure 2 [Belytschko and Leviathan, 1994], has been widely used to benchmark the performance of various shell elements. Due to the symmetrical conditions, only oneeighth of the cylinder is modeled, meshed and analyzed. The corresponding results are shown in Table 1. In this problem, three types of shell elements are compared, namely the 4node quadrilateral element S4 in ABAQUS [Dassault Syst¨¨mes Simulia Corp., 2009], MITC4 in OpenSees [Dvorkin and Bathe, 1984], and the proposed flat shell element NLDKGQ. For all three types of shell elements, the geometric nonlinearity is not considered. Table 1 shows that the accuracy of the NLDKGQ element is higher than that of the other two shell elements, especially when the mesh is coarse. This advantage of the NLDKGQ element is very critical for the simulation of tall buildings consisting of thousands of structural components, because a very fine mesh requiring large amount of computational workload is impractical and unacceptable.
Table 1. Normalized results for the pinched cylinder problem
^{a} The theoretical solution is 1.8541¡Á10^{5} [Belytschko and Leviathan, 1994] 3.2 MacNeal¡¯s beam problemThe MacNeal¡¯s beam problem [MacNeal and Harder, 1985] is a classical benchmark to test the mesh distortion sensitivity of 4node quadrilateral elements. Three different meshes are adopted: (a) rectangular, (b) parallelogram and (c) trapezoid, as shown in Figure 3. Pure bending and transverse loading cases are analyzed. Geometric nonlinearity is not considered in this case. In addition to the NLDKGQ elements, the normalized results obtained from other types of elements are also given for comparison, as summarized in Table 2. Both the NLDKGQ and the S4 elements are free from mesh distortion and membrane locking problems and their predictions are of similar accuracy; whereas serious membrane locking problems are encountered by the MITC4 elements in OpenSees and the Q4 [Shi and Voyiadjis, 1991]. Note that a similar phenomenon associated with the MITC4 element was also found by Cardoso et al. [2006]. This comparison confirms the advantage of NLDKGQ in simulating the beams with irregular meshes.
3.3 Large deformation problem of a cantilever beamA cantilever beam subjected to a pure bending load (outofplane) is analyzed to validate the capacity of the NLDKGQ element in modeling large deformations (geometric nonlinearity), as shown in Figure 4. The commonly used mesh of 1´10 is adopted [Horrigmoe and Bergan, 1978; Park et al., 1995]. Figure 5(a) shows the relationship between the normalized moment (k=M/M_{max}) and the horizontal and vertical displacements at the loading point. Figure 5(b) shows the deformed shape of the cantilever beam under different bending moments. The results show that the NLDKGQ element can simulate the large deformation and rotation problems with good accuracy, which is similar to the S4 element in ABAQUS. Such a large deformation capacity makes the NLDKGQ element highly suitable for the collapse simulation of tall buildings. Note that the MITC4 element in OpenSees fails to simulate this problem accurately even when k is small, due to its limited capacity in modeling large deformations. It should also be noted that the same model using the MITC4 element fails to converge when k reaches 0.4485. Table 2. Normalized results for the MacNeal¡¯s slender beam problem
^{a} The theoretical value is 0.1081. ^{b} The theoretical value is 0.0054.

4 Validation through analysis of RC specimens4.1 RC shear wall experimentsTo investigate the performance of the NLDKGQ element in simulating RC specimens, the multilayered shell section model [Lu et al., 2015b] is used herein. Details of this section model in OpenSees are provided by Lu et al. [2015b] and the corresponding command lines are given in Appendix B. In general, concrete is assumed to be in a state of planar stress in the multilayered shell element. The analytical model for concrete is based on the concept of damage mechanics and the smeared crack model. The reinforcing bars are simulated as smeared steel layers of equivalent thicknesses. The material model for reinforcing steel, designated as PlateRebar, is derived from the uniaxial materials for rebars in conjunction with their corresponding angles with respect to the Cartesian coordinate. Based on the multilayered shell section model [Lu et al., 2015b] and the NLDKGQ element, four shear wall specimens are analyzed using OpenSees to investigate the performance of the NLDKGQ element in simulating the hysteretic behavior of the RC shear walls. The test specimens include two rectangular walls [Zhang, 2007], one flanged wall [Chen and Qian, 2005], and one coupled wall [Chen and Lu, 2003]. The element mesh recommended by Lu et al. [2015b] is used. The simulation results predicted by the MITC4 element are also presented for comparison. The lateral force versus top displacement hysteretic curves shown in Figure 6 indicate that both the NLDKGQ and the MITC4 elements perform satisfactorily against the test results.
4.2 A pseudostatic collapse experiment of an RC columnMost RC shear wall experiments are limited to small deformations. To validate the performance of the NLDKGQ element for RC components under large deformations, a pseudostatic collapse experiment of an RC column [Xie et al., 2015] is simulated herein. Column B tested by Xie et al. [2015] is modeled with the MITC4 and NLDKGQ elements using the same mesh strategy and same multilayered shell section. Figure 7 shows the loading condition of Column B. The column is meshed with 15 (along height) ´ 4 (along width) elements. Figure 8 shows the comparison of the hysteretic curves. It is evident that, the MITC4 element fails to simulate the collapse characteristics of Column B because of its locking problems and limited capacity to simulate large deformations. In contrast, the NLDKGQ element can accurately simulate the strength degradation due to large deformations.

5 IDA of an RC framecore tube tall buildingAs described in Sections 3 and 4, the accuracy and reliability of the proposed NLDKGQ element have been validated. In this section, the proposed element is used for the collapse simulation of a typical tall building. In this work, a 42story RC framecore tube building (Figure 9) with a height of 141.9 m (denoted as Building 2N by Lu et al. [2015a]) is simulated using OpenSees. More details about this building can be found in Lu et al. [2015a]. Note that the element mesh of this building is further refined compared with that in Lu et al. [2015a] to better represent the local failure of shear walls during the collapse procedure. The IDA of this tall building is performed to demonstrate the collapse simulation capacity of the NLDKGQ element. Two FE models are established for comparison, i.e., ¡°Building 2NNLDKGQ¡± (the shear wall is meshed with the NLDKGQ element) and ¡°Building 2NMITC4¡± (the shear wall is meshed with the MITC4 element). In general, it is rather challenging to establish a complicated tall building model in OpenSees directly. Hence a conversion program, named MSC.Marc to OpenSees [Lu et al., 2015b], is used herein. The modeling strategy is as follows. The DispBeamColumnelement in OpenSees is adopted to simulate the beams and columns using fiber sections. The multilayered shell elements (i.e. ShellMITC4 and ShellNLDKGQ) are used to simulate the shear walls, core tubes and coupling beams. The rigidDiaphragm is used to simulate the slabs. The Rayleigh damping consisting of mass and current stiffness matrix is adopted. The analysis domain recommended by Lu et al. [2015b] is adopted, as given in Table 3. A comparison between the fundamental periods predicted by using the MITC4 and NLDKGQ elements is presented in Table 4. The comparison indicates that both models exhibit similar basic dynamic properties at the elastic stage.
Table 3. Analysis domain used for Building 2N
Table 4. Comparison of the periods predicted by the MITC4 and NLDKGQ elements
Note that the results of modal analysis in Table 4 are lower than the results of a similar US design (denoted as ¡°Building 2A¡±) provided by Moehle et al. [2011]. There are two main reasons for the shorter design period of Building 2N: (1) Moehle et al. [2011] used effective component stiffness values to predict the vibration period; while the elastic stiffness values are used herein following the Chinese design codes. (2) The basement is included in the analysis model of Building 2A [Moehle et al., 2011], which also lengthened the design periods. Regardless of the above two reasons, if the periods of Building 2A are calculated using the same method for Building 2N without including the basement, the elastic fundamental period of Building 2A is approximately 2.9 s. This value is similar to the translational period of the models established in this work. More details have been discussed by Lu et al. [2015a]. For tall buildings, highorder vibration modes in general contribute largely to the seismic responses. Therefore, performing a pushover analysis for tall buildings may not be considered appropriate [Krawinkler and Seneviratna, 1998; Fajfar, 2000; Chopra and Goel, 2002; Jan et al., 2004]. The nonlinear time history analysis, which is more suitable for tall buildings, is thus conducted in this study. Note that the PDelta effect is also considered in the analysis. To validate the models in OpenSees, a careful comparison of the simulation results predicted by OpenSees and a well validated commercial software MSC.Marc [MSC Software, 2007] is conducted, and the results are shown in Figure 10. In the time history analysis, the widely used ElCentro 1940 ground motion is input to the structure along the Xdirection after the gravity analysis is completed. Two peak ground accelerations (PGA) are chosen here, including 110 cm/s^{2} and 510 cm/s^{2} (63.2% and 2% probability of exceedance at the 50year hazard level of this building defined in the Chinese code [CMC, 2010]). It can be seen that, the results agree well at both elastic and nonlinear stages, confirming the reliability of the models being implemented in OpenSees. In the IDA, the ElCentro 1940 ground motion is also used as input to the structure along the Xaxis after the gravity analysis is completed. Several typical PGAs are considered, including 510 cm/s^{2}, 2000 cm/s^{2}, 4000 cm/s^{2}, 5000 cm/s^{2}, and 6000 cm/s^{2}. Note that such an IDA for tall buildings is very time consuming if it is implemented with the CPUbased matrix solvers in OpenSees. Thus, a graphics processing unit (GPU)based matrix solver (i.e., solver CuSP) is developed to accelerate the simulation process. Details of the GPUbased solver are presented in Appendix A. It should also be noted that, in this study, elemental deactivation technology is not used. In the collapse simulation, the stress and stiffness of the materials (steel and concrete) are considered to decrease with increased deformation after the peak point, following the material constitutive law. In consequence, when the deformation is large enough, the stress and stiffness of the materials are approaching zero. For instance, the concrete strength reaches zero when the strain of concrete fibers/layers is greater than the ultimate strain. Such modeling technique has proven to successfully simulate the collapse behavior together with the geometric nonlinearity of the finite element model. The envelopes of the interstory drifts at different PGAs are shown in Figure 11 for Building 2NNLDKGQ and Building 2NMITC4. The relationship between the PGA and the corresponding maximum interstory drift ratio is illustrated in Figure 12. With respect to the interstory drift ratios, Villaverde [2007] summarized a number of methods in assessing the seismic collapse capacity of buildings and concluded that in the incremental dynamic analysis (IDA), ¡°the global collapse capacity is considered reached when the curve in this plot becomes flat¡±. In general, 20% of the initial slope is adopted as the criterion for defining the ¡°flat¡± curve [FEMA, 2000; Jalayer, 2003; Riahi et al., 2015]. From the IDA curve shown in Figure 12, this criterion is satisfied when the interstory drift ratio is greater than 7%. The corresponding collapse mode of Building 2NNLDKGQ at the PGA level of 6000 cm/s^{2} is displayed in Figure 13.
The comparisons given in Figures 11 and 12 reveal that Building 2NNLDKGQ and Building 2NMITC4 present similar results when the displacement is small (i.e., the maximum interstory drift is smaller than 2% when PGA is smaller than 4000 cm/s^{2}). However, a significant difference can be observed when the structures undergo large deformations and when collapse occurs under the strong ground motion (i.e., PGA is greater than 4000 cm/s^{2}). Considerable deformations taking place at the bottom 5 stories of Building 2NNLDKGQ (Figure 13) signify the collapse of the building under large deformations. Such a collapse behavior cannot be predicted by the MITC4 element, as indicated by the continuously increasing IDA curve (Figure 12) of Building 2NMITC4, even under the excessive ground motion. This outcome is no doubt unreasonable.

6 ConclusionsIn this work, a highperformance quadrilateral flat shell element NLDKGQ is proposed based on the theory of generalized conforming element and updated Lagrangian formulation. Through the analysis of the classical benchmark problems and the typical experiments of RC specimens, the performance of the proposed element is validated. Furthermore, the collapse simulation of an RC framecore tube building with a height of 141.9 m is successfully conducted using the proposed element. The following conclusions are drawn. (1) The NLDKGQ element is perfectly adaptive to the framework of OpenSees and can be conveniently added to the software. (2) The NLDKGQ element is insensitive to mesh distortion and free from locking problems. In addition, it can concurrently consider the geometric and material nonlinearities, and produce accurate solutions to geometric and material nonlinear problems. (3) In combination with the multilayered shell section model, the NLDKGQ element can accurately simulate the nonlinear behavior of RC specimens and perform better than the MITC4 element in OpenSees when the RC specimens undergo large deformations. (4) The NLDKGQ element can be used to simulate large deformations and collapse behavior of typical RC tall buildings. This simulation cannot be performed in OpenSees without the NLDKGQ elements. The highperformance quadrilateral flat shell element NLDKGQ proposed in this work is expected to provide an effective tool and a useful reference for further indepth investigation of the collapse mechanisms of tall and supertall buildings using OpenSees. 
AcknowledgementsThe authors would like to acknowledge the financial supports of the National Natural Science Foundation of China (Nos. 11272181, 51261120377), the Specialized Research Fund for the Doctoral Program of Higher Education of China (No. 20120002110080), and the Tsinghua University Initiative Scientific Research Program (Project No. 2014z09099). The authors are also grateful for Prof. Quan Gu for his help in adding the source code of this new element into OpenSees. ReferencesBao, Y. and Kunnath, S. K. [2010] ¡°Simplified progressive collapse simulation of RC framewall structures,¡± Engineering Structures 32(10), 31533162. Batoz, J. L. and Tahar, M. B. [1982] ¡°Evaluation of a new quadrilateral thin plate bending element,¡± International Journal for Numerical Methods in Engineering 18(11), 16551677. Belytschko, T. and Leviathan, I. [1994] ¡°Physical stabilization of the 4node shell element with one point quadrature,¡± Computer Methods in Applied Mechanics and Engineering 113(3), 321350. Belytschko, T., Liu, W. K., Moran, B. and Elkhodary, K. [2013] Nonlinear finite elements for continua and structures. John Wiley & Sons. Cardoso, R. P. R., Yoon, J. W., and Valente, R. A. F. [2006] ¡°A new approach to reduce membrane and transverse shear locking for onepoint quadrature shell elements: linear formulation,¡± International Journal for Numerical Methods in Engineering 66(2), 214249. Chen, Q. and Qian, J. R. [2005] ¡°Static elastoplastic analysis of RC shear walls with one row opening,¡± Chinese Journal of Computational Mechanics 22(1), 1319. [in Chinese] Chen, X. W., Han, X. L., Jack, C., Lin, S.Y. and Mao, G.N. [2008] ¡°Dynamic inelastic numerical simulation for a shaking table test of a full scale steel moment frame structure based on OpenSEES,¡± Proc. of the 14th World Conference on Earthquake Engineering, Beijing, China. Chen, Y. T. and Lu, X. L. [2003] ¡°Seismic behavior of coupled shear wallsexperiment and theoretical analysis,¡± Journal of Building Structures 24(4), 2534. [in Chinese] Chopra, A. K., and Goel, R. K. [2002] ¡°A modal pushover analysis procedure for estimating seismic demands for buildings,¡± Earthquake Engineering & Structural Dynamics 31(3), 561582. CMC. [2010] Code for Seismic Design of Buildings (GB500112010), China Ministry of Construction; China Architecture and Building Press, Beijing, China. [in Chinese] Dassault Syst¨¨mes Simulia Corp.. [2009] Abaqus Theory Manual. Abaqus 6.9 HTML Documentation. http://abaqusdoc.ucalgary.ca/v6.9/books/stm/default.htm Dvorkin, E. N. and Bathe, K. J. [1984] ¡°A continuum mechanics based fournode shell element for general nonlinear analysis,¡± Engineering computations 1(1), 7788. Fajfar, P. [2000] ¡°A nonlinear analysis method for performancebased seismic design,¡± Earthquake spectra 16(3), 573592. Fan, H., Li, Q. S., Tuan, A. Y. and Xu, L. [2009] ¡°Seismic analysis of the world¡¯s tallest building,¡± Journal of Constructional Steel Research 65, 12061215. FEMA. [2000] ¡°Recommended seismic design criteria for new steel moment frame buildings,¡± Publ. No. 350, Washington, D.C. Fischinger, M., Rejec, K. and Isaković, T. [2012] ¡°Modeling inelastic shear response of RC walls,¡± Proc. of the 15th World Conference on Earthquake Engineering, Lisboa. Horrigmoe, G. and Bergan, P. G. [1978] ¡°Nonlinear analysis of freeform shells by flat finite elements,¡± Computer Methods in Applied Mechanics and Engineering 16, 1135. Jan, T. S., Liu, M. W., and Kao, Y. C. [2004] ¡°An upperbound pushover analysis procedure for estimating the seismic demands of highrise buildings,¡± Engineering Structures 26(1), 117128. Jalayer, F. [2003] ¡°Direct probabilistic seismic analysis: implementing nonlinear dynamic assessments,¡± PhD. Dissertation, Department of Civil Engineering, Stanford University. Krawinkler, H., and Seneviratna, G. [1998] ¡°Pros and cons of a pushover analysis of seismic performance evaluation,¡± Engineering structures 20(4), 452464. Lignos, D. G., Chung, Y., Nagae, T. and Nakashima, M. [2011] ¡°Numerical and experimental evaluation of seismic capacity of highrise steel buildings subjected to long duration earthquakes,¡± Computers & Structures 89(11), 959967. Long, Y. Q., Cen, S. and Long, Z. F. [2009] Advanced finite element method in structural engineering, Tsinghua University Press. Lin, K. Q., Xie, L. L., Lu, X. Z. and Ye, L. P. [2015] ¡°Earthquakeinduced collapse simulation of a superlong span cablestayed bridge based on an open source FE program,¡± Proc. of IABSE ConferenceStructural Engineering: Providing Solutions to Global Challenges, Geneva, Switzerland. Lu, X. Z., Lu, X., Guan, H., Zhang, W.K. and Ye, L. P. [2013b] ¡°Earthquakeinduced collapse simulation of a supertall megabraced framecore tube building,¡± Journal of Constructional Steel Research 82, 5971. Lu, X., Lu, X. Z., Guan, H., and Xie, L. L. [2016] ¡°Application of earthquakeinduced collapse analysis in design optimization of a supertall building,¡± The Structural Design of Tall and Special Buildings 25(17), 926¨C946. DOI: 10.1002/tal.1291. Lu, Y., Panagiotou, M. and Koutromanos, I. [2014] ¡°Threedimensional beamtruss model for reinforcedconcrete walls and slabs subjected to cyclic static or dynamic loading,¡± Pacific Earthquake Engineering Research Center, Los Angeles, California. MacNeal, R. H. and Harder, R. L. [1985] ¡°A proposed standard set of problems to test finite element accuracy,¡± Finite Elements in Analysis and Design 1(1), 320. Madenci, E. and Barut, A. [1994] ¡°A freeformulationbased flat shell element for nonlinear analysis of thin composite structures,¡± International Journal for Numerical Methods in Engineering 37(22), 38253842. Moehle, J., Bozorgnia, Y., Jayaram, N., Jones, P., Rahnama, M., Shome, N., Tuna, Z., Wallace, J., Yang, T., and Zareian, F. [2011] ¡°Case Studies of the Seismic Performance of Tall Buildings Designed by Alternative Means,¡± Pacific Earthquake Engineering Research Center, University of California, Berkeley, California, USA, July. MSC Software. [2007] Marc 2007 r1 User's Guide, MSC Software, Santa Ana, CA. NVIDIA. [2014a] CuSP Home Page. http://cusplibrary.github.io/. NVIDIA. [2014b] CUDA C programming guide, http://docs.NVIDIA.com/cuda/pdf/CUDA_C_Programming_Guide.pdf; Date accessed: May 2016. Park, H. C., Cho, C. and Lee, S. W. [1995] ¡°An efficient assumed strain element model with six DOF per node for geometrically nonlinear shells,¡± International Journal for Numerical Methods in Engineering 38(24), 41014122. PEER. [2014] Subversion Repositories, http://opensees.berkeley.edu/WebSVN/listing.php?repname=OpenSees&path=%2Ftrunk%2FSRC; Date accessed: May 2016. PodioGuidugli P. [1989] ¡°An exact derivation of the thin plate equation,¡± Journal of Elasticity 22(23), 121133. Riahi, H.T., Amouzegar, H., and Falsafioun, M. [2015] ¡°Seismic collapse assessment of reinforced concrete moment frames using endurance time analysis,¡± The Structural Design of Tall and Special Buildings 24(4), 300315. Shi, G. and Voyiadjis, G. Z. [1991] ¡°Efficient and accurate fournode quadrilateral C0 plate bending element based on assumed strain fields,¡± International Journal for Numerical Methods in Engineering 32(5), 10411055. Villaverde, R. [2007] ¡°Methods to assess the seismic collapse capacity of building structures: State of the art,¡± Journal of Structural Engineering, ASCE, 133(1), 5766. Xu, Y. and Long, Y. Q. [1993] ¡°Quadrilateral membrane element with vertex rotational freedom from generalized compatible condition,¡± Engineering Mechanics 10(3), 2736. [in Chinese] Zhang, H. M. [2007] ¡°Study on the performancebased seismic design method for shear wall structures,¡± Ph.D. thesis, Tongji University, Shanghai, China. [in Chinese] Zhang, Q., Lu, M. and Kuang, W. Q. [1998] ¡°Geometric nonlinear analysis of space shell structures using generalized conforming flat shell elementsfor space shell structures,¡± Communications in Numerical Methods in Engineering 14(10), 941957. Zhang, Y. X. and Yang, C. H. [2006] ¡°A family of simple and robust finite elements for linear and geometrically nonlinear analysis of laminated composite plates,¡± Composite Structures 75(1), 545552. Appendix A. GPUbased Matrix Solvers in OpenSeesA new GPUbased parallel solver for the sparse systems of equations (SOEs) is developed and added to OpenSees, namely CuSP solver [Xu et al., 2016]. The matrix algorithm of CuSP solver is based on CuSP [NVIDIA, 2014a], a GPUaccelerated numerical library for the sparse SOEs developed by Compute Unified Device Architecture (CUDA) [NVIDIA, 2014b]. The CuSP solver is developed based on the following rules: (1) The iterative algorithms (i.e., the conjugate gradient (CG) algorithm, BiCG algorithm and generalized minimal residual (GMRES) algorithm) are used to maximize the parallel computing ability of the GPU. (2) The Solver class (corresponding to the LinearSOESolver class of OpenSees) and the solving function are designed separately. The Solver class is designed inheriting the LinearSOESolver class, to improve the compatibility of the solver. (3) The kernel algorithm of CuSP solver is packaged in the form of a dynamiclink library (DLL), which is convenient to be updated without changing the framework of OpenSees. Figure A.1 illustrates the specific solving process of CuSP solvers. The matrix is integrated in CPU, and then copied into the graphics memory to perform the parallel computation. Finally, the results are returned into CPU for subsequent computation. The corresponding source codes of CuSP solver are available at the website of OpenSees [PEER, 2014].
Appendix B. Command lines to implement the multilayered shell element and NLDKGQ elementThe twodimensional concrete material is introduced to the model through the following command lines as
The reinforcing steel material for the multilayered shell element is introduced to the model through the following command line as
The multilayered shell section is defined through the following command line as
The NLDKGQ shell element is defined through the following command line as
* Address correspondence to Xinzheng Lu, Department of Civil Engineering, Tsinghua University, Beijing, P.R. China. Email: luxz@tsinghua.edu.cn 