A non-linear response history model for the seismic analysis of high-rise framed buildings

Computers and Structures 84 (2006) 318–329 www.elsevier.com/locate/compstruc

A non-linear response history model for the seismic analysis of high-rise framed buildings S.M. Wilkinson *, R.A. Hiley School of Civil Engineering and Geosciences, University of Newcastle, Drummond Building, Newcastle upon Tyne NE1 7RU, United Kingdom Received 27 October 2004; accepted 29 September 2005

Abstract A materially non-linear plane-frame model is presented that is capable of analysing high-rise buildings subjected to earthquake forces. The model represents each storey of the building by an assembly of vertical and horizontal beam elements The model introduces yield hinges with ideal plastic properties in a regular plane frame. The displacements are described by the translation (sway) of each floor and the rotation of all beam–column intersections. The mass is only associated with the translations, and thus the analysis can be carried out as a static condensation of the rotations, combined with integration of the dynamic equations for the translations. The dynamic integration is here carried out by use of the Runge–Kutta scheme. This approach allows a building to be modelled by m(n + 2) degrees of freedom (where m is the number of storeys and n is the number of bays). The rank of the condensed stiffness matrix is only m. Its construction, which requires the inversion of the rotational, rank m(n + 1), stiffness matrix, is required only at time-steps where the pattern of yielding has altered from the previous time-step. This model is particularly attractive for non-linear response history analysis of high-rise buildings as it is efficient, allows each storey to have multiple redundancies, and each connection to be modelled with any suitable moment–rotation relationship. Three verification examples are presented and the results from static push-over analysis are compared with time–history results from the simplified model. The results verify that the model is capable of performing non-linear response history analysis on regular high rise buildings.
2005 Elsevier Ltd. All rights reserved. Keywords: Dynamic analysis; Earthquake engineering; Multi-storey buildings; Response history analysis; Push-over analysis; Non-linear analysis

1. Introduction Because of the nature of earthquakes, a dual design philosophy has been adopted for the design of buildings in earthquake prone regions. The first design criterion is to ensure that little or no damage is suffered during an earthquake that can reasonably be expected to occur during the life-time of the structure. The second is that the building does not collapse during the most severe probable earthquake that could occur at that site. The corollary of this is that if the building is to remain cost effective the second

*

Corresponding author. Fax: +44 191 222 6613. E-mail address: [email protected] (S.M. Wilkinson).

0045-7949/$ - see front matter
2005 Elsevier Ltd. All rights reserved. doi:10.1016/j.compstruc.2005.09.021

criterion will make it necessary to design the building inelastically. It is for this reason that all buildings designed in regions where earthquakes pose a serious threat to infrastructure are in some way designed inelastically. At present three main methods are used to analyse buildings subjected to earthquakes. These are: 1. Response history analysis. 2. Response spectrum analysis. 3. Quasi-static method. Response history analysis is potentially the most accurate but there are two problems associated with it. The first is that it can be difficult to choose an appropriate earthquake to use as the loading, while the second is that it is

S.M. Wilkinson, R.A. Hiley / Computers and Structures 84 (2006) 318–329

generally too computer-intensive to be practical—especially if inelastic analysis is considered (see e.g. [1]). The computer resources required to perform a response history analysis on a detailed inelastic finite element model are generally considered prohibitive. The most commonly employed method is the quasistatic method, as it is the simplest, requires only static analysis, and estimates the response of the structure for an ensemble of earthquakes. The response spectrum method is identical to the quasistatic method except that it considers more than just the fundamental mode of vibration. Most codes require that enough modes of vibration are considered to account for 90% of the modal mass (EC8: 1994 [2], IBC: 2000 [3], SANZ: 1992 [4]). For the quasi-static method and the response spectrum method the earthquake forces are divided by a behaviour factor (also known as a structural response factor or response modification coefficient). This factor accounts for the reserve strength of the building after the formation of the first plastic hinge and allows a pseudo inelastic design to be achieved without complicating the analysis. The only extra requirement to account for inelastic behaviour is for the designer to choose an appropriate building behaviour factor. Typically, this is done by choosing a value from a table in a relevant earthquake code. This is simple and reasonably effective but it is overly conservative. The various ductility factors have been arrived at empirically based on past experience of structural behaviour during earthquakes and based on generalised analysis of simple models of various building types. Recently EC8:1994 [2] and the New Zealand earthquake code SANZ:1992 [4] has tried to improve the selection of the behaviour factor by introducing Ôcapacity designÕ. The philosophy adopted in capacity design is best outlined by Park [5] ‘‘In the capacity design of structures, appropriate regions of the primary lateral earthquake force resisting structural system are chosen and suitably designed and detailed for adequate strength for a severe earthquake. All other regions of the structural system and other possible failure modes, are then provided with sufficient strength to ensure that the chosen means for achieving ductility can be maintained throughout the post-elastic deformations that may occur.’’ With certain limitations, this procedure could allow designers to choose whatever behaviour factor they want as long as they then ensure that the elements and connections have sufficient rotational capacity to redistribute the forces and hence achieve the ultimate load. However, because of the simplifications used in the analysis it is not a straightforward task to relate the actual ductility of a structure to a building behaviour factor. The ductility of a structure can be defined as the ratio of ultimate load to the load producing the first plastic hinge and is generally referred to as the ductility factor, which is different from a behaviour factor. Various researchers have tried to relate behaviour factors to ductility factors, notably Uang [6] and Whittaker et al. [7], but there is still no exact way to do this.

