A methodology for deriving analytical fragility curves for masonry buildings based on stochastic nonlinear analyses

Engineering Structures 32 (2010) 1312–1323

Contents lists available at ScienceDirect

Engineering Structures journal homepage: www.elsevier.com/locate/engstruct

A methodology for deriving analytical fragility curves for masonry buildings based on stochastic nonlinear analyses M. Rota a,b,∗ , A. Penna b , G. Magenes a,b a

University of Pavia, Department of Structural Mechanics, via Ferrata 1, Pavia, Italy

b

European Centre for Training and Research in Earthquake Engineering, via Ferrata 1, Pavia, Italy

article

info

Article history: Received 31 July 2009 Received in revised form 11 December 2009 Accepted 13 January 2010 Available online 1 February 2010 Keywords: Fragility curves Seismic vulnerability Masonry buildings Nonlinear analysis Monte Carlo simulation

abstract A new analytical approach for the derivation of fragility curves for masonry buildings is proposed. The methodology is based on nonlinear stochastic analyses of building prototypes. Since such structures are assumed to be representative of wider typologies, the mechanical properties of the prototypes are considered as random variables, assumed to vary within appropriate ranges of values. Monte Carlo simulations are then used to generate input variables from the probability density functions of mechanical parameters. The model is defined and nonlinear analyses are performed. In particular, nonlinear static (pushover) analyses are used to define the probability distributions of each damage state whilst nonlinear dynamic analyses allow to determine the probability density function of the displacement demand corresponding to different levels of ground motion. Convolution of the complementary cumulative distribution of demand and the probability density function of each damage state allows to derive fragility curves. © 2010 Elsevier Ltd. All rights reserved.

1. Introduction The evaluation of seismic vulnerability of existing buildings has become really relevant in the last decades due to the frequent occurrence of earthquakes, which have demonstrated that the number of victims and the amount of economic losses depend significantly on the seismic behaviour of structures. Masonry structures are a much diffused type of construction which can be built rapidly, cheaply and often without any plan or particular technical competence. Indeed, this is one of the most common housing types built across the world and will continue to be so in the foreseeable future. Masonry also represents the structural type of a large architectural heritage that needs to be preserved. Moreover, old unreinforced masonry buildings constitute the large majority of most urban aggregates in several seismicprone countries. The main characteristics of masonry buildings are high rigidity, low tensile and shear strength, low ductility and low capacity of bearing reverse loading. These are the main reasons for the frequent collapse of masonry buildings during earthquakes, often responsible for a considerable number of casualties.

∗ Corresponding author at: University of Pavia, Department of Structural Mechanics, via Ferrata 1, Pavia, Italy. Tel.: +39 0382 516950; fax: +39 0382 529131. E-mail addresses: [email protected] (M. Rota), [email protected] (A. Penna), [email protected] (G. Magenes). 0141-0296/$ – see front matter © 2010 Elsevier Ltd. All rights reserved. doi:10.1016/j.engstruct.2010.01.009

It should be clear hence that earthquake risk reduction should mainly concentrate on masonry structures. Nevertheless, the efforts put into the study and improvement of this type of construction are very limited since masonry buildings tend to be regarded as out-of-date, even though they are still diffusely built even in seismic areas and, in any case, they still make up today a very large proportion of the world’s existing building stock. In particular, the knowledge of vulnerability functions for what concerns masonry buildings is quite limited. Recently, some studies on large scale vulnerability of masonry buildings [1–4] have been made, but the available research is still far behind that concerning reinforced concrete buildings. However, the reactions and mechanisms of failure of rigid unframed buildings to a seismic event are quite different from those of flexible framed buildings and hence specific vulnerability studies are required. Although linear analyses are not adequate for masonry buildings (e.g. [5]), the lack of reliable models for calculating the inelastic response to input accelerograms of load-bearing masonry structures has limited, at least in the past, the execution of nonlinear analyses. Most of the few recent studies available in the literature concerning the derivation of analytical fragility curves for masonry buildings (e.g. [6–12]) are based on more or less approximate structural models and/or simplified analysis methods. This is also due to the inherent complexity and relative unfamiliarity with nonlinear dynamic analyses, which are more computationally demanding than any other procedure. In this study, a new advanced methodology for obtaining analytical fragility curves for classes of buildings is proposed.

M. Rota et al. / Engineering Structures 32 (2010) 1312–1323

This approach allows to take rigorously into account the different sources of uncertainty involved in the problem, by derivation of the probability distributions of both capacity and demand through 3D nonlinear analyses of entire structures. To illustrate the methodology, nonlinear structural analyses have been performed on a prototype building representative of a structural typology that is common in Italy. Fragility curves for such typology have been derived. The approach takes advantage of the capabilities of the program TREMURI, a frame-type macro-element global analysis program developed at the University of Genoa [13] which is able to perform nonlinear pushover and time history analyses on masonry buildings [14]. Since in most real cases building properties are not well known and since, on the other hand, it is assumed that the prototype building is representative of a class of buildings with similar structural characteristics, mechanical parameters have been considered as random variables. Hence Monte Carlo simulations have been carried out based on realistic ranges of the parameters derived from experimental tests and other sources of information. First of all, the assumptions adopted for identifying mechanical damage states are presented. Four damage levels have been considered, two of which are defined from the single pier capacity curve while the other two are identified on the global pushover curve of the building. The prototype building is hence described and the assumptions used for defining the model explained. After that, the results of the different sets of analyses performed are presented. In particular nonlinear static pushover analyses have been carried out with the aim of defining the probability distributions associated to each damage state. On the other hand, time history analyses were necessary to define the probability density function (pdf) of the displacement demand imposed on the building by different levels of ground motion, taking into account the significant variability of natural records. Convolution of the complementary cumulative distribution of demand and the probability density functions of each damage state allows to derive fragility points which will then be fitted by lognormal distributions in order to obtain analytical fragility curves. 2. Mechanical definition of damage states In order to derive analytical fragility curves it is necessary to define damage states in terms of a mechanical parameter that can be directly obtained from the analysis. Three limit states are usually defined for masonry buildings (e.g. [15,16]):

• Elastic (cracking) limit, where the wall displacement exceeds the elastic limit and the first significant crack forms in the wall, which changes the initial stiffness; • Maximum resistance, determined by the shear and displacement corresponding to the attainment of the building maximum resistance; • Ultimate state, where the building resistance deteriorates below an acceptable limit, e.g. 80% of the maximum resistance [17]. Notice that a real masonry building usually does not collapse at the so-defined ultimate state. However, although the actual collapse takes place at a lateral displacement larger than that corresponding to the ultimate state, at this latter displacement the structure is already damaged beyond repair [16]. In this work, four mechanical damage states have been considered. Two of them can be identified from the response of a single masonry pier while the other two are found from the global response of the building.

1313

Fig. 1. Identification of the yield, cracking and ultimate drifts on the pushover curve (envelope of force–displacement curve from cyclic shear testing) of a single pier and its bilinear approximation.

The first damage state, DS1, is identified by the attainment of the yield displacement δY of the bilinear approximation to the capacity curve of a single masonry pier, as indicated in Fig. 1. The bilinear approximation is obtained by fixing the initial stiffness as being secant to the point corresponding to 70% of the maximum resistance and the equivalent resistance so that the area below the bilinear curve up to the ultimate limit state coincides with the area below the capacity curve up to the same point [17]. It should be noticed that δY is always larger than the drift corresponding to 70% of the shear resistance and hence it is likely that some cracks have already occurred at such deformation value. The second damage state, DS2, is identified by the drift corresponding to the first shear cracking of the pier, δS , and is also shown in Fig. 1. This value can be derived from experimental test results, whenever available. In this work, values of δS have been obtained directly from the experimental test report, where they were explicitly indicated. Sometimes they are not explicitly reported since the first cracking is not easy to detect during experimental tests, especially if the wall surface is not painted or covered by plaster. For this reason the experimental values associated to δS are typically larger than those associated to δY . The damage levels DS3 and DS4 have been derived from global pushover curves of the building. In particular, DS3 is assumed to correspond to the attainment of the maximum shear resistance, while DS4 corresponds to the attainment of 80% of that value. Notice that a damage level corresponding to complete collapse of the structure cannot be analytically defined since the total collapse condition cannot be numerically represented. Moreover it would occur in correspondence of very high values of drift, at which the structure or some parts of it loose their vertical loadbearing capacity. The considered damage states are represented in Fig. 2, where they are identified on the global pushover curve of the structure. 3. Stochastic Monte Carlo simulations The proposed methodology for deriving analytical fragility curves requires analyses of a population of buildings subjected to different levels of ground motion. For each building typology, a prototype is identified and the population of buildings belonging to that typology is generated through the treatment of selected structural parameters as random variables. This allows to take into account both the fact that in most real cases building properties are not well known and also the variability intrinsic in the definition of building typologies grouping together buildings with different characteristics and mechanical properties.

1314

M. Rota et al. / Engineering Structures 32 (2010) 1312–1323

Fig. 2. Identification of the considered damage states on a building pushover curve.

A probability density function (typically gaussian) is then assigned to each random variable based on realistic ranges of variation derived from experimental tests and other sources of information. Values of the parameters are then extracted from the distributions using appropriate sampling techniques. Finally, the sampled values are combined to define a series of structures with different characteristics, all nominally representing the same building. Clearly, consideration of all possible uncertainties in the global and local structural characteristics and ground motion yields an extremely large number of permutations for the analysis. This is neither practical nor justified, considering the general nature of the application of the final vulnerability curves. In the proposed approach, the geometry of the building is treated deterministically, meaning by this the plan and elevation configuration of the structure, the wall size and position and so on. On the other hand, all the mechanical properties of masonry necessary for the definition of the building model have been considered as random variables. This is justified by the assumption that, as stated also in [18] concerning reinforced concrete structures, material variability can be regarded as the most important source of uncertainty in the determination of structural response, with all other sources either deriving directly from its effects or being insignificant in size compared to it. It is noteworthy that this type of uncertainty is also considered in the design process, where characteristic values for material strengths are used and are further reduced using partial safety factors to ensure structural safety. The Monte Carlo algorithm, which is a very diffused and well known method, has been used for sampling the uncertain parameters from their associated probability density functions (e.g. [19,20] and references therein). Monte Carlo analyses have been performed using a software called STAC, developed at the CIMNE

(International Centre for Numerical Methods in Engineering) of Barcelona [21]. This software allows to perform stochastic simulations using any type of analysis tool. The user is required to define appropriate probability distributions of the input variables of interest. The software then calls the structural program which runs the analysis and returns the statistics of the output variables. A number of Monte Carlo simulations is performed until convergence of both the input and output variables to their mean is reached and the value of standard deviation becomes stable. Two sets of Monte Carlo analyses have been carried out. In the first set, values of the mechanical properties have been varied within the selected ranges and randomly assigned uniformly over the entire building model. In the second set, an additional randomness was introduced with the creation of a library of randomly defined materials which are then randomly assigned to the different structural elements of the model. 4. Selected building case study 4.1. Building description The selected building prototype consists in a three-storey masonry building located in Benevento (southern Italy), in the neighbourhood called Rione Libertà. This building is typical of the construction typologies of the entire Rione Libertà. It has been constructed in 1952 and has plane dimensions of 17.7 × 14.3 m and a height of 11.05 m. The masonry bearing structure is entirely realised with tuff units, while floors are made of reinforced concrete. Reinforced concrete tie beams guarantee the connection between floors and masonry walls. A picture of the selected building is reported in the left part of Fig. 3. The geometry of the structure has been considered deterministic. The building shows plan regularity and symmetry around the axis of the stairs. The typical plan of the building (which is practically identical at all storeys) is reported in the right part of Fig. 3. For what concerns the material characteristics, since this building is a prototype considered to be representative of a class of similar buildings, the characteristics of the materials have been assumed to be random variables, defined as already discussed in the previous section. 4.2. Building model and analysis tool For what concerns modelling and analysis of masonry structures, either very sophisticated finite element models or extremely simplified methods are commonly adopted. The analysis tool used