319

One of the ways of judging whether a building has sufficient ductility to resist the ultimate earthquake is to perform a push-over analysis. Push-over analysis subjects an inelastic building model to a static loading. This load is increased until a collapse mechanism is formed. The rotations of all the plastic hinges are calculated and the joints detailed to ensure that these rotations can be achieved. Typically designers do not actually perform an inelastic analysis but rather use the superposition of successive elastic analyses, placing pinned connections at the positions of the plastic hinges. Although effective, for a multi-storey building this procedure can be tedious and some doubt has been raised about its accuracy [8]. The inelastic capacity of a building will depend on the distribution of statically equivalent earthquake forces. These forces are applied up the building in proportion to the mass of each storey and the fundamental mode shape of the building. When designing high-rise buildings it is often necessary to consider more modes than just the fundamental in order to account for 90% of the modal mass. In such cases the inelastic capacity of a building may be misrepresented by a push-over analysis based only on the fundamental mode. The validity of inelastic response spectrum methods is still in doubt. Berg [9] quotes an example of non-conservative error of 300% when response spectrum analysis is used. Part of the reason for this is that the response spectrum method uses mode superposition to combine the responses of the individual modes. Mode superposition is only applicable to linear analysis and so the method of using building response factors is not strictly correct. At present this is not a great a problem as the building response factors are crude and conservative; however, as capacity design becomes more widely used the effects of inelastic behaviour may become more important. A recent improvement to this technique, known as modal push-over analysis (which accounts in an approximate way for the effects of yielding) has been developed and evaluated by Chintanapakdee and Chopra [10]. The combination of modal responses remains problematic, and results suggest that significant errors may arise in the analysis of tall and/or reduced-strength frames. The method is evaluated by comparing results with those of (non-linear) response history analysis, which evidently remains our most reliable analytical tool. Other recent developments of the push-over technique include that of Kim et al. [1] where the procedure is enhanced by considering more than just the fundamental mode and recalculating mode shapes whenever yielding occurs. 2. Model development Many simple models have been developed to perform response history analysis in the inelastic range. To combine the best of all models some researchers have looked at stick models. These models replace each storey of a building

320

S.M. Wilkinson, R.A. Hiley / Computers and Structures 84 (2006) 318–329

with a single element. This makes the problem sufficiently simple to allow efficient response history analysis. As plastic hinges are formed and unformed in the building the stiffness of each element can be altered to reflect the resulting change in stiffness. Thambiratnam and Thevendran [11] proposed a 5 degree of freedom model to analyse highrise buildings. Wilkinson and Thambiratnam [12] extended this model to consider more accurately the rotations at the columns. The stiffness of each floor was modified by a factor to account for the contribution of horizontal elements and the configuration and variation of columns. The model was extended to inelastic analysis, by Wilkinson and Hiley [13], by adding a torsional spring (with stiffness Kb) to the top and bottom of each column before assembling the stiffness matrix, instead of a modification factor being added to the final stiffness matrix. These torsional springs were successively modified as progressive yielding occurred. The problem with these models is that they lump the behaviour of all elements in a floor into a single element which approximates the ensemble behaviour of individual elements. For a small reduction in computational efficiency it is possible to provide a (static) degree of freedom for each connection and thus model the hysteretic behaviour of each connection in detail, individually. Krenk et al. [14] developed an efficient collapse analysis procedure for use in the offshore structures which incorporates yield hinges at the end of beam elements, thus permitting elimination of the additional plastic rotations from the global equations. A similar approach is used here in the context of high-rise buildings. For the case of high-rise buildings, the floors of the buildings are assumed to be rigid diaphragms enabling the lateral displacement of each floor to be modelled by only one degree of freedom. In the model presented in this paper, the idealised frame is represented as an assemblage of beam finite elements. Static condensation of the structure stiffness matrix, whereby static equilibrium is imposed for the rotational degrees of freedom, together with neglect of axial deformations allows the equations of dynamic equilibrium to be expressed in terms of one degree of freedom per storey. The model thus retains the computational efficiency of stick models. The advantage of the new method is that the inelastic behaviour of each individual element may be modelled to any degree of complexity. The static condensation is based on an instantaneous linearization of the element moment–rotation relationships, and must be performed whenever there is a change in that linearization (i.e. in the elasto-plastic case, whenever there is a change in the pattern of yielded joints). Only at these discrete instants does the full (rank m(n + 2)) stiffness matrix enter the calculations; otherwise the Runge–Kutta solution proceeds efficiently and effectively without any need for matrix inversion. The numerical integration scheme adopted is also amenable to error-based adaptive time-step control. At present the model is 2D, but could be extended to three dimensions. The simpler two dimensional case is investi-

gated in this paper to demonstrate that deficiencies occur in the present analysis methods even for simple 2D frames. 3. Structural model The idealized structure is a plane frame, with m floors and n bays, as shown (for m = n = 2) in Fig. 1. Column spacing is invariant throughout, and there are (n + 1) columns per floor. Axial degrees of freedom (DOFs) are neglected in both beams and columns. Thus the model has (n + 2) DOFs per floor. In the sequel: i = 1, . . . , m is the floor/storey index; j = 0, . . . , n is the column (or bay, if j = 1, . . . , n) index; k = 0, 1 is the beam-end index; and A(1,2) refers to element (1, 2), for example, of any matrix A. 3.1. Members The stiffness matrix for a column (superscript ÔCÕ), subject to a known axial compression force P (due in our case to gravity), is given by Krenk [15]. 2 3 12ui;j =L2i;
6=Li;
12ui;j =L2i;
6=Li; 6 7 3 þ u0i;j 6=Li; 3
u0i;j 7 wi;j EI Ci;j 6
6=Li; C 6 7; Ki;j ¼ 2 Li; 6 12ui;j =L2i; 6=Li; 7 4
12ui;j =Li; 6=Li; 5 0 0
6=Li; 3
ui;j 6=Li; 3 þ ui;j ð1Þ where EI is flexural rigidity; L is member length; u  a cot a and w  13 a2 =ð1
uÞ are the symmetric and anti-symmetric pffiffiffiffiffiffiffiffiffiffiffi bending stiffness coefficients; u 0  u/w; and a  12 L P =EI . The bending stiffness coefficients tend to unity as P tends to zero, and are defined as such. The corresponding vector of displacements is ½ ui
1 hi
1;j ui hi;j T , where u is storey ceiling-level displacement (relative to the ground) and h is column upper end-rotation. The elastic stiffness matrix for a beam (superscript ÔBÕ), neglecting shear in accordance with neglect of column axial DOFs, is given by

IB21 IC21

IB22 IC22

IB11 IC11

LC2

IC13

LC1

IB12 IC12

LB1

IC23

LB2 Fig. 1. Structural model.

S.M. Wilkinson, R.A. Hiley / Computers and Structures 84 (2006) 318–329

KBi;j

  2EI Bi;j 2 1 ¼ ; L;j 1 2

h with corresponding vector of displacements hBi;j;0

ð2Þ iT hBi;j;1 ,

where hBi;j;k is the rotation at end k of beam (i, j). 3.2. Yielding The moment rotation relationship is elasto-plastic. Inelastic behaviour of beams is modelled by allowing plastic deformation of the beam/column connections (i.e. the beam ends): the connection at end k of beam (i, j) yields whenever the magnitude of the end moment exceeds its prescribed yield value. The elastic end-moments are, from Eq. (2), given by M Bi;j;k ¼ ð4EI Bi;j =L;j Þ ~ hi;j;k , where ~ h is defined B B 1 ~i;j;k  ðh þ h by h Þ. The yield criterion therefore i;j;k 2 i;j;1
k amounts to a constraint on the absolute value of ~h. This is achieved by introducing the plastic end-rotation hP, which represents deformation of the beam–column connection, shown in Eq. (3) and Fig. 2: hBi;j;k ¼ hi;jþk
1
hPi;j;k .

ð3Þ

The plastic end-rotation remains constant under elastic loading/unloading; on the other hand, during yielding, changes in the effective rotation are constrained by the moment–rotation relationship (e.g. in the elastoplastic model hP is constrained to vary in such a way that ~ h remains constant). The response to a given set of incremental column end-rotations is calculated in terms of ~ hi;j;k , following which the plastic deformations are recovered by applying the inverse transformation:   4 ~ 1~ P hi;j;k
hi;j;1
k . hi;j;k ¼ hi;jþk
1
ð4Þ 3 2 The two end-rotations of any given beam exhibit a coupled response. This is realised computationally through local

θ

θP θB

Fig. 2. Beam rotations.

321