Fig. 3. Prototype building located in Benevento: A picture (left) and the plan of a storey (right).

M. Rota et al. / Engineering Structures 32 (2010) 1312–1323

1315

Table 1 Parameters derived from calibration based on experimental results.

T1-1 T1-2 T1-3 T1-4 T2-5 T2-6 Min Max

Fig. 4. 3D view of the building model.

in this work is based on the effective macro-element approach which allows an accurate modelling strategy of masonry buildings at a reasonable computational burden. This model is implemented in the program TREMURI. The software and the algorithms embedded in it are described in detail in several literature works (e.g. [14,22–24]). The program is based on the nonlinear macro-element model proposed by Gambarotta and Lagomarsino [25], further modified by Penna [14], representative of a whole masonry panel (pier or spandrel beam). The software allows to perform nonlinear seismic analyses of unreinforced masonry buildings using a set of analysis procedures: incremental static analysis (Newton–Raphson) both with force or displacement control, 3 pushover analysis and 3D time history dynamic analysis using the Newmark integration method and Rayleigh viscous damping. A view of the 3D model of the building is shown in Fig. 4. 3D modelling of the building is based on the identification within the construction of the seismically resistant structure, constituted by walls and floors. The walls are the bearing elements, while the floors, apart from sharing vertical loads to the walls, are considered as planar stiffening elements (flexible diaphragm modelled by orthotropic 3–4 nodes membrane elements) on which the horizontal action distribution between the walls depends. The model is based on the hypothesis that the seismic response of the building is governed by a global box-type behaviour, assuming that local mechanisms (mainly out-of-plane) are prevented by appropriate structural details and/or connecting devices (e.g. tie rods or tie beams). With this assumption, the local flexural behaviour of the floors and the local wall out-of-plane response are not computed because they are considered negligible with respect to the global building response, which is governed by their inplane behaviour. This is acceptable for this type of building, where stiff diaphragms and low height/thickness ratios for the walls render out-of-plane response a secondary phenomenon. Notice that a global seismic response is possible only if vertical and horizontal elements are properly connected. A frame-type representation of the in-plane behaviour of masonry walls is adopted: each wall of the building is subdivided into piers and lintels (2-nodes macro-elements) connected by rigid areas (panel nodes). 4.3. Identification of mechanical parameters’ ranges The first step for implementing a Monte Carlo simulation is identifying realistic ranges of variation for the various parameters needed for the numerical model. The building under study is made of tuff masonry. For this type of material, several experimental tests have been carried out and

G (MPa)

fv o (MPa)

µ (–)

G c (– )

β (–)

δ(%)

500 500 500 500 700 750 500 750

0.105 0.105 0.16 0.15 0.13 0.2 0.105 0.2

0.06 0.06 0.05 0.05 0.05 0.08 0.05 0.08

8 10 4 6 8 6 4 10

0.3 0.2 0.3 0.2 0.4 0.3 0.2 0.4

0.78 0.72 0.52 0.6 0.56 0.6 0.52 0.78

the results are available in the literature. Among others, the tests carried out by Faella et al. [26] have been selected. They consist in in-plane cyclic shear-compression tests carried out on specimens made of cement mortar and tuff units obtained from demolished buildings erected in Naples in the last 2 centuries. Tests have been performed on three types of masonry piers, all with 1250 × 1300 × 500 mm (h × l × t) dimensions, differing for the applied axial load and the constructive details:

• The first type (T1) is a double-wythe wall made of roughly dressed tuff blocks, without through stones connecting the external wythes and with incoherent filling material in between. • The second type (T2) consists of a monolithic wall made of tuff units similar to those of type one. • The third type (T3) consists of panels of the first type, injected with cement mixture and will not be considered in this work. The tests, carried out at the ISMES laboratory of Bergamo, have been performed in displacement control by imposing a cyclic horizontal displacement at the top of the wall and measuring the horizontal force. The vertical axial force was either equal to 130 kN (σ = 0.2 MPa) or 325 kN (σ = 0.5 MPa), while the imposed displacement has been incremented by 0.5 mm every three cycles. Each test has been stopped when the wall has shown evidence of critical damage. All the tests carried out on the first two types of walls have been simulated with the program TREMURI in order to identify meaningful ranges of variation for the different model parameters. With the aim of obtaining the desired parameter ranges, different cyclic pushover analyses were carried out, playing with the model parameters until a good fit of the experimental tests’ results was obtained. A plot of the cyclic force–displacement curves obtained form the experimental tests with the nonlinear model results superimposed is shown in Figs. 5 and 6 for the different tested specimens. Based on the comparison with experimental results, it has been possible to calibrate only the parameters related to shear failure modes, since all experimental specimens failed in shear. In particular, the following parameters have been obtained: shear modulus G, initial shear resistance for zero compression (cohesion) fv o , friction coefficient µ, nonlinear shear deformability ratio Gc , softening parameter β , ultimate shear drift δv . The results of the nonlinear identification are summarised in Table 1. The missing parameters, which cannot be directly obtained from cyclic tests, have been derived from the indications reported in Annex 11.D of the OPCM 3274 [17]. For tuff masonry with good quality mortar and lack of through stones the following ranges of values for the different parameters are suggested:

• Compressive strength of masonry: 1.2 < fm < 1.8 MPa • Young modulus in compression: 1350 < E < 1890 MPa • Specific weight: w = 16 kN/m3 . In case of transversal connections, the compressive strength range can be corrected in:

1316

M. Rota et al. / Engineering Structures 32 (2010) 1312–1323

Fig. 5. Comparison of experimental test results (grey thin curve) and numerical simulations (black thick curve) for the different tested specimens: T1–1, T1–2, T1–3 from top to bottom. Table 2 Parameters’ ranges adopted. E (MPa) G (MPa) fm (MPa) fv o (MPa) µ (–) Min 1350 Max 1890

500 750

1.2 2.7

0.105 0.2

0.05 0.08

Gc (–)

β (–) δv (%) δf (%)

4 10

0.2 0.4

0.52 0.78

– –

• Compressive strength of masonry: 1.8 < fm < 2.7 MPa. In addition, for what concerns the ultimate flexural drift δf , a value of 0.8% has been adopted as the mean value (as also suggested in [27]), with an assumed coefficient of variation of 10%. In conclusion, the ranges reported in Table 2 have been adopted. In the numerical model, normal distributions have been assumed for all the mechanical parameters, with a mean value corresponding to the central value of the interval and a standard deviation such that 95% of the data lie within the interval. The parameters of the normal distribution of each mechanical quantity needed for the numerical analysis are summarised in Table 3.

Fig. 6. Comparison of experimental test results (grey thin curve) and numerical simulations (black thick curve) for the different tested specimens: T1–4, T2–5, T2–6 from top to bottom.

In order to define the damage states, a probability distribution of values of the drift corresponding to first cracking is also required. The interpretation of the experimental tests provides two sets of drift values: the first (δs ) is obtained directly from the tests and is defined as corresponding to the first evidence of shear cracking; the second (δy ) is derived from the bilinear approximation of the envelope capacity curve and represents the drift associated to the elastic limit of the equivalent bilinear curve. The values of these two parameters for the different tests are reported in Table 4, while the parameters of the probability distribution adopted in the damage state definition are summarised in Table 5. 5. Results of pushover analyses A significant number of pushover analyses has been performed on the selected building. The analyses have been carried out only in the x direction, since the building appears to be more vulnerable in x than in y and the first vibration mode is along the x direction.

M. Rota et al. / Engineering Structures 32 (2010) 1312–1323

1317

Table 3 Parameters used in the stochastic analyses.

Mean Std. dev.

E (MPa)

G (MPa)

fm (MPa)

fv o (MPa)

µ (–)

Gc (–)

β (– )

δv (%)

δf (%)

1620 135

625 62.5

1.95 0.375

0.1525 0.02375

0.065 0.0075

7 1.5

0.3 0.05

0.65 0.065

0.8 0.08

Table 4 First cracking parameters derived from experimental tests.

T1-1 T1-2 T1-3 T1-4 T2-5 T2-6 Min Max

δs (%)

δy (%)

1.72 1.99 1.62 1.60 2.05 1.69 1.60 2.05

0.65 1 1.07 1.16 0.68 0.85 0.65 1.16

Table 5 Parameters used in the stochastic analyses.

Mean Std. dev.

δs (%)

δy (%)

1.825 0.1125

0.905 0.1275

The force distribution proportional to the first vibration mode has been selected. If a uniform force distribution is applied, it is very likely that damage would be concentrated at the first storey of the building. The assumption of a modal distribution, i.e. inverse triangular forces, should be more coherent with the results of time history analyses. Another assumption embedded in the performed analyses consists in neglecting damage in the spandrel beams for the individuation of the limit states of interest. The reason is that the model has been calibrated on experimental tests on masonry piers, while results concerning spandrel beams were not available (as it is usually the case). It should be noticed that in any case this procedure has a methodological purpose: in order to obtain more definitive results all the assumptions should be better tested and the other options of force distribution and direction of analysis should be explored. In some cases, results of pushover analyses may be strongly affected by the choice of the control node; however, this is not an issue for this building, since floors are rigid and hence any node at the top storey can be assumed as the control one [28]. Two different sets of stochastic pushover analyses have been carried out, each one consisting of 1000 analyses. In the first set the geometry was assigned, the materials were assigned to each macro-element, but the mechanical properties were considered as random variables varying within the ranges defined at the previous section. In the second set, not only are the mechanical parameters random variables, but also the material associated to each macroelement is randomly selected from a predefined library of materials. In particular, 30 different materials have been defined for the 165 structural elements of the model. Each material is characterised by values of the mechanical parameters randomly selected from the intervals. In each analysis, materials are redefined. The number of analyses performed for each of the two sets (1000) was considered sufficient, since convergence of each input and output variable to the corresponding mean value and stabilisation of the standard deviations were observed. The comparison of the pushover curves obtained for the two sets of stochastic analyses shows that the mean curve is not so different in the two cases. However, the dispersion in the base shear is more reduced in the second case, up to the failure range of the curves; after this point the dispersion is larger than in the first case. Moreover, the dispersion associated to the maximum force limit state is significantly reduced in the second set of analyses.

Fig. 7. Results of the second set of stochastic pushover analyses, with identification of the mean and mean plus or minus one standard deviation curves.

Fig. 8. Identification of global damage states based on the second set of pushover analyses.

However, the dispersion associated to the ultimate limit state (80% of the maximum shear) is slightly increased. Fig. 7 shows the pushover curves obtained from the second set of analyses, just to give an idea of the variability observed in the results. Also the mean pushover curve and the mean plus or minus one standard deviation are shown in the same figure. From each pushover curve it is possible to identify two global limit states, one corresponding to the attainment of the maximum base shear resistance, the other corresponding to the ultimate condition, identified by a value of base shear deterioration corresponding to 20% of the maximum base shear. Since a set of 1000 analyses is available, it is possible to define a distribution of values for these two limit states. Both the mean value and the mean plus or minus one standard deviation values obtained form the analyses are identified in Fig. 8 by vertical lines. In order to derive a probability density function associated to damage state DS2, it has been necessary to define a relationship (obviously in probabilistic terms) between the global displacement of the structure and the corresponding maximum drift observed in the elements. The results for the performed analyses are reported in Fig. 9, in which also the mean and the mean plus or minus one