(i.e. at the level of a beam) subdivision of increments into linear sub-increments, so that the overall response of a beam is mathematically consistent with the prescribed moment–rotation relationships of its two ends, that is, the response to a particular sequence of increments is identical to that for any equivalent sequence, regardless of step size. In the applications reported here, end moments remain constant during yielding and there is no degradation of elastic tangent moduli; thus, the relationship between a beam end-moment and its associated column end-rotation is elastoplastic. For an individual connection this representation is not particularly accurate, but: (i) it is the simplest of the inelastic constitutive models, and is adequate for a prototypical demonstration; (ii) the overall behaviour of a storey with multiple columns will nevertheless be multilinear and will include Bauschinger effects; (iii) other piecewise-linear relationships are in any case easily incorporated into the model and (iv) it allows comparison with pushover analysis—where it is common practice to place pinned connections at the location of plastic hinges. Calculations involving ultimate failure of connections, wherein an endmoment becomes permanently zero on attainment of its prescribed ultimate effective rotation, have also been performed (though not reported here) and present no apparent difficulty to the numerical integration scheme. 3.3. Statics The contributions of all columns (represented by Eq. (1)) and beams (represented by Eq. (2)) associated with a given storey combine to give a storey stiffness matrix, written in block form as " # K0;1 K0;0 i i Ki ¼ ; ð5Þ B K1;0 K1;1 i i þ Ki where the displacement vector is T ½ ui
1 hi
1;0    hi
1;n ui hi;0    hi;n  , and the upper left block for example, which gives, for the columns of storey i, the end-forces at level i
1 corresponding to displacements at that level, is formed by combining the upper left quadrants (i.e. elements (1, 1), (1, 2), (2, 1) and (2, 2)) of the individual column stiffness matrices of Eq. (1), thus: 3 2 n P C C C K K       K i;0ð1;2Þ i;nð1;2Þ 7 6 j¼0 i;jð1;1Þ 7 6 7 6 C C 7 6 Ki;0 K 0    0 i;0ð2;2Þ ð2;1Þ 7 6 7 6 0;0 ð6Þ Ki ¼ 6 .. 7. .. .. .. 7 6 . . . 0 . 7 6 7 6 .. .. .. .. 6 . 0 7 . . . 5 4 KCi;nð2;1Þ 0    0 KCi;nð2;2Þ The upper right block, Ki0;1 , has exactly the same form: it gives the end-forces at level i
1 corresponding to displacements at level i, and is a combination of the upper right quadrants of the column stiffness matrices. The lower

322

S.M. Wilkinson, R.A. Hiley / Computers and Structures 84 (2006) 318–329

blocks, K1;0 and K1;1 i i , similarly, give end-forces at level i, and are formed from the lower left and right quadrants, respectively, of the column stiffness matrices. The beam reactions are included by combining the individual beam stiffness matrices of Eq. (2) into a block that augments K1;1 i . In the initial elastic case: 2

0

0

0

 .. .

 .. .

0 .. .

3

6 7 6 0 KB 7 KBi;1ð1;2Þ i;1ð1;1Þ 6 7 6 7 6 7 . . . B B .. .. .. 7 6 0 KB K þ K 6 7 i;1ð2;1Þ i;1ð2;2Þ i;2ð1;1Þ 6 7 KBi ¼ 6 . 7. . . 6 .. B B .. . Ki;2ð2;1Þ . Ki;n
1ð1;2Þ 0 7 6 7 6 7 6. 7 . . . 6. 7 B B B . . . . . . Ki;n
1ð2;2Þ þ Ki;nð1;1Þ Ki;nð1;2Þ 7 6. 4 5 0   0 KBi;nð2;1Þ KBi;nð2;2Þ

ð7Þ To accommodate yielding: elements of the member stiffness matrices are scaled appropriately to reflect the variation of hB with h; and additional force terms are generated to compensate for plastic deformations (see Eq. (3)). The structure stiffness matrix is now constructed from the storey stiffness matrices: neglecting the ground-level forces and displacements, which are determined by boundary conditions, we have 2

3 0;0 B K0;1 0   0 K1;1 1 þ K1 þ K2 2 6 7 .. .. .. .. 6 7 . . . . K1;0 6 7 2 6 7 6 7 .. .. 0;1 K0 ¼ 6 7. 0 . . K 0 6 7 m
1 6 7 .. 6 7 0;1 1;1 B 0;0 0;1 4 5 . K2 Km
1 þ Km
1 þ Km Km 1;0 1;1 B 0   0 Km Km þ Km

ð8Þ With the stiffness matrix in this form, the corresponding vector of displacements is D 0 ¼ ½ u1

h1;0

. . . h1;n

. . . . . . um

hm;0

. . . hm;n T . ð9Þ