1318

M. Rota et al. / Engineering Structures 32 (2010) 1312–1323

the figure. Such areas correspond to the integral of the joint probability distributions. The probability density functions associated to each of the four considered damage states have hence been determined. They are summarised in Fig. 11. 6. Incremental dynamic analyses 6.1. Selection of accelerograms

Fig. 9. Relationship between global structural displacement and maximum element drift. The thick continuous line represents the mean curve while the thick dashed lines represent the mean plus or minus one standard deviation curves.

standard deviation curves are highlighted. It can be noticed that the results are not much scattered from the mean curve and that the curve dispersion is not symmetrical around the mean curve. At this point, for each value of global displacement, the probability density function of the maximum value of drift in the elements is known. What is needed is the probability that this value of drift exceeds the drift threshold associated to the considered damage states. As indicated in a previous section, the drift threshold of each damage level is also defined in terms of a probability density function. In order to obtain the final probability distribution relating DS1 and DS2 to the global structural displacement, it is necessary to perform the convolution of the two corresponding probability distributions. Fig. 10 qualitatively describes the procedure followed for such convolution: for each global displacement value, the drift demand on the elements is evaluated and represented by a probability density function of values (right part of Fig. 10). The complementary cumulative distribution function of such values (dashed line in the left part of the figure) is then convolved with the probability density function of the two limit state thresholds (two bell curves in the left part of the figure). The result of the convolution for each DS is represented by the shaded areas in

In order to carry out incremental dynamic analyses, an appropriate set of acceleration time histories is required. As well known, accelerograms may be selected from data banks of real accelerograms or alternatively they can be generated synthetically. The advantage of real records over artificial or synthetic accelerograms is that genuine records of ground motion shaking are more realistic since they carry all the ground motion characteristics (amplitude, frequency, energy content, duration, number of cycles and phase) and reflect all the factors that influence accelerograms (source, path and site). In this study, the accelerograms adopted for the dynamic analyses have been selected from online strong motion record databases (www.isesd.cv.ic.ac.uk/, peer.berkeley.edu/smcat/, db. cosmos-eq.org), with the constraint of the spectrum-compatibility with a target response spectrum. A single-record time history analysis cannot fully capture the behaviour of a building, since results may be strongly dependent on the record chosen (e.g. [29]). So, a minimum of 7 accelerograms must be applied to the structure to be allowed to use average results instead of the most unfavourable ones, as suggested by several modern seismic codes [17,30,31] and by research works (e.g. [32]). The 7 real spectrum-compatible accelerograms used for time history analyses are plotted in Fig. 12. They have been selected through an algorithm based on a Monte Carlo random selection of groups of accelerograms: the program automatically combines the records downloaded from the strong motion databases and identifies the set best reproducing the target response spectrum. The approach followed is described in detail in [33]. Finally, the selected accelerograms have to be scaled to the target PGA in order to have a good fitting of the mean response

Fig. 10. Scheme of the procedure used for the identification of drift-dependent limit state probabilities.

M. Rota et al. / Engineering Structures 32 (2010) 1312–1323

1319

6.2. Results of the incremental dynamic analyses

Fig. 11. Probability density functions of the different considered damage states.

spectrum with respect to the target spectrum. There is widespread concern in the engineering community regarding the practice of scaling records. For a detailed discussion on the methods adopted in the literature for scaling accelerograms and their limitations, the reader is referred for example to [34]. In this work, records have been scaled linearly in order to match their PGA with the target PGA of the selected response spectrum. Notice that the scale factors applied to the accelerograms were all very close to one, within 0.6 and 1.1. Since the concern usually has something to do with ‘‘weaker’’ records not being representative of ‘‘stronger’’ ones [29] and since the maximum PGA value considered is 0.3 g (not much larger than the spectrum PGA of 0.25 g), it is believed that the considered case is not one of the most critical ones. A more detailed study of this issue will be definitely a further development of the proposed methodology, for example by selecting a different set of ground motion records for each PGA level, without scaling them. Fig. 13 shows the comparison of the average response spectrum of the selected accelerograms and the target code spectrum. The code spectrum is the one of the EC8 (type 1) for a PGA of 0.25 g.

Each time history analysis with a given record provides a value of the selected demand parameter. If the analysis is repeated using a different record, another value of this parameter is obtained. As the number of analyses increases, it is possible to define a probability distribution of values of this demand measure. If this process is repeated for a range of ground motion values, for example obtained by scaling each record to different values of PGA, it is possible to obtain a structural response function indicating the probable distribution of the engineering demand parameter for different levels of ground motion. This is the so-called incremental dynamic analysis [29]. A set of incremental time history analyses has been performed by applying to the structural model the accelerograms defined in the previous section. A value of Rayleigh critical viscous damping equal to 2% has been assumed for all the analyses, since most of the energy dissipation was provided by the hysteretic response. The aim of this study is to derive fragility curves for a typology of structures, identified by the selected building prototype, for which the ground motion and material properties are prone to variability. In order to take into account both the variability due to different ground motions (7 selected accelerograms) and the variability of mechanical parameters, a very large number of analyses would be required, wherein, for each time history, a number of material properties is extracted from properly defined probability density functions by stochastic Monte Carlo simulations. Since each nonlinear time history analysis requires a significant computational time, it is necessary to reduce the number of analyses by making some assumptions. A first set of time history analyses has been carried out using the mean value of each mechanical parameter, as reported in Table 3, and applying all the 7 accelerograms scaled so that the mean PGA matches 0.05, 0.1, 0.15, 0.2, 0.25 and 0.3 g. These analyses will be referred in the following as the deterministic set. This set corresponds to the assumption that the dispersion

Fig. 12. Spectrum-compatible real accelerograms used for time history analyses.

1320

M. Rota et al. / Engineering Structures 32 (2010) 1312–1323

Fig. 14. Comparison of standard deviation values obtained with mean mechanical parameters and 7 accelerograms, with stochastic parameters and only 1 accelerogram and combining the two cases.

Fig. 13. Comparison of the mean response spectrum of the selected records and the code spectrum.

in displacement demand due to the variability of mechanical parameters is negligible with respect to the dispersion introduced by the use of different ground motions. A second set of analyses has been carried out taking into account the variability due to material parameters and neglecting the variability induced by the different ground motion records. In this case, 100 realisations of each material parameter have been generated through Monte Carlo simulations and a corresponding number of analyses has been carried out for the values of PGA = 0.1, 0.2 and 0.3 g and applying only the first accelerogram. In order to check the assumption that the variability related to ground motion has a stronger influence on the displacement demand than the one due to mechanical parameters, for each considered PGA level the total standard deviation associated to displacement demand has been calculated as the square root of the sum of the squares of the two standard deviations associated respectively to ground motion and mechanical parameter variability. A similar approach has been followed in several works for combining the uncertainty associated to capacity and demand in order to obtain the total damage state variability (e.g. [2,35]), since the latter was not directly available from the analyses as it is instead the case in this work. The results are summarised in Fig. 14, where the total standard deviation obtained is compared to the standard deviations corresponding to the deterministic case (mean mechanical parameters and 7 accelerograms) and the probabilistic case (stochastic mechanical parameters and a single accelerogram). Observation of the plot shows that, as expected, the standard deviation of the displacement demand increases with the severity of ground motion due to the increase of nonlinear behaviour, which modifies the structural response of the building. Moreover, the figure demonstrates that, for this prototype and with the selected ground motions, the variability due to mechanical properties is significantly smaller than the variability associated to different ground motions, although the number of analyses in the former case (100 realisations of each material parameter for each PGA) is significantly larger than in the latter (7 accelerograms for each PGA). Hence, it seems acceptable to neglect the variability associated to different values of material parameters and work with the deterministic hypothesis. The results reported from now on will refer to this case. However, it has to be noticed that this assumption has been tested only with reference to the first accelerogram. Further analyses, with different ground motion records, would be needed to confirm the results obtained. This assumption is actually confirmed by the results of many researchers who identified the randomness in earthquake excitation as the most significant

Fig. 15. Complementary cumulative distribution functions of the displacement demand at various PGA levels.

Fig. 16. Identification of maximum displacement demand (mean and mean plus or minus one standard deviation) for a PGA of 0.1 g.

source of uncertainty in the probabilistic assessment of structures (e.g. [36] and references therein). At this point, the probability density functions of displacement demand at different values of PGA are needed, since they will be convolved with the distributions of each damage state in order to derive fragility curves. In particular, the complementary cumulative distributions reported in Fig. 15 will be used for such purpose. Such complementary cumulative distributions of maximum displacement demand have been obtained from time history analyses with the 7 accelerograms at each value of PGA. From the results of such analyses, the parameters of the corresponding lognormal distributions have been obtained. As an example, Fig. 16 shows the hysteretic cycles obtained from all the time history analyses with the 7 accelerograms scaled to a PGA of 0.1 g: the mean and mean plus or minus one standard deviation value of the maximum displacement demand obtained from the different analyses are indicated by vertical lines. 7. Derivation of fragility curves Once the probability density functions of the different considered damage states have been defined from pushover analyses and

M. Rota et al. / Engineering Structures 32 (2010) 1312–1323

1321

Fig. 17. Derivation of the fragility points for a PGA = 0.25 g, from convolution of pdf of limit states and complementary CDF of demand.

the displacement demand imposed on the building by different levels of ground motion has been evaluated from time history analyses, it is possible to convolve these two results to obtain fragility curves. As reported for example by [20], the probability of failure associated to a predefined limit state can be defined as: Pf (PGA) =

Table 6 Parameters of the lognormal distributions fitting the analytical fragility points. DS

µ

σ

DS1 DS2 DS3 DS4

−2.026 −1.645 −1.351 −1.169

0.362 0.273 0.218 0.175

∞

Z

1 − FD|PGA (α) fC (α)dα





(1)

0

where FD|PGA is the cumulative density function (cdf) of the demand for a given level of PGA and fc is the probability density function (pdf) of the capacity, i.e. of the damage state under consideration. The procedure is graphically illustrated in Fig. 17. The top part of the figure shows the complementary cumulative distribution function of the demand (dash-dotted line) for a PGA of 0.25 g and the probability density functions of the various damage states (continuous bell curves). These curves are convolved and the obtained curves are shown, for each damage state, with a dashed bell curve under the corresponding continuous bell curve. The integral of the area below the curve representing the result of the convolution gives a value of probability, i.e. the probability of exceedance of that damage state for a PGA of 0.25 g. These values are then reported as points in the lower part of the figure and represent the fragility points corresponding to this PGA level. The procedure described above is then repeated for all the considered PGA values. The integral of the areas below the curves obtained from the convolution of capacity and demand for the different PGA levels provide the fragility points representing, for each PGA level,

the probability of exceeding each damage state. The fragility points so obtained have been fitted using a lognormal probability density function, which is very often adopted to describe fragility curves. Obviously other types of probability distributions could be used to fit the fragility points, possibly leading to even better approximations. Also, additional points corresponding to intermediate values of PGA could be added although the five considered levels appear to be sufficiently well fitted by the lognormal probability distribution. The points obtained and the lognormal curves fitting them are plotted in Fig. 18, while the parameters (µ and σ ) of the lognormal curves are summarised in Table 6. 8. Comments and conclusions A methodology is presented for obtaining analytical fragility curves for classes of masonry buildings from detailed 3D nonlinear dynamic analyses of entire structures. This work was motivated by the scarcity of research specifically addressed to the derivation of analytical fragility curves for masonry buildings using detailed models and rigorous analysis techniques. The existing literature

1322

M. Rota et al. / Engineering Structures 32 (2010) 1312–1323

Fig. 18. Lognormal fragility curves fitting the analytically derived points.