It is useful in what follows to introduce the superscript notation whereby u denotes the set of elements or rows or columns corresponding to the positions of the ui in Eq. (9), and h denotes the remainder, i.e. the set corresponding to the hi,j. 3.4. Dynamics The lumped-mass idealization is used, together with the constraint (reasonable for several types of floor system [16, Section 9.2.4] that each floor diaphragm is rigid in its own plane but flexible in bending. The resulting mass matrix lumps the mass, mi, of storey i onto the diagonal element corresponding to the deflection DOF for that storey: using notation similar to that for the storey stiffness matrix in Eq. (5):

 Mi ¼



0

0

0

M1;1 i

;

ð10Þ

1;1 are zero. where M1;1 ið1;1Þ ¼ mi and all other elements of Mi The mass matrix for the structure, M0, is formed from these sub-matrices in the same way as the structure stiffness matrix in Eq. (8). The system of equations of motion of the structure, omitting damping forces, may now be written:

€ 0 ðtÞ þ K0 D0 ðtÞ ¼ P0 ðtÞ; M0 D

ð11Þ

where P0(t) is the vector of time-dependent external forces. As we are concerned here only with the dynamics induced by horizontal ground motion, the effective earthquake force P0(t) is given by € g ðtÞ. P0 ðtÞ ¼
M0 D

ð12Þ

Dg(t) is the rigid-body displacement vector corresponding to the ground motion, ug(t): it is given by iug(t), where i is the influence vector, defined (using the superscript notation introduced above) by iu = 1 and ih = 0. The fact that all but m of the elements of M0 and of P0 are zero, allows the set of equations of motion to be conveniently separated into two coupled sets: one dynamic, but smaller and hence more efficiently solved than the full set; the other static, expressing rotational equilibrium (i.e. static condensation, such as that used by Chopra [16, Section 9.3] can be employed). First, Eq. (11) is rewritten as u

uu u uh h u € Muu 0 D0 þ K0 D0 þ K0 D0 ¼ P0 ;

ð13Þ

u Khu 0 D0

ð14Þ

þ

h Khh 0 D0

¼ 0.

Eq. (14) gives the rotations as a linear combination of the translations, so that the former may be eliminated from Eq. (13) to give u

 u u € Muu 0 D0 þ K0 D0 ¼ P0 ;

ð15Þ

where the condensed stiffness matrix is  hh 
1 hu uh K0 ¼ Kuu K0 . 0
K0 K0

ð16Þ

3.5. Damping Modal damping has been included in the model and this was set at 5% for all modes as is usual practice for high-rise buildings. For simplicity, the damping matrix was formulated using the properties of the unyielded structure and was kept constant during the analysis. Free vibration analysis, for the initial elastic system, entails computation of the characteristic values and vectors of Eq. (15) with no external forcing. Time history analysis was achieved by solving Eq. (15) using a fourth-order Runge–Kutta numerical integration scheme. Accuracy was controlled by comparing displacements and velocities at each time step with those calculated using two half time steps.

S.M. Wilkinson, R.A. Hiley / Computers and Structures 84 (2006) 318–329

4. Verification examples

P3

323

I

M

4.1. Example 1 I

To validate the proposed model, two examples are presented. In the first example the validity of the simplifying assumptions are investigated. This is achieved by comparing the results of an analysis of the three-storey plane frame, shown in Fig. 3, obtained with a commercially available structural analysis package, to the results obtained using the proposed model. All structural elements are beam elements and their section and material properties are shown in Fig. 3. The frame was subjected to successive elastic static analyses (push-over analysis) as would be typical for a seismic resistant design. After each analysis the locations of plastic hinges and the associated loads required to develop them were calculated. The base shear was distributed up the building in proportion with the first mode shape of the unyielded structure as per common practice. The collapse sequence and final collapse mechanism are displayed in Fig. 4 together with the natural period of the unyielded structure. The natural periods calculated by the proposed model are in close agreement with those calculated by the finite element analysis and demonstrate that the model can predict accurately the free vibration response of the building.

Period T1 T2 T3

P2

I

I

P1

h

M

I

I

h

I

h

M

I

L Fig. 3. Example 1—I = 5.208 · 10
3 m4, M = 62,500 kg, h = 3 m, L = 10 m, E = 200 · 106 kN/m2, P1 = 134 kN, P2 = 350 kN, P3 = 516 kN, yield moment = 1000 kNm.

In this paper, the definition of collapse is considered to be when sufficient plastic hinges exist to form a mechanism. Since the material behaviour of the beams is assumed to be

FEM Simplified Model 0.375 s 0.0860 s 0.0367 s

I

0.374 s 0.0860 s 0.0367 s

Fig. 4. Collapse mechanisms.

324

S.M. Wilkinson, R.A. Hiley / Computers and Structures 84 (2006) 318–329

perfectly elastic perfectly plastic and gravity forces (P–d effects) are neglected, the model structure will always gain static equilibrium after cessation of the earthquake loading. Collapse, for the model, can occur only if ultimate failure of connection is allowed or if gravity is included, in which case the external loading remains non-zero as long as there is lateral deformation of the structure.

In addition to the push-over analysis by the commercially available structural analysis package (where hinges were placed at the point of plastic moment) a response history analysis was performed using the simplified model. To compare the two analyses, it is necessary to assume perfectly elastic-perfectly plastic material behaviour. While this is not the most appropriate relationship for the

900 800

Base Shear (kN)

700 600 500 storey 1 push-over storey 1 time-history storey 2 push-over storey 2 time-history storey 3 push-over storey 3 time-history

400 300 200 100 0 0

0.01

0.02

0.03

0.04

0.05

Displacements (mm)

(a) Base Shear Vs Deflection 900 800

Base Shear (kN)

700 600 storey 1 push-over

500

storey 2 push-over

400

storey 3 push-over

300

storey 1 time-history

200

storey 2 time-history storey 3 time-history

100 0 -0.007

-0.006

-0.005

-0.004

-0.003

-0.002

-0.001

0

Rotations (rads)

(b) Base Shear Vs Rotation 900 800

Base Shear (kN)

700 600 500 400 storey 1 push-over storey 2 push-over storey 3 push-over storey 1 time-history storey 2 time-history storey 3 time-history

300 200 100 0 0

200

400

600

800

1000

Moment (kN-m)

(c) Base Shear Vs Beam Moments Fig. 5. Comparison of push-over analysis with time–history analysis.

1200

S.M. Wilkinson, R.A. Hiley / Computers and Structures 84 (2006) 318–329

analysis of structures subjected to earthquake loading it is convenient for verifying the accuracy of the model. It is also convenient to use this relationship to study the appropriateness of push-over analysis when higher order modes of vibration become significant in the overall response of the structure. For more realistic non-linear behaviour, models such as Takeda et al. [17] could be used for reinforced concrete, or Krawinkler et al. [18] for steel frames. The loading chosen for this analysis was a cosine-wave displaced by half its amplitude. To simulate static loading a cosine-wave loading with a period that was very long (10 s) compared to the natural period of the structure (0.375 s) was chosen. The shape of the wave was chosen also to minimise the dynamic effects. The derivative of the cosine-wave (which is proportional to velocity in this analysis) at the start of the analysis is zero. The comparison of push-over analysis with time–history analysis is shown in Fig. 5. Fig. 5a shows the base shear versus deflection, Fig. 5b shows the base shear versus rotation and Fig. 5c shows the base shear versus the beam moments. The response history analysis was terminated when the plastic hinge in the top storey was reached (this represents collapse of the structure and occurred at V = 846.8 kN compared to 839.6 kN for the push-over analysis). From Fig. 5, it can be seen that the two analyses are in very good agreement. The collapse loads calculated by the two analyses are within 1% of each other. Furthermore the time–history of the second analysis follows the static analysis very closely. The only parts of the graphs that do not match very closely are those located directly after yielding (indicated with markers). The reason for this is that because of the perfectly plastic behaviour assumed for the analysis, after yielding there is a sudden loss of stiffness. This loss in stiffness results in the frame accelerating in the direction of loading which in turn produces inertia

1

1

forces. These forces result in the small differences between the graphs. 4.2. Example 2 In the first example the two columns were identical and therefore had identical rotations. If the columns have different stiffnesses or the frame has more than one bay, then the rotations at the top of each column will not be equal. To demonstrate the methodÕs ability to deal with buildings that have columns whose rotations differ from each other, the 2 storey frame with 2 bays, shown in Fig. 6, is analysed. All structural elements are beam elements and their section and material properties are the same as for the first example. Again a push-over analysis was performed using a commercially available structural analysis package and the results compared to a time–history analysis using the simplified model. The load was similar to that used in Example 1. The collapse sequence and final collapse mechanism are displayed in Fig. 7 together with the loads

1

2 2

V= 840 kN

I

P2

M

I

I

I I

P1

M

I

L

Fig. 6. 2 Storey 2 bay example.

1

2 2

1

V= 640 kN 3

3

4 4

3

1

1

2 2

1

V = 856 kN

Fig. 7. Collapse sequence and collapse loads.

I

h

I

h

I

I

L

V = 604 kN 3

325

326

S.M. Wilkinson, R.A. Hiley / Computers and Structures 84 (2006) 318–329

required to form the plastic hinges. In Fig. 8a and b, the time–history analysis and push-over analysis are compared. Fig. 8a plots the storey deflections against the base shears, while Fig. 8b plots the rotation of each joint against base shear. The figures show excellent agreement between the two analyses and demonstrate the modelÕs ability to extract rotation information for each column in a storey. 4.3. Example 3 Finally, to demonstrate differences between a push-over analysis and materially non-linear time–history analysis, the 15 storey example shown in Fig. 9 was analysed. Again the properties of the structure are the same as in the previous two examples, except that in this example the model was subjected to a real earthquake load and the yield moment was reduced to a value that allowed the formation of a mechanism (which represents collapse as defined previously). The earthquake load chosen was El Centro

1940 N–S. For the push-over analysis the base shear of V = 1969 kN was calculated by multiplying the relevant pseudo acceleration from the response spectrum, by the total mass of the structure. Again this base shear was distributed up the building in proportion to the first mode shape. The first and second mode shapes are also given in Fig. 9, along with the first three natural periods of the unyielded structure. In Fig. 10, the collapse sequence and collapse load is presented for each type of analysis. The possibility of premature local failure of individual elements before the collapse mechanism developed was excluded by monitoring the rotations at the ends of the beams as well as inter-storey drifts. The maximum value of inter-storey drift for the example was 0.86%. This is greater than the serviceability allowance of Eurocode 8 [2] (0.6% for structures where non-structural elements are unlikely to be damaged by excessive deflections) but less than the IBC [3] ultimate limit (1.5% for structural use II). The inter-storey drift ratios themselves will not result in collapse because the joint rotations remained below

900 800

Base Shear (kN)

700 600 500 Storey 1 push-over 400

Storey 1 time-history Storey 2 push-over

300

Storey 2 time-history 200 100 0 0.00E+00

2.00E-03

4.00E-03

6.00E-03

8.00E-03

1.00E-02

1.20E-02

1.40E-02

Deflection (mm)

(a) Base Shear vs Deflection 900 800 700 600 500 400 300 200 100 0 -0.0025

-0.002

-0.0015

-0.001

-0.0005

Rotation (rads.)

(b) Base Shear vs Rotation Fig. 8. Deflections and rotations—Example 2.

0

Base Shear (kN)

joint 1 time-history joint 1 pushover joint 2 time history joint 2 pushover joint 3 time history joint 3 pushover joint 4 time history joint 4 pushover

S.M. Wilkinson, R.A. Hiley / Computers and Structures 84 (2006) 318–329

15 Storeys @h = 15h

Is

M

Ic

L V

327

1st mode 2nd mode 1.000 -1.000 0.989 -0.901 0.971 -0.741 0.944 -0.523 0.909 -0.264 0.866 0.017 0.815 0.297 0.757 0.553 0.691 0.764 0.620 0.912 0.542 0.985 0.459 0.975 0.371 0.882 0.275 0.705 0.163 0.435 Mode shapes Natural Periods T1 = 1.49s, T2 = 0.48s, T3 = 0.28s.

Fig. 9. 15 Storey example.

15 14 13 12 11 10 9 8 7 6 5 4 3 2 1

15 14 13 12 11 10 9 8 7 6 5 4 3 2 1

14 12 7 6 3 2 1 9 10 10 15 13 7 5 4

12 8 8 8 7 8 14 14 12 6 5 4 3 2 1

14 12 7 6 3 2 1 9 10 10 15 13 7 5 4

Mp = 1933 kN-m

Mp = 56 kN-m

(a) Pushover analysis

(b) RHA

12 8 8 8 7 8 14 14 12 6 5 4 3 2 1

Mp = 185 kN-m (c) Geometrically non-linear RHA

Fig. 10. Collapse sequence.

0.01 radians—and this figure should be easily achieved with a good earthquake resistant connection. As can be seen the collapse load and collapse sequences for the two analyses are different. The push-over analysis predicts a plastic moment capacity of 1913 kN-m to prevent collapse with a collapse sequence where plastic hinges initiate at the bottom of the building and propagate to the top (i.e. the collapse sequence is storey 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15). For the response history analysis the required plastic moment capacity to prevent collapse is 56 kN-m. If the yield capacity of the beam is less than this value, collapse will occur with plastic hinges forming in the following sequence: namely storey 9, 10, 11, 1, 2, 12, (3, 13), 8, (6, 7), 14, 4, 15, 5. In the response history analysis, the col-

lapse sequence starts in the middle of the structure. The reason for this is that the second mode of vibration is significant in this example. A modal analysis was performed and the results showed that the first mode of vibration was responsible for 86% of the modal mass while the second mode of vibration was responsible for 9%. In Fig. 11, the response spectrum of El Centro 1940 is presented. In Fig. 11, it can be seen that the acceleration for the first mode is 2.1 m/s/s while the acceleration for the second mode is 9.9 m/s/s resulting in base shears of 1693 kN and 835 kN, respectively. This relatively high contribution of the second mode to the response results in plastic hinges forming higher up the building than the push-over analysis predicted.

328

S.M. Wilkinson, R.A. Hiley / Computers and Structures 84 (2006) 318–329 14 12

A (m/s/s)

10 8 6 4 2 T1

T2

0 0

0.5

1

1.5

2

2.5

Tn (sec.)

Fig. 11. Response spectrum of El Centro.

The required plastic moment to prevent collapse calculated using response history analysis is significantly lower than that predicted by the push-over analysis. The reasons for this are, firstly that the base shear for the push-over analysis was based on an elastic response, while the response history analysis assumed a perfectly elastic-perfectly plastic moment rotation relationship. Since thirty plastic hinges are required to create a collapse mechanism the model can absorb a large amount of energy through hysteretic damping. For the equivalent static analysis the base shear would be divided by the response modification factor. For a steel moment resisting frame the value would 8. As discussed previously response modification factors are empirical and as such are not particularly relevant to a theoretical example such as this. Of more relevance is the displacement ductility factor—which is a measure of the ratio of ultimate displacement to yield displacement. Investigation of the displacement histories of the first storey of this example revealed a displacement ductility factor of 45. Chopra [16, Section 7.10], constructs constant ductility response spectra for elasto-plastic systems with various specified values of ductility factor. While there is no precise relationship, in general, for buildings with natural periods over 0.6 s (as is the case with this example) the required yield strength (strength demand) is roughly proportional to the inverse of the ductility factor. The very large ductility factor results from choosing a low yield moment to create a collapse mechanism. In reality stiffness degradation and P–d effects would result in collapse at a lower load. This analysis has been performed using the matrix of Eq. (1). For this example, the vertical load was calculated by multiplying the storey masses by the acceleration due to gravity. For simplicity and so the analysis is compatible with the geometrically linear example, the vertical loads have only been used for the purpose of P–d effects and are distributed to each column based on tributary area. The results for geometrically non-linear analysis are: the required plastic moment capacity to prevent collapse is 185 kNm. If the section capacity is less than this value collapse will occur with plastic hinges forming in the following sequence: namely storey 1, 2, 3, 4, 5, 6, 11, (10, 12, 13, 14), (7, 15), (8, 9)

as shown in Fig. 10c. Notwithstanding the simplifications used in the model, the results of Example 3, show that push-over analysis may result in an incorrect sequence of plastic hinge formation and consequently the incorrect collapse load. 5. Conclusion An efficient simplified model for the analysis of high-rise buildings has been presented. This model can analyse nominally symmetric structures using only v + 1 degrees of freedom per floor (where v is the number of vertical elements in a floor). Furthermore prior to solution the stiffness matrix is partitioned and condensed into two smaller matrices. The rank of the condensed stiffness matrix is only m. Its construction, which requires the inversion of the rotational, rank mv, stiffness matrix, is required only at time-steps where the pattern of yielding has altered from the previous time-step. This makes it particularly attractive for non-linear response history analysis. The model presented incorporates perfectly elastic perfectly plastic moment rotation relationship for the beam to column connections, but any suitable relationship could be incorporated. Three verification examples have been presented. The model was shown to be capable of analysing simple structures to within less than 1% of the results obtained by finite element analysis. The model accurately predicts higher modes of vibration and therefore can be used to consider the influence of these on the collapse of buildings. Using a simple push over analysis for structures with significant 2nd modal masses and accelerations may produce the wrong collapse sequence and therefore the incorrect collapse load. To accurately determine the collapse loads of structures, P–d effects need to be considered as does the true moment rotation relationship of the connections especially if stiffness degrades over successive cycles. Acknowledgement The research presented in this paper was funded by the Engineering and Physical Sciences Research Council of the

S.M. Wilkinson, R.A. Hiley / Computers and Structures 84 (2006) 318–329

UK (EPSRC) under Grant GR/R40937/01. This financial support is gratefully acknowledged. References [1] Kim SD, Hong WK, Ju YK. A modified dynamic inelastic analysis of tall buildings considering changes of dynamic characteristics. Struct Des Tall Buildings 1999;8:57–73. [2] Eurocode No 8, EC8:1994, Common Unified Rules for Structures in Seismic Regions, Commission of European Communities; 1994. [3] International Code Council, International Building Code 2000, Joint UBC, BOCA, SCCI, Whittier, California; 2000. [4] SANZ—Standards Association of New Zealand, New Zealand Standard Code of Practice for General Structural Design and Design Loadings for Buildings, 1992, NZS 4203; 1992. [5] Park R. Capacity design of ductile RC building structures for earthquake resistance. Struct Eng 1992;70(16/18):279–89. [6] Uang CM. Establishing R and Cd factors for building seismic provisions. J Struct Eng 1991;117(1):19–28. [7] Whittaker A, Hart G, Rojahn C. Seismic response modification factors. J Struct Eng 1999;125(4):438–43. [8] Krawinkler H, Seneviratna G. Pros and cons of a push-over analysis of seismic performance evaluation. Eng Struct 1998;20(4–6):452–64.

329

[9] Berg GV. Elements of structural dynamics. Englewood Cliffs, New Jersey: Prentice Hall; 1989. [10] Chintanapakadee, Chopra. Evaluation of modal pushover analysis using generic frames. Earthquake Eng Struct Dyn 2003;32:417–42. [11] Thambiratnam DP, Thevendran V. Simplified analysis of asymmetric buildings subjected to lateral loads. Comput Struct 1992;34(5): 873–80. [12] Wilkinson S, Thambiratnam D. Development and application of a simplified procedure for seismic analysis of asymmetric buildings. Comput Struct 2001;79(32):2833–45. [13] Wilkinson SM, Hiley R. Use of advanced analysis techniques to determine the collapse loads of buildings. In: Proc ASRANET Int Colloquium, 8–10 July, 2002, Glasgow. [14] Krenk S, Vissing-Jørgensen C, Thesbjerg L. Efficient collapse analysis techniques for framed structures. Comput Struct 1999;72:481–96. [15] Krenk S. Mechanics and analysis of beams columns and cables. 2nd ed. Denmark: Springer-Verlag; 2001. Section 3.4. [16] Chopra AK. Dynamics of structures. Englewood Cliffs, New Jersey: Prentice Hall; 1995. [17] Takeda T, Sozen M, Nielsen N. Reinforced concrete response to simulated earthquakes. J Struct Eng (ASCE) 1970:2257–573. [18] Krawinkler H, Bertero V, Popov E. Inelastic behaviour of steel beam to column subassembleges. Report No EERC 71-7, University of California, Berkley; 1971.