approaches are indeed often based on simplified models of the buildings and approximate analysis types (linear or nonlinear static analyses). The proposed approach is based on the results of nonlinear stochastic analyses of a prototype building, considered to be representative of a building typology. All the mechanical properties of the structure are assumed to be random variables, to which normal probability distributions are associated based on realistic ranges of variation. Even the association of mechanical parameters to each structural element is treated as a random variable. Monte Carlo simulations are then used to extract the input parameters from such distributions and then advanced nonlinear analyses of the whole structure are performed. The probability density functions of selected damage thresholds are determined based on pushover analyses and then convolved with the probability distribution of displacement demand obtained from nonlinear time history analyses. Lognormal fragility curves are hence obtained by fitting the fragility points representing the probability of exceeding different damage levels for discrete levels of PGA. The main advantages of this approach consist in:

• The definition of the probability distributions of global displacement limit states through nonlinear static analyses on stochastically defined structural models; • The identification of the variability of material characteristics based on the analytical simulation of the results of cyclic experimental tests; • The nonlinear dynamic capability of the structural software, which allows to carry out accurate analyses of the whole structure with a reasonable computational effort; • The definition of the probability distribution of displacement demand through the use of a set of real spectrum-compatible accelerograms for performing time history analyses, which allow to include all the ground motion characteristics (amplitude, frequency, energy content, duration, number of cycles and phase) and reflect all the factors influencing accelerograms (source, path and site). Due to the definition of the damage state DS1, the obtained fragility curves predict a null damage for values of PGA smaller than 0.05 g. This is not always realistic, since in many cases damage has been observed for lower values of ground motion. One of the reasons of this lower estimated vulnerability is related to the definition of the damage state DS1 based on visual interpretation of experimental tests, in which typically the external surface of the tested specimens is not painted nor plastered and hence low levels of damage are not easy to detect. On the other hand, in

existing masonry buildings, some level of damage can usually be noticed even prior to the occurrence of an earthquake, due to poor maintenance conditions [1]. Dynamic time history analyses are all performed with a value of 2% equivalent viscous damping. This assumption should be verified more in detail: analytical models, indeed, often tend to underestimate the effective dissipation capacity of the structure due to the shape and type of the implemented hysteretic model, particularly when damage propagates into the structure. The assumption of 2% equivalent viscous damping, if appropriate for low values of PGA, may be too low for higher intensities of ground motion, where the behaviour of the building becomes highly nonlinear and the hysteretic model may no longer be able to adequately reproduce the dissipation capacity of the building. The consequence of this may be an overestimation of the related vulnerability. Although this effect is normally on the conservative side if analyses are used for the aim of assessing or designing a structure, in this case the overestimation of vulnerability may determine an incorrect evaluation of seismic risk. To compensate for this numerical effect, a higher value of equivalent viscous damping may be appropriate. Moreover the local flexural out-of-plane behaviour of walls is not modelled, since a global box-type behaviour of the building is assumed. In this case the wall out-of-plane flexural contribution to global stiffness and resistance is generally negligible with respect to the in-plane one. This assumption represents a good approximation for low levels of in-plane damage, but for severely damaged structures in which the wall in-plane capacity is strongly deteriorated, it may induce a non-negligible underestimation of the seismic building performance. The proposed approach can be considered appropriate for structural typologies with good connections between orthogonal walls and between walls and floors and with rigid diaphragms (such as the prototype building considered in this study) whose behaviour is dominated by a global response governed by in-plane mechanisms. In case of buildings with inappropriate connections and lack of any specific device preventing local collapses (e.g. tie rods, tie beams, etc.), consideration of such local failure modes (mainly associated to out-of-plane response) needs to be incorporated in the assessment procedure. It is important to note that time history analyses have been performed using accelerograms compatible with the response spectrum for rock soil and hence they refer to a predefined spectral shape (EC8, type 1) and they do not represent all the possible stratigraphic and topographic conditions occurring in real cases. A possible improvement of the proposed incremental time history analysis would consider the selection of a wider database of natural records able to cover a significant ground motion severity range without any need for scaling. Also different soil conditions would be taken into account. The quantitative results obtained in this study are only applicable to structures belonging to the same typology as the prototype building. However, the proposed procedure has a general validity and can be easily extended to other types of buildings, with the limitations previously discussed concerning out-of-plane behaviour. Clearly the presented approach is still at a methodological stage and needs to be further tested and improved by applying it to different building prototypes representative of the same typology and also to structures belonging to different types. Moreover, for each building prototype, several geometrical configurations could be taken into consideration by treating the variability of some geometrical characteristics by means of additional stochastic variables. Acknowledgements The authors would like to thank Prof. Alex Barbat of the CIMNE (Centro Internacional de Métodos Numéricos en la Ingeniería), Universitat Politècnica de Catalunya, Barcelona (Spain), for kindly

M. Rota et al. / Engineering Structures 32 (2010) 1312–1323

providing the program STAC used to perform Monte Carlo simulations. The financial support of the Reluis Project (Line 10) ‘‘Definition and development of databases for risk evaluation and emergency management and planning’’ is also acknowledged. Finally, the manuscript has benefited from the thorough reviews of two anonymous reviewers. References [1] Rota M, Penna A, Strobbia CL. Processing Italian damage data to derive typological fragility curves. Soil Dyn Earth Eng 2008;28:933–47. [2] RISK-UE. The European risk-Ue project: An advanced approach to earthquake risk scenarios (2001–2004). 2004. www.risk-ue.net. [3] Lagomarsino S, Giovinazzi S. Macroseismic and mechanical models for the vulnerability and damage assessment of current buildings. Bull Earthq Eng 2006;4:415–43. [4] Restrepo-Velez LF, Magenes G. Simplified procedure for the seismic risk assessment of unreinforced masonry buildings. In: Proc. 13th WCEE. 2004. Paper No. 2561. [5] Magenes G. Masonry building design in seismic areas: Recent experiences and prospects from a European standpoint. In: Proc. 1st ECEES. 2006. [6] D’Ayala DF. Force and displacement based vulnerability assessment for traditional buildings. Bull Earthq Eng 2005;3:235–65. [7] Cattari S, Frumento S, Lagomarsino S, Parodi S, Resemini S. Multi-level procedure for the seismic vulnerability assessment of masonry buildings: The case of Sanremo (north-western italy). In: Proc. 1st ECEES. 2006. [8] Kappos AJ, Panagopoulos G, Panagiotopoulos C, Penelis G. A hybrid method for the vulnerability assessment of R/C and URM buildings. Bull Earthq Eng 2006;4:391–413. [9] Bothara JK, Mander JB, Dhakal RP, Khare RK, Maniyar MM. Seismic performance and financial risk of masonry houses. ISET J Earthq Technol 2007; 44(3–4):421–44. [10] Erberik MA. Generation of fragility curves for Turkish masonry buildings considering in-plane failure modes. Earthq Eng Struct Dyn 2008;37(3):387–405. [11] Pasticier l, Amadio C, Fragiacomo M. Non-linear seismic analysis and vulnerability evaluation of a masonry building by means of the sap2000 v.10 code. Earthq Eng Struct Dyn 2008;37(3):467–85. [12] Park J, Towashiraporn P, Craig JI, Goodno BJ. Seismic fragility analysis of low-rise unreinforced masonry structures. Eng Struct 2009;31(1):125–37. [13] Galasco A, Lagomarsino S, Penna A. TREMURI program: Seismic analyses of 3d masonry buildings. Italy: University of Genoa; 2009. [14] Penna A. A macro-element procedure for the non-linear dynamic analysis of masonry buildings. Ph.D. dissertation. Italy: Politecnico di Milano; 2002 [in Italian]. [15] Tomaževi˘c M, Klemenc I. Seismic behaviour of confined masonry walls. Earthq Eng Struct Dyn 1997;26:1059–71. [16] Tomaževi˘c M. In: Bernardini A, editor. Seismic damage to masonry buildings. Correlation between damage and seismic resistance of masonry walls and buildings. Rotterdam: Balkema; 1999.

1323

[17] OPCM 3274. Ordinanza del Presidente del Consiglio dei Ministri n. 3274 del 20 Marzo 2003: Primi elementi in materia di criteri generali per la classificazione sismica del territorio nazionale e di normative tecniche per le costruzioni in zona sismica. GU n. 72 del 8-5-2003, with further modifications in OPCM 3431. 2005 [in Italian]. [18] Rossetto T. Vulnerability curves for the seismic assessment of reinforced concrete building populations. Ph.D. dissertation. London: Imperial College; 2004. [19] Rubinstein RY. Simulation and the Monte Carlo method. New York: Wiley; 1981. [20] Pinto PE, Giannini R, Franchin P. Methods for seismic reliability analysis of structures. Pavia (Italy): IUSS Press; 2004. [21] Zárate F, Hurtado J, Oñate E, Rodríguez J. STAC program: Stochastic analysis computational tool, CIMNE, Centro Internacional de Métodos Numéricos en la Ingeniería. Barcelona, España: 2002. [22] Lagomarsino S, Penna A. A nonlinear model for pushover and dynamic analysis of masonry buildings. In: Proc. int. conf. computational and experimental engineering and sciences. 2003. [23] Penna A, Cattari S, Galasco A, Lagomarsino S. Seismic assessment of masonry structures by nonlinear macro-element analysis. In: Proc. 4th int. seminar SAHC. 2004. [24] Lagomarsino S, Galasco A, Penna A. Nonlinear macro-element dynamic analysis of masonry buildings. In: Proc. ECCOMAS conf. computational methods in structural dynamics and earthquake engineering. 2007. [25] Gambarotta L, Lagomarsino S. On dynamic response of masonry panels’. Gambarotta L, editor. Proc. Nat. Conf. La meccanica delle murature tra teoria e progetto. 1996 [in Italian]. [26] Faella G, Manfredi G, Realfonzo R. Comportamento sperimentale di pannelli in muratura di tufo sottoposti ad azioni orizzontali di tipo ciclico. In: Proc. 5th it. conf. earthquake engineering. 1991 [in Italian]. [27] ENV 1998-3 2005 Eurocode 8: Design of structures for earthquake resistance — Part 3: Assessment and retrofitting of buildings. [28] Galasco A, Lagomarsino S, Penna A. On the use of pushover analysis for existing masonry buildings. In: Proc. 1st ECEES. 2006. Paper No. 1080. [29] Vamvatsikos D, Cornell CA. Incremental dynamic analysis. Earthq Eng Struct Dyn 2002;31:491–514. [30] UBC Uniform Building Code. Structural engineering design provisions. In: International conference of building officials, vol. 2. 1997. [31] ENV 1998-1 2005 Eurocode 8: Design of structures for earthquake resistance — Part 1: General rules. seismic actions and rules for buildings. [32] Bommer JJ, Acevedo AB, Douglas J. The selection and scaling of real earthquake accelerograms for use in seismic design and assessment. In: Proc. ACI Int. conf. seismic bridge design and retrofit. La Jolla(California): American Concrete Institute; 2003. [33] Dall’Ara A, Lai CG, Strobbia C. Selection of spectrum-compatible real accelerograms for seismic response analyses of soil deposits. In: Proc. 1st ECEES. 2006. Paper No. 1240. [34] Bommer J, Acevedo AB. On the use of real earthquake accelerograms as input to dynamic analysis. J Earthq Eng 2004;8(Special Issue 1):43–91. [35] FEMA. HAZUS earthquakes loss estimation methodology. Washington: US federal emergency management agency; 1999. [36] Dymiotis C, Kappos AJ, Chryssanthopoulos MK. Seismic reliability of RC frames with uncertain drift and member capacity. ASCE J Struct Eng 1999;1038–47.