Effects of ground-motion sequences on fragility and vulnerability of case-study reinforced concrete frames

Abstract

This study investigates the effects of ground-motion sequences on fragility and vulnerability of reinforced concrete (RC) moment-resisting frames (MRFs). Two four-storey, four-bay RC MRFs are selected as case studies. These structures, which share the same geometry, are representative of distinct vulnerability classes in the Mediterranean region and are characterized by different material properties, cross-section dimensions, and detailing. The first case study is a ductile MRF designed according to Eurocode 8 (i.e., a special-code frame), while the second is a non-ductile MRF designed to sustain only gravity loads (i.e., a pre-code frame). The influence of masonry infills on their seismic performance is also investigated. Advanced numerical models are developed to perform cloud-based sequential nonlinear time history analyses using ground-motion sequences assembled by randomly pairing two real records via Latin hypercube sampling. Different structure-specific damage states are considered to derive fragility curves for the undamaged structures, when subjected to a single ground-motion record, and state-dependent fragility curves by considering the additional damage induced by a second ground-motion record within the sequence. Damage-to-loss models are then used to derive mean vulnerability relationships. Results of the analysis show the importance of considering the effect of damage accumulation in the pre-code frames. Moreover, the presence of infills shows an overall positive contribution to the seismic performance of both frame types. Vector-valued vulnerability relationships accounting for the damaging effect of two ground-motion records are finally presented in the form of mean vulnerability surfaces.

1 Introduction and motivation

Earthquakes typically occur in sequences characterized by a seismic event with a larger magnitude than all others. Such major events constitute mainshocks, generally followed by several smaller-magnitude aftershocks clustered in space and time. This may not strictly apply to triggered earthquakes where the energy/stress due to the seismic waves generated by a fault rupture triggers a distant earthquake with a different rupture zone (and potentially similar or even larger magnitude). These sequences can cause vast losses due to direct repair costs, business interruption and casualties, especially if the affected structures/infrastructure are left unrepaired upon initial damage due to the first event (often due to the short time between events).

Recent worldwide earthquakes have demonstrated the devastating effects of such sequential earthquake-induced ground shaking. On September 2010, the city of Christchurch in New Zealand, was hit by a mainshock with a moment magnitude (M_w) of 7.1, which was followed, in February 2011, by an aftershock with M_w of 6.3. This seismic sequence resulted in 185 casualties and nearly $15 billion financial losses (Parker and Steenkamp 2012). After the September mainshock, 90% of reinforced concrete (RC) structures located in the central business district of Christchurch were classified as safe to re-occupy with none or minor damage (Kam et al. 2010). Following the aftershock of February, only 48% of those structures were classified as safe to re-occupy, whilst the portion of unsafe (i.e., red-tagged) ones raised to 19% and approximately 33% had limited accessibility or restricted use (Kam et al. 2011).

A similar situation was observed in the 2016–2017 Central Italy earthquake sequence, which consisted of several moderate-to-high magnitude earthquakes, each centered in a different location and with its own sequences of aftershocks spanning several months (e.g., Zimmaro and Stewart 2017). The first event on the 24th of August 2016 with M_w of 6.1 was followed by two seismic events in October 2016 with M_w of 5.9 and 6.5 respectively. The sequence hit a large sector of the Central Apennines of Italy, particularly affecting Marche and Umbria regions, and significantly damaging several towns such as Amatrice, Norcia and Accumoli (e.g., Stewart et al. 2018). For instance, after the August 2016 event, 72% of the buildings in Accumoli town showed an acceptable performance and remained undamaged; 24% had damages ranging from negligible to heavy; whilst only 4% collapsed. Following the two events of October 2016, the percentage of collapsed buildings raised to 65%, as reported in Fig. 1.

Fig. 1

Building damages in Accumoli after the first event of August and after the two events of October 2016 (modified after Zimmaro and Stewart 2017)

These examples, among many others, highlight the devastating effects of damage accumulation on building performance due to seismic sequences. They also show the potential of aftershocks (or triggered events) to lead to higher damage or even collapse of mainshock-damaged buildings. Many studies in the last few years focused on these aspects; a comprehensive overview of them is outside the scope of this paper, but some key findings are briefly reviewed here.

A number of authors investigated the issue of sequential ground shaking adopting single-degree-of-freedom (SDoF) systems (e.g., Hatzigeorgiou and Beskos 2009; Goda and Taylor 2012; Di Sarno 2013; Goda et al. 2015). Amongst others, Hatzigeorgiou and Beskos (2009) quantified the inelastic displacements due to ground-motion records repeated multiple times (i.e., back-to-back analysis), whilst Di Sarno (2013) used 15 as-recorded sequences. Both studies found that the deformation demand due to multiple ground motions increases significantly (40–400%) in many analysis cases compared to those associated with single events.

Additional studies focused on RC structures modelled as multi-degrees-of-freedom (MDoF) systems (e.g., Hatzigeorgiou and Liolios 2010; Raghunandan et al. 2015; Jeon et al. 2015; d’Aragona et al. 2017). Jeon et al. (2015) considered RC frames designed using an old seismic code and performed incremental dynamic analysis (IDA) (Vamvatsikos and Cornell 2002) using 100 mainshocks to bring the structures to prescribed initial damage states (DSs). Next, the authors randomly chose unscaled records applied as artificial aftershocks via cloud analysis. Similarly, Raghunandan et al. (2015) analyzed ductile RC frames via IDA and used artificial mainshock–aftershock (MS–AS) sequences. Both studies showed that the initial mainshock-induced damage make the considered structures more fragile as indicated by the leftward shift of the derived aftershock fragility curves. Specifically, Jeon et al. (2015) showed that the median collapse capacity of mainshock-damaged frames, i.e., ground-motion intensity causing collapse, is reduced up to 34% compared to the undamaged (i.e., as-built) frames when the initial mainshock DS was moderate to severe; Raghunandan et al. (2015) found a reduction up to 46%. Both studies showed that high capacity reduction generally occurs when the mainshock DS is significant, whilst no notable reduction is generally observed for low mainshock damage.

Many studies also addressed the influence of infills on the performance of RC frames subjected to multiple excitations (e.g., Tesfamariam et al. 2015; Burton and Sharma 2017; Furtado et al. 2018). Amongst others, Tesfamariam et al. (2015) assessed the inelastic response of seismically-designed bare and infilled frames using real MS–AS sequences scaled to match a seismic-hazard target intensity. They showed that uniformly distributed infills improve the seismic performance of the considered case-study structures, i.e., the infilled frame can achieve the immediate occupancy performance objective in most analysis cases, whereas the bare one reaches either life safety or collapse prevention in 90% of the cases. The study also found that ≤ 10% of the aftershocks push the frames from immediate occupancy to higher damage levels. Burton and Sharma (2017) assessed seismically-designed infilled RC frames via IDA using artificial MS–AS sequences. They found that, when the mainshocks are scaled up to 50% of the median collapse capacity of the undamaged structures, the collapse fragility of the mainshock-damaged structures does not notably increase. The opposite is observed when the mainshocks were scaled to higher intensity levels and a reduction of 26–63% is observed.

Although many of the studies discussed above addressed fragility analysis against seismic sequences, none of them explicitly investigated vulnerability relationships (i.e., likelihood of loss levels as a function of a hazard intensity measure) characterizing the increase in expected losses due to sequential seismic excitations. Moreover, most of these studies derived aftershock fragility curves accounting for (near) collapse DS only without considering other DSs that also contribute to seismic losses in probabilistic seismic risk assessment exercises.

Other studies focused on different structural configurations and materials (e.g., Li et al. 2014; Ruiz-García and Aguilar 2015; Di Trapani and Malavisi 2019; Zhang et al. 2020). For example, Ruiz-García and Aguilar (2015) analyzed a ductile seismically designed steel moment-resisting frame (MRF) via IDA and adopted real MS–AS sequences considering the influence of mainshock-induced residual drift. The study resulted in similar observations to those highlighted above in terms of collapse capacity reduction. Zhang et al. (2020) developed state-dependent fragility curves for wood-frame houses subjected to real MS–AS sequences scaled with moderate scaling factors. They used residual drift to measure mainshock DS and maximum first-storey drift for aftershock DS. The authors showed that the median values of aftershock fragility curves are reduced by 14–46% compared to the undamaged state and the reduction is much more significant when the initial mainshock DS is large. Some studies have also focused on different aspects such as MS–AS record/sequence selection, with special focus on hazard consistency (e.g., Goda 2015; Papadopoulos et al. 2020).

In summary, the seismic performance and fragility of structures subjected to seismic sequences has been largely investigated in the last decade. Past research efforts focused mainly on seismic performance assessment and changes in collapse fragility due to such sequences. Only a few studies developed state-dependent fragility curves for seismically damaged structures considering different DSs with structure-specific DS thresholds. Moreover, as pointed out above, very limited attention has been given to the development of vulnerability relationships that could be practically implemented in seismic risk assessment of building portfolios against seismic sequences.

Accordingly, this paper develops state-dependent fragility and vulnerability relationships for ground-motion sequences. Direct monetary losses (i.e., in terms of repair costs) is the considered loss metric. Two case-study structures, representing widespread vulnerability classes in Europe and particularly Italy, are selected: (1) ductile RC MRFs (i.e., special-code design level) designed as per the seismic provisions of Eurocode 8 (EC8) (EN 1998-1 2004); (2) non-ductile (i.e., a pre-code design) RC MRFs designed for gravity loads only. Both frames are investigated in both bare and uniformly-infilled configurations to understand the effects of infills on seismic-demand estimates, damage and vulnerability. Advanced nonlinear models are developed using OpenSees (Mazzoni et al. 2009) and artificial ground-motion sequences are assembled via a randomized approach and then used to perform sequential cloud-based nonlinear time history analyses (NLTHAs). The results are used to develop fragility curves for both the undamaged and seismically-damaged structures considering three structure-specific DSs. The derived fragility curves are coupled with damage-to-loss models, which are used in the literature for generic Italian buildings, to generate vulnerability relationships.

This paper is organized as follows: in Sect. 2 the general methodology of the work is presented with detailed discussion of some important aspects such as ground-motion selection and derivation of fragility and vulnerability relationships; Sect. 3 describes the four case-study structures, nonlinear modelling strategies and DSs thresholds mapping; Sect. 4 addresses the performance assessment and fragility curves derivation considering undamaged and seismically-damaged structures. Finally, vulnerability relationships are derived with an attempt to generate vector-valued vulnerability relationships (i.e., surfaces); Sect. 5 provides the main conclusions drawn from this study.

2 Methodology

The methodology implemented in the current study to perform fragility and vulnerability assessment against seismic sequences consists of the following steps:

• 1. Identification of the layout, geometric features, material characteristics and detailing of case-study structures representative of the considered vulnerability classes for the geographic area of interest;

• 2. Development of advanced numerical models using OpenSees to simulate the nonlinear hysteretic response of the case-study structures under seismic sequences;

• 3. Definition of structure-specific DSs, which are used to describe the damage conditions (i.e., the performance objectives) of the considered structures, and of DS thresholds in terms of global engineering demand parameters (EDPs), i.e., the maximum interstorey drift ratio (MIDR) among all storeys during the entire NLTHA. Three different DSs, i.e., DS1, DS2 and DS3, are considered in this study to account for increasing levels of damage. The thresholds of these DSs are identified based on pushover analyses and local EDPs, e.g., rotation, strength of cross-sections, and material strains (e.g., Freddi et al. 2017);

• 4. Definition of artificial ground-motion sequences by randomly pairing two real records via Latin hypercube sampling (LHS) (McKay et al. 1979). The sequences are formed of two ground motions each, named GM1 and GM2 hereinafter;

• 5. Performing cloud-based sequential NLTHAs (e.g., Jalayer et al. 2014) to derive seismic demand estimates in terms of the EDPs of interest;

• 6. Derivation of fragility and vulnerability relationships for the case-study structures in the undamaged state, i.e., when subjected to the initial ground shaking (i.e., GM1) only, and then derivation of state-dependent fragility curves considering the additional damage induced by the following record in the sequence (i.e., GM2);

• 7. Derivation of vector-valued vulnerability relationships accounting for the damaging effects of the two ground motions and their combined effect on the expected losses.

The following subsections highlight further details regarding cloud analysis, fragility curve derivation, ground-motion selection, and derivation of vulnerability relationships.

2.1 Cloud analysis and fragility curve derivation

The geometric mean of the 5%-damped pseudo-spectral acceleration (avgSa) over a selected range of periods is the selected ground-motion intensity measure (IM) in this study. This IM is used to indirectly account for the higher-mode effects and period elongation due to the damage of non-structural components (i.e., infills) and stiffness and strength degradation of the structure (Baker and Cornell 2006; Eads et al. 2015; Kazantzi and Vamvatsikos 2015; Kohrangi et al. 2017). When compared to other IMs more widely used in the literature (e.g., spectral acceleration at the fundamental period of the structure, ${\text{S}}_{{\rm a}}\text{(}{\text{T}}_{{\rm 1}}\text{))}$, avgSa shows better efficiency, minimizing the response variability conditioned to a given IM level (Kohrangi et al. 2017), and it has higher (relative) sufficiency (e.g., Minas and Galasso 2019), yet being hazard computable.

To account for damage accumulation due to sequential ground shaking, EDPs monotonically increasing with respect to the duration of the excitation must be used, such as the hysteretic energy or the cumulative ductility demand (e.g., Morfuni et al. 2019; Gentile and Galasso 2020). However, it is challenging to correlate such EDPs with (observed) structural and non-structural damage to define DS thresholds. Moreover, the hysteretic energy dissipation characteristics depend on the specific plastic mechanism and number of structural components contributing to the nonlinear behavior, which makes it difficult to generalize the results. Consequently, this study uses MIDR as the considered EDP even though it does not explicitly quantify damage accumulation (e.g., related to cyclic strength degradation). Despite this disadvantage, MIDR remains a reliable and widely-used proxy to account for structural and non-structural damage (Raghunandan et al. 2015) and to define DS thresholds. Moreover, when damage accumulation is properly incorporated in the nonlinear model through plastic hinges, hysteretic dissipation and cyclic degradation modelling, the MIDR values resulting from the analysis can be used to indirectly reflect damage accumulation and characterize global damage.

The fragility curves in this study are initially derived considering the structures in the undamaged (i.e., as-built) state by fitting an appropriate Probabilistic Seismic Demand Model (PSDM) via a power-law regression (Cornell et al. 2002) to the IM vs. EDP clouds resulting from subjecting the structures to GM1 records only. Successively, state-dependent fragility curves are derived for the structures already damaged by GM1, which provide the likelihood of exceeding different DSs over a range of IM levels of GM2 conditioned on the initial DS of the structure caused by GM1. This process can be performed by first separating the analysis results based on the initial DS, which can be either DS1 or DS2 as shown later. Next, additional clouds are generated using the IM vs. EDP resulting from applying GM2 provided that each cloud includes the analysis cases sharing the same initial DS caused by GM1. After fitting an appropriate PSDM to each cloud, state-dependent fragility curves can be generated respectively for DS2 and DS3 given that the structure is initially in DS1 (i.e., DS2|DS1 and DS3|DS1) and for DS3 given that the structure is initially in DS2 (i.e., DS3|DS2). The few analysis cases resulting in a global collapse identified by excessive MIDR values larger than 10% and/or by non-convergence (e.g., Vamvatsikos and Cornell 2002; Raghunandan et al. 2015; Jalayer et al. 2017) are not incorporated in the fragility derivation for simplicity. It is worth mentioning that the identified collapse cases represent a very limited portion of analyses and they are not expected to introduce a significant bias in the fitted fragility functions.

2.2 Ground-motion selection and randomized sequence assembling

Ground-motion sequences used in the literature are generally in the form of real or artificial MS–AS sequences. Choosing either of them might influence the analysis results due, for instance, to the differences in ground-motion features (e.g., frequency content, spectral shape, and duration) between real and artificial aftershocks within a sequence (Song et al. 2014). Although some past studies adopted real MS–AS sequences (e.g., Ruiz-García 2012; Tesfamariam et al. 2015; Zhang et al. 2020), the available records with high aftershock intensities are generally scarce, which might limit the observations of significant damage increase due to sequential excitations. This can prevent in turn from developing statistically robust models for the derivation of aftershock fragility/vulnerability relationships conditioned on initial mainshock DS. In fact, studies utilizing such real sequences focused mainly on general assessment of seismic performance (e.g., median seismic demands) under such sequences (rather than developing fragility curves) and on understanding the effects of aftershock ground-motion characteristics on the seismic performance of case-study structures. Real sequences are usually scaled (e.g., in IDA, implementing linear scaling) if the purpose is developing aftershock fragility/vulnerability relationships.

To overcome this, some studies adopted artificial sequences via back-to-back analysis, in which the same mainshock record is applied as both a mainshock and as an aftershock (e.g., Hatzigeorgiou and Beskos 2009). This approach may yield conservative results as it assumes that the ground-motion features of aftershock-induced records are identical to those of mainshock-induced ground motions within each sequence (Li et al. 2014). Other studies tackled this issue by combining each mainshock with other mainshocks (e.g., Raghunandan et al. 2015). Another way to reduce the computational effort associated with the previous approach is to randomly select mainshocks to be seeded as aftershocks (e.g., Jeon et al. 2015; Burton and Sharma 2017). However, most of these studies (i.e., Raghunandan et al. 2015; Burton and Sharma 2017) apply the developed artificial MS–AS sequences in conjunction with the IDA approach, implying excessive scaling on the records and related issues (e.g., Miano et al. 2018). More recently, Papadopoulos et al. (2020) developed an advanced strategy combining spectrally compatible ground motions into MS-AS sequences that are representative of potential rupture scenarios, by implementing the appropriate triggering mechanisms and recently-developed MS–AS spectral acceleration correlation models.

The current study uses a randomized approach to assemble artificial ground-motion sequences representing multiple ground shaking that may occur in the form of typical aftershocks or triggered events. Almost all the records used to build the sequences are unscaled (i.e., as-recorded) with a few exceptions as specified later. These records are selected from three different databases: (1) the 2012 KKiKSK Ground-Motion Database, developed by Goda et al. (2015); (2) the database developed by Goda and Taylor (2012); and (3) SIMBAD Database developed by Smerzini et al. (2014). The issue of hazard consistency is not specifically addressed in the ground-motion selection performed here, and this can introduce some bias in the analysis results. However, the considered record selection procedure is consistent with the adopted analysis strategy (i.e., cloud analysis) that relies on selecting unscaled records covering a wide range of IM levels. This approach is deemed appropriate for applications related to seismic risk assessment of building portfolios as opposed to site- and building-specific performance-based assessment exercises. Moreover, the IM selected in this study (i.e., avgSa) has proven to be more sufficient (and efficient) than other conventional IMs (e.g., Kazantzi and Vamvatsikos 2015). This means that the analysis results will be less sensitive to seismological features (e.g., magnitude, source-to-site distance, among others), which are explicitly accounted for when performing a site-specific, hazard-consistent record selection.

The 2012 KKiKSK Database includes real MS–AS sequences developed based on the ground motions recorded by the Japanese national strong-motion network K-NET/KiK-net¹. This database includes, among others, 112 mainshocks (with two horizontal components) pertaining to crustal earthquakes. Those are characterized by: M_w ≥ 5.90, peak ground acceleration PGA ≥ 0.075 g and source-to-site distance R ≤ 300 km. The second database was developed by Goda and Taylor (2012) based on the NGA-West Database (Chiou et al. 2008), including 122 real mainshocks (with two horizontal components) related to crustal events. Those are characterized by M_w ≥ 5.0, PGA ≥ 0.04 g, peak ground velocity ≥ 1.0 cm/s, average shear wave velocity in the uppermost 30 m of soil, ${\text{V}}_{{\rm S30}}$, between 100 and 1000 m/s, lowest usable frequency of 1.0 Hz or less. The related recording stations are all free field. It should be noted that the previous two databases include real MS–AS sequences that were implemented during the initial stages of the current study. However, real sequences with aftershocks characterized by high IM values are scarce, which reduced the number of damage observations and affected the robustness of the derived state-dependent fragility curves, especially for high DSs (Aljawhari et al. 2019). Therefore, as stated above, it was decided to form artificial sequences, such that only mainshocks from these databases are considered and separated from their original sequences. Considering both horizontal components of the mainshocks from the two databases, a total of 468 records are obtained. However, a notable portion of the selected ground motions are weak, practically irrelevant for this study. After calculating avgSa between 0.17 and 1.38 s (see Sect. 3.4 for a justification to this period range), the number of records is reduced to 335 in accordance with the following criteria: (1) 30% of the sequences with avgSa smaller than 0.5 g are removed uniformly at random; (2) ground motions with avgSa smaller than 0.1 g are removed. These conditions are arbitrarily set to respond to the trade-off of having a fairly-small ground-motion set size and having records of practical engineering interest.

Finally, to include a sufficient number of strong records for capturing high levels of damage, further ground-motion records are added from the SIMBAD database (Smerzini et al. 2014). SIMBAD incorporates 467 three-component records, related to shallow crustal earthquakes with 5.0 ≤ M_w ≤ 7.3 and epicentral distances R < 30 km (Smerzini et al. 2014). The 20 records with the highest avgSa are selected from this database. Those are also arbitrarily linearly-scaled in amplitude by a factor of 1.5. As a result, the total number of selected records for this study is equal to 355.

Artificial sequences are assembled by considering combinations of two records (herein called GM1 and GM2), also allowing 40 s of free vibration in between. Considering all the possible pairs of records would lead to 126,025 assembled sequences (Fig. 2a). Therefore, to reduce the computational effort needed to perform NLTHAs, 500 random pairs from all records are selected adopting the LHS approach. The LHS is used to generate pairs of integer numbers representing the ground motion indices for the GM1 and GM2. To do so, the interval [0, 355], where 355 is the total number of ground motions, is first subdivided in 500 equal-length bins. Based on a uniform distribution, one non-integer is extracted (without replacement) from each bin and it is rounded up to obtain a ground motion index. The process is repeated twice and the results are randomly combined to generate random pairs defining the sequences. It is worth noting that this is not a full LHS design, strictly speaking, because ground motions are assembled via their index number rather than their IM levels and also some records will appear twice. However, each pair of assembled ground motions remains unique. Figure 2b illustrates the values of avgSa for GM1 (avgSa_GM1) vs. GM2 (avgSa_GM2) in each sequence, which shows a wide variety of low and high IM levels for both GM1 and GM2 in each sequence.

Fig. 2

Scatter of the avgSa values for GM1 versus GM2 considering the a all possible assemblies (126,025 sequences); b randomly selected 500 sequences via the Latin Hypercube Sampling (LHS) approach

2.3 Vulnerability relationships

A building-level approach is adopted in this study to develop vulnerability relationships, for which direct economic loss is the selected loss metric. This approach is deemed appropriate for seismic risk assessment of building portfolios, as opposed to more refined, component-based loss estimation procedures, which are more suitable for building-specific applications (e.g., FEMA P-58, Federal Emergency Management Agency 2018). Vulnerability relationships are herein expressed in terms of expected loss ratio (LR), i.e., the expected repair-to-replacement cost ratio of the building, conditional on the level of ground-shaking IM. Such functions are derived according to Eq. (1), using the total probability theorem:

$$\text{LR}\left(\text{IM}\right)=\sum _{{\rm i}=1}{\text{DLR}}_{{\rm i}} \text{P}\left(\text{DS}={\text{ds}}_{{\rm i}}|\text{IM}\right)$$

(1)

where DLR_i is the expected damage-to-loss ratio for the ith damage level DS_i, while $\text{P(DS}={\text{ds}}_{{\rm i}}\text{|IM)}$ is the probability that the DS is equal to ${\text{ds}}_{{\rm i}}$. Vulnerability relationships in this study are presented in two different ways. The first one consists of plotting vulnerability relationships providing the expected LR as a function of the ground-motion IM for the structures in the undamaged state, or the state-dependent vulnerability relationships that consider the structures with the initial DS caused by GM1. The second approach consists of developing vulnerability surfaces combining vulnerability relationships of both the undamaged state and the structures with an initial damage into a vector-valued relationship, which allows for estimating the total expected LR due to the sequence made by GM1 and GM2.

3 Case study structures

3.1 Layout, geometry and material characteristics

A four-storey, four-by-four-bays RC MRF building located in Torre del Greco (Naples, Italy) is considered for illustrative purposes. The total height is equal to 13.5 m with a first storey of 4.5 m, upper storeys of 3 m and a bay width of 4.5 m in both directions (Minas and Galasso 2019).

Four structures are derived using the same geometry. The first two case-study MRFs, termed as special-code frames, are designed and detailed according to the seismic provisions of EC8 for high ductility class. These provisions include capacity design, various requirements in terms of cross-sectional dimensions and seismic detailing to ensure a ductile global performance and preventing the formation of localized brittle failure modes. Some sample cross-sectional details are shown in Fig. 3c. For the considered site, the design seismic action for the ultimate limit state (i.e., probability of exceedance of 10% in 50 years) is characterized by the Type 1 elastic response spectrum, with a PGA equal to 0.16 g (Stucchi et al. 2011) and ground type B. The upper limit of the behavior factor (q) equal to 5.85 is used in the design. One of the two special-code frames is considered in the bare configuration, whereas the other one is considered in a uniformly-infilled configuration.

Fig. 3

Elevation of the a special-code (SB); b pre-code bare (PB) frame; and c sample cross-sectional details (cross-sections dimensions in cm)

The other two frames, termed as pre-code frames, are designed for gravity loads only as per the Royal Decree n. 2239 of 1939 (Consigliore dei ministri 1939) that regulated the structural design in Italy until 1974. Therefore, these frames do not conform to modern seismic requirements and they are characterized by a non-ductile overall behavior due to the lack of strength hierarchy rules, i.e., capacity design, sufficient amount of longitudinal reinforcement (e.g., Fig. 3c shows a column with less than 1% longitudinal reinforcement), poor confinement, and susceptibility to develop brittle failure mechanisms (e.g., shear failure). Similarly to the special-code frames, the first considered pre-code frame is a bare frame, whilst the other is a uniformly-infilled frame and both of them share the same geometry and detailing.

According to this distinction, the four case-study frames include: a special-code bare (SB) frame; a special-code infilled (SI) frame; a pre-code bare (PB) frame and a pre-code infilled (PI) frame. Figure 3 provides the frames layout including cross-sectional dimensions for the SB and PB respectively.

The average compressive strength of the concrete (f_cm) and average yield strength of the reinforcement (f_ym) for the special-code frames are respectively 37 MPa and 490 MPa. Typical average values of these properties for the pre-code frames are significantly lower (Verderame et al. 2001; Verderame et al. 2011); in particular, f_cm and f_ym are assumed respectively equal to 19 MPa and 360 MPa. Different infill properties, as summarized in Table 1, are also selected for the pre-code and special-code frames. For the former, the average mechanical properties assigned to the infills are representative of those adopted in Italy, and some regions in southern Europe, in frames designed according to old seismic codes (Liberatore and Mollaioli 2015). On the other hand, the (average) mechanical properties of the infills used for the special-code frames are based on the study of Noh Mohammad et al. (2017), which implemented the same properties in investigating the infill response of a highly reinforced and ductile experimental frames. It can be noticed here that the thickness of the infills in the SI frame is less than that in the PI frame. In fact, the SI frame is representative of modern RC frames with weak infills according to the classification of Hak et al. (2013). Such weak infills have been increasingly used and they are more recommended to ensure better seismic performance since they minimize the effects of frame-infill interaction on structural elements. Nevertheless, thicker infills are also used in modern frames, but they are not explicitly investigated in the current study.

Table 1 Properties of the infill walls for the pre-code and special-code frames

3.2 Numerical modelling

The response of the case-study structures is simulated by 2D numerical models developed via OpenSees (Mazzoni et al. 2009). The nonlinear behavior of the structural components is modelled through a lumped plasticity approach for both beams and columns using zero-length rotational springs. Shear springs are also added in series to the rotational ones in the pre-code frames only to account for possible shear failure modes. Both beam-column joints and floor diaphragms are modelled as rigid. Gravity loads and the masses are respectively uniformly distributed on the beams and concentrated at the beam-column intersections. Damping sources other than the hysteretic energy dissipation are modelled through the Rayleigh damping matrix as per Zareian and Medina (2010), using a 5% damping ratio for the first two vibration modes.

The Ibarra–Medina–Krawinkler hysteretic model (Ibarra et al. 2005) is implemented in OpenSees to describe the M-θ hysteretic behavior of the rotational springs. Amongst others, this M-θ hysteretic behavior includes both strength and stiffness degradation during the loading and unloading phases which is essential to correctly evaluate the seismic response during ground-motion sequences. The main parameters of the M-θ backbone curve, as illustrated in Fig. 4a, are the effective initial stiffness K_e, the yielding bending moment M_y, the yielding rotation θ_y, the hardening stiffness K_s, the maximum bending moment M_c, the rotation at the onset of capping θ_cap, the capping negative stiffness K_c and the post-capping rotation θ_pc. The cyclic deterioration is quantified based on energy dissipation using the normalized hysteretic energy dissipation parameter λ. The values of M_y and θ_y are defined according to Panagiotakos and Fardis (2001) while the other parameters are defined according to Haselton et al. (2016). The former study is based on the results of 1,000 tests, mainly cyclic, for different types of RC members, whilst the latter provides formulations based on cyclic and monotonic tests of 255 RC columns failing in flexural or in flexural-shear modes.

The Setzler and Sezen (2008) model is implemented in OpenSees to define the force-deformation relationship of the shear springs as shown in Fig. 4b. It is characterized by the maximum shear strength V_n calculated according to Sezen and Moehle (2004), shear deformation at the onset peak shear strength Δ_v,n, the shear deformation at the beginning of shear failure Δ_v,u and the shear deformation at the axial load failure Δ_v,f.

Fig. 4

Backbone curve of the a moment rotation relationship and b force-deformation for shear

Masonry infills are modelled as equivalent struts and the force-deformation relationship (H–u), i.e., backbone curve, developed by Liberatore and Decanini (2011) is implemented. The main parameters of this backbone curve, as illustrated in Fig. 5a, are the uncracked infill stiffness K_mf, the strength at the first crack H_mf, the displacement at the first crack u_f, the stiffness at complete cracking K_mfc, the maximum strength H_mfc, the full crack displacement u_fc, the residual strength H_mr and the residual displacement u_r. According to Liberatore and Decanini (2011), the parameter H_mfc was considered the minimum strength of four possible failure modes, which are the diagonal tension, sliding shear, corner crushing and diagonal compression. H_mf and H_mr were assumed 80% and 5% of H_mfc respectively. These parameters are used as input to the ‘Hysteretic’ material available in OpenSees. This material is assigned then to the equivalent struts simulating the nonlinear behavior of infills. The parameters of the ‘Hysteretic’ material that determine the shape of the hysteretic cycles, degradation and pinching effect are adopted from Noh Mohammad et al. (2017), who calibrated them experimentally.

Fig. 5

Modelling of the infills equivalent struts: a force-deformation relationship; b modelling approach for special-code frame; c modelling approach for pre-code frame

Single diagonal struts (Fig. 5b) connecting the nodes at the beam-column intersections are modelled in the special-code frames since their columns are characterized by high shear strength and ductility (i.e., capacity design, over-strength factors, high confinement) and, therefore, they are not expected to experience shear failure even if shear demands are increased due to localized interaction with the infills.

Conversely, infills in the pre-code frames are modelled using the double strut approach (Fig. 5c) suggested by Burton and Deierlein (2014). In this case, one strut is diagonal and it connects the nodes at the beam-column intersections while the other is an off-diagonal strut connecting the shear springs in columns. According to Burton and Deierlein (2014), 75% of the total strut strength and stiffness is assigned to the diagonal strut, whilst 25% is assigned to the off-diagonal one. Although this modelling strategy does not consider the full distribution of column shear increase at the contact with the infill, it allows capturing the increase in the columns shear demand due to the infill-frame interaction, and thus accounting for possible changes in the plastic mechanism of the overall frame.

3.3 Definition of DSs and mapping their thresholds by pushover analysis

Structure-specific MIDR thresholds for three different DSs are calibrated via pushover analyses by assessing multiple measurable criteria, which are based on Rossetto et al. (2016) and Eurocode 8 Part 3 (EC8-3) (EN 1998-3 2005) as described in Table 2. The parameters θ_y and θ_u in Table 2 (see Fig. 4) denote the yield and the ultimate rotations calculated according to Panagiotakos and Fardis (2001) and EC8-3 respectively. It should be noted that the MIDR threshold for each DS is mapped based on the first occurrence of any of the selected criteria.

Table 2 Criteria used for mapping of DSs thresholds

Pushover analyses are performed for the four case-study frames to derive the MIDR DSs thresholds based on the damage criteria defined in Table 2. The pushover incremental load patterns are defined according to the first mode shape as indicated in the EC8. Accordingly, Fig. 6 shows IDR vs. storey shear curves of the first storey for all case-study frames. Despite the special-code frames showing a better distribution of the inelastic demands among the stories, the first storey is the one that experiences the largest drift demands compared to the upper ones in all frames. The mapping of MIDR thresholds of the DSs is entirely governed by the drift demands associated with this floor.

Fig. 6

Pushover curves and DS thresholds for the a special-code frames; b pre-code frames

Figure 6a shows that the SB frame undergoes a ductile nonlinear behavior without a significant loss in lateral strength even at high levels of storey drifts. On the other hand, the presence of infills provide the SI frame with a significant increase in its lateral stiffness and strength. However, such increase in strength and stiffness drops down quickly as the infills undergo high levels of deformation and experience damage until the pushover curve of the SI frame matches the one of the SB frame. Modern seismic code provisions also provide the frames with regularity in elevation, which allows the SB and SI frames to undergo uniform nonlinear drifts at different storeys (not illustrated here for brevity) to avoid deformation demand concentration. A global plastic mechanism is registered for both the SB and SI frames.

Conversely, the PB frame shows much lower ductility as indicated in Fig. 6b compared to the special-code frames. It also starts losing lateral strength rapidly after reaching DS2. The contribution of infills in the PI frame is also obvious with respect to increasing the frame strength and stiffness. Moreover, the PI frame reaches DS3 even before the complete damage of the infills. It is important to highlight that the PB frame experiences a significant concentration of drift demands in the first storey, indicating a soft storey mechanism. A similar mechanism is recorded in the PI frame, indicating that there are no detrimental effects due to the local infill-frame interaction, i.e., no column shear failure is observed.

Additional relevant information on the dynamic behavior of the case-study buildings such as, the first-mode vibration periods and the corresponding mass participation ratios are summarized in Table 3. The MIDR thresholds for the DSs are reported in Fig. 6 and summarized in Table 4. As expected, the special-code (i.e., SB and SI) have significantly higher DS thresholds with respect to the pre-code frames (i.e., PB and PI).

Table 3 First-mode vibration periods and corresponding mass participation ratios

Table 4 MIDR thresholds of DSs

3.4 Average spectral acceleration (avgSa)

The avgSa, selected here as IM, is conventionally calculated by considering a range of 10 equally-spaced periods spanning approximately from a lower bound of $0.2{\text{T}}_{{1}}$ and an upper bound of $1.5\text{T}_{{1}}$, where ${\text{T}}_{{1}}$ is the fundamental period of the structure (e.g., Kohrangi et al. 2016). Accordingly, the four considered case studies should be characterized by different definitions of the IM but this would not allow a consistent comparison between the fragility curves. In order to overcome this issue, the same range of periods is considered for all case-study frames using the ${\text{T}}_{{1}}$ of the PB frame only (i.e., 0.17–1.38 s). This range of periods is deemed appropriate for the SB frame as well since special-code frames are ductile and can undergo as large as ${{3}\text{T}}_{{1}}$ of period elongation due to high levels of nonlinearity at high seismic intensities (Baker and Cornell 2006; ASCE 7–16). Furthermore, the 0.17–1.38 s period range is also deemed appropriate for the infilled frames (i.e., SI and PI) for two reasons: (1) the higher modes are practically negligible (Table 3); (2) their expected period elongation can reach 3–4${\text{T}}_{{1}}$ (O’Reilly et al. 2018), which is close enough to the adopted upper bound (1.38 s).

3.5 Implemented damage-to-loss ratios

DLRs are usually estimated empirically through post-earthquake reconnaissance or by the means of expert judgment. These ratios are region-specific and building-type-specific and must be carefully selected while developing vulnerability relationships (e.g., Rossetto and Elnashai 2003). For this study, the DLRs suggested by Di Pasquale et al. (2005), which are implemented in seismic risk analysis of Italian buildings are considered. Table 5 shows the DLRs for all the DSs used in this study. It should be noted that the number of DSs used by Di Pasquale et al. (2005) is larger than the number of DSs used in this study. Therefore, to be able to couple DLRs with fragility curves, different DSs are combined as indicated by Table 5. It is worth mentioning that such DLRs were defined under the assumption that the repair process reinstates the original conditions of buildings before the occurrence of seismic events. No considerations are made to account for possible increase in the replacement cost due to retrofitting and upgrading of the building seismic resistance. Moreover, the DLRs used here are only mean values, so the uncertainty around them is neglected in order to derive mean vulnerability relationships. The DLRs derived by Martins et al. (2016) could also be implemented in the present study. Those DLRs are in good agreement with those by Di Pasquale et al. (2005). It is worth highlighting that the DLR for each DS considered in this study is assumed to be the same for both infilled and bare frame configurations and that other strategies, e.g., component-by-component loss assessment procedures, could be implemented to provide better estimates of the contribution of the infill panels to losses at low DSs. However, a detailed definition of DLRs is outside the scope of the current study and the use of global building-level DLRs is a widely-used and recognized procedure, especially for risk assessment of building portfolios.

Table 5 Damage-to-loss ratios for the case-study structures

4 Results and discussion

4.1 Performance assessment of undamaged structures

The case-study frames are initially subjected to the first ground motion of all sequences (i.e., GM1) to assess their performance in the undamaged state. Figure 7 shows the results of NLTHA in the form of avgSa_GM1 versus MIDR clouds, in addition to the PSDM fitted as a line (in the log-log space) via a least-squares approach. The DSs are also shown for each analysis case based on the MIDR thresholds defined in Sect. 3.3. The special-code frames, as expected, show an overall satisfactory performance as reported in Fig. 7a, b, which can be summarized as follows:

• Both configurations (bare and infilled) remain undamaged (ND) in most analysis cases, particularly in 58.5% for the SB frame and 50.3% of analysis cases for the SI frame;

• The SB frame experienced DS1 in 32% of the cases, while the combined DS0 and DS1 cases for the SI frame are approximately 42% of the total analysis cases;

• Only a small portion of cases resulted in high DSs (i.e., DS2 and DS3). Specifically, 4.5% of the cases represent DS2 in the SB frame and around 3.5% led to DS2 in the SI frame;

• The near-collapse (i.e., DS3) is observed in less than 4% of the cases for both SB and SI frames. It is also noticed that the presence of the infills reduced the number of cases with high DSs compared to SB frame.

• For both special-code frames, the collapse cases are observed in less than 1% of the analysis cases.

Fig. 7

Cloud analysis results for GM1 only for: a special-code bare (SB); b special-code infilled (SI); c pre-code bare (PB); d pre-code infilled (PI) frames

In contrast, the pre-code frames show a significantly weaker performance and experience comparatively-higher DSs in much more analysis cases, as shown in Fig. 7c, d and summarized in the following points:

• The pre-code frames remained undamaged in only 19% and 37% for the PB and PI frames respectively;

• The PB and PI frames reached either DS1 or DS0 in 47% and 49.5% of the analysis cases, respectively;

• The PB frame entered the DS2 range in 17% of the analysis cases, whereas the PI frame experienced the same DS in only 6.5% of the analysis cases;

• DS3 is reached in approximately 9.2% and 5.8% of the analysis cases for the PB and PI frames, respectively. Similarly to the SI frame, an overall positive contribution of the infills in reducing the DSs can be observed.

• The collapse cases represent 7.8% and 1.2% of the analysis cases for the PB and PI frames respectively.

4.2 Fragility curves for the undamaged structures

Figure 8 shows the fragility curves derived for the undamaged case-study structures. The median IM values (μ_DS) and the dispersion ($\beta$) for each DS are given in Table 8. As expected, the fragility of the pre-code frames (Fig. 8c, d) is higher than the special-code ones (Fig. 8a, b). The significant (and beneficial) effect of the infills can also be noticed from the figure/table, even if the infills develop cracking at low IM levels. Specifically, the μ_DS values of the fragility curves for DS1, DS2 and DS3 for the SI frame are respectively 83%, 37% and 27% higher than those for the SB frame. The effect of infills on the pre-code frame is even more prominent since the μ_DS for DS1, DS2 and DS3 for the PI frame are respectively 133%, 104%, and 77% higher than those for the PB frame. This can be attributed to the fact that, being thicker, the infills in the PI frame show higher stiffness and strength when compared to those used in the SI frame.

Fig. 8

Undamaged state fragility curves for a special-code bare (SB); b special-code infilled (SI); c pre-code bare (PB); d pre-code infilled (PI) frames

It is interesting to note that the fragility curves for DS1 and DS2 are more separated from each other in the infilled frames compared to the bare frames. This means that the IM levels needed to cause DS2 are much higher than those causing DS1 in the infilled frames, unlike the bare ones. For example, the difference in μ_DS between DS1 and DS2 is equal to 1.08 g for the SI frame, but it is around 0.90 g for the SB frame. Such observation is more prominent in the PI frame as the difference in μ_DS is equal to 0.88 g, whereas is only 0.45 g for the PB frame. This can be explained as follows: when the infilled frames are in DS1, the infills are still significantly contributing to the lateral strength and energy dissipation even if they have some initial cracks and damage at that stage. This observation is more dominant in the PI frame compared to the SI one due to the different properties and behavior of infills in the two frames. In fact, as shown by the pushover curves in Fig. 6, it is noticed that when the SI frame achieves DS1, the infills have already reached their peak strength and they begin degradation after that as the lateral displacement increases. Conversely, when the PI is in DS1, the infills are still undergoing hardening and they have not reached their peak strength yet.

4.3 State-dependent probabilistic seismic demand models

The considered frames are now subjected to the full ground-motion sequences to investigate their performance under sequential excitations. Even though the case-study structures might remain elastic without showing any damage for a few weak sequences, damage accumulation is expected to take place in the vast majority of the considered sequences, and the adopted nonlinear modelling strategy allows to explicitly capture such an effect. However, since the MIDR can be only used to indirectly capture damage accumulation, it results in some limitations. Indeed, many analysis cases show that the MIDRs resulting from GM2 (MIDR_GM2) are lower than those obtained during GM1 (MIDR_GM1), which mainly occurs in sequences with low IM values for GM2. This is not an indicator of a decreased damage level, but rather related to other reasons, including: (1) the IM values corresponding to GM2; (2) the effects of ground-motion polarity; and (3) the change of dynamic properties and structural response of the frames due to the initial damage caused by GM1. Therefore, any MIDR-based inference on the increased damage level is hardly achievable in the analysis cases with MIDR_GM2 < MIDR_GM1. Moreover, incorporating those analysis cases in fitting linear PSDMs can result in crossover between the fragility curves of the undamaged structure and the state-dependent ones for the damaged structures, which provides an erroneous impression that the damaged structures are less fragile than the undamaged ones.

Bilinear PSDMs (e.g., Ramamoorthy et al. 2006; Jeon et al. 2015; Tubaldi et al. 2016; Freddi et al. 2017) can be more appropriately fitted to the IM vs. EDP clouds of the structures having initial DSs to overcome this issue and reflect the physics of the problem under investigation. To shed further light on this, Fig. 9 shows the clouds of the PI frame resulting from GM2 when the frame is initially in DS1 (Fig. 9a) and DS2 (Fig. 9b) respectively, in addition to the PSDMs for both the undamaged and damaged cases. It can be observed that the fitted bilinear PSDMs have two branches partitioning each cloud into two groups of points and each branch line is fitted independently, which allows for different dispersions (e.g., Tubaldi et al. 2016). The two branches intersect at the break point (i.e., IM_o or avgSa_o) that can be defined by minimizing the sum of the square of residuals (i.e., least-square) both before and after the break point (e.g., Jeon et al. 2015). The values of IM_o in Fig. 9a, b are equal to 0.45 g and 0.29 g respectively. The number of analysis cases used in the fitting of the bilinear PSDMs are as follows: 148 for the SB, 75 for the SI, 193 for the PB, and 72 cases for the PI when all these frames are initially in DS1; 72 for the PB and 29 cases for the PI when they are initially in DS2. In theory, fragility curves developed using bilinear PSDMs have two branches (i.e., two lognormal approximations). One branch is derived based on the fragility parameters resulting from the line fitted prior to the breaking point (i.e., IM_o), whilst the other branch relies on the line fitted after IM_o. Both branches are then joined together at IM_o, as shown for instance in Tubaldi et al. (2016). For simplicity, the state-dependent fragility curves in this study are derived using only the line fitted after IM_o (i.e., the steeper) over the whole range of IMs. This simplification is deemed appropriate since the probability of exceeding different DSs is particularly low for IMs lower than IM_o (e.g., Jeon et al. 2015).

Fig. 9

IM versus EDP clouds of the pre-code infilled (PI) frame when it has an initial damage state of a DS1; and b DS2 resulting from GM1

It is interesting to note that in most analysis cases of which the DS due to GM2 is higher than that caused by GM1, the values of avgSa_GM2 are also higher than avgSa_GM1 as shown in Fig. 10a, b for the SB and PI frames respectively. This is an indicator that generally, if the MIDR is used as an EDP, achieving a higher DS in a structure with an initial DS requires an IM level for GM2 that is higher than or at least comparable to that for GM1. On the contrary, GM2 with low IM values are likely to result in MIDR values lower than those from GM1, as the MIDR is not monotonically increasing with the length of ground shaking. Therefore, using bilinear PSDMs in such situations is beneficial to derive appropriate and statistically-consistent fragility curves.

Fig. 10

Values of avgSa_GM1 versus avgSa_GM2 with the final DSs due to GM2 for a special-code bare (SB) frame; b pre-code infilled (PI) frame

4.4 State-dependent fragility curves

State-dependent fragility curves are derived here for the structures with initial DSs due to GM1 by using the analysis results of the entire sequences (i.e., both GM1 and GM2) as shown in Fig. 11. The x-axis here represents the avgSa_GM1 for the fragility curves of the undamaged structures, whilst the avgSa_GM2 should be used for the state-dependent fragility curves. The resulting fragility parameters are also summarized in Tables 7, 8, 9 and 10. The reduction in the μ_DS values of the state-dependent fragility curves with respect to the undamaged state is also shown.

It should be noted that fragility curves conditioned on DS0 for the infilled frames are not derived as this DS is very close to DS1. Moreover, DS0 is only related to the initial cracking of infills and at that stage all the structural elements are still undamaged. In fact, μ_DS values for the fragility curves conditioned on DS0 are less than 8% different from those for the curves conditioned on DS1. In addition to that, the fragility curves of DS3 conditioned to DS2 (i.e., DS3|DS2) for the special-code frames are not defined due to the lack of relevant analysis cases. In summary, the following observations can be highlighted based on Fig. 11:

• It is clearly shown that the presence of initial damage due to GM1 increases the fragility of all the structures;

• The reduction in the μ_DS values of the state-dependent fragility curves compared to the undamaged state is less than 21.6% for all frames when they are in DS1 due to the initial shaking caused by GM1;

• When the initial DS caused by GM1 reaches DS2, the μ_DS values of DS3 (i.e., DS3|DS2) are significantly reduced, i.e. 31.1% for the PB frame and 41.8% for the PI frame compared to the undamaged state. This observation could be only made to the pre-code frames since the DS3|DS2 fragility curves are not derived for the special-code frames as explained earlier.

Fig. 11

State-dependent fragility curves for the a special-code bare (SB) frame; b special-code infilled (SI) frame; c pre-code bare (PB) frame; and d pre-code infilled (PI) frame

The aforementioned results show that, in general, when the case-study frames are subjected to slight damages (or DS1 here) due to the initial ground shaking by GM1, the increase in their fragility against the final shaking (i.e., GM2) is not very considerable, especially for the near-collapse (i.e., DS3) fragility curves. Conversely, if the initial DS due to GM1 leads to a partial or complete plastic mechanisms (DS2), the increase in their fragility against GM2 is very significant. It is worth mentioning that for each set of fragility curves, the $\beta$ value is constant across different DSs as appears in Tables 6, 7, 8,9 and 10. This represents a limitation of cloud analysis as a single power-law model is fitted via the least-square regression over the entire IM versus EDP cloud regardless of the DS, resulting in constant $\beta$ values, i.e., homoscedasticity assumption.

Table 6 Median and standard deviation values of fragility curves for the undamaged state

Table 7 Special-code bare (SB) frame

Table 8 Special-code infilled (SI) frame

Table 9 Pre-code bare (PB) frame

Table 10 Pre-code infilled (PI) frame

4.5 Derivation of vulnerability relationships

Vulnerability relationships are finally derived for the case-study frames accounting first for the undamaged state and then considering the structures with the initial DS due to GM1 (Fig. 12). As explained earlier, these curves are generated by coupling the DLRs defined in Sect. 3.5 with the DS probabilities as per Eq. (1). The x-axis in Fig. 12 represents the avgSa_GM1 when using the vulnerability relationships of the undamaged structures, whereas it represents avgSa_GM2 if the state-dependent vulnerability relationships are used. It should be also noted that the LR values in the state-dependent vulnerability relationships at avgSa_GM2 = 0 represent the DLR values of the DSs caused by GM1.

Fig. 12

Undamaged and state-dependent vulnerability relationships for the a special-code bare frame (SB); b special-code infilled frame (SI); c pre-code bare frame (PB); d pre-code infilled frame (PI)

As expected, the vulnerability relationships for the special-code frames are significantly different from those of the pre-code ones. High differences are also registered between the bare and infilled frames. LR values also vary depending on the level of initial DS, in line with the fragility curves developed earlier. Considering that the case-study structures are in the undamaged state, the following observations can be highlighted regarding their vulnerability:

• LR values in the special-code frames are significantly lower than those for the pre-code ones. For example, the expected LR in the PB frame (Fig. 12c) reaches approximately 1.00 at an avgSa_GM1 of 1.60 g, whilst the LR in the SB frame (Fig. 12a) is only 0.58 at the same level of IM;

• Similarly, when the PI frame (Fig. 12d) reaches an almost 1.00 LR at an avgSa_GM1 value of 3.0 g, the LR in the SI frame (Fig. 12b) is expected to be only 0.71 at the same level of avgSa_GM1;

• The positive contribution of infills in reducing losses is also notable. For instance, the PB frame has a LR of 1.00 at an avgSaGM1 of 1.60 g, whereas the same level of avgSaGM1 causes a LR of 0.60 only in the PI frame. However, this observation cannot be generalized as it is limited to cases in which infills do not cause localized failure due to demand concentration in surrounding columns, which might trigger collapse;

• Similarly, the SB frame achieves an expected LR of 1.00 at an avgSa_GM1 value of 3.23 g, while the LR of the SI is expected to be 0.76 only at the same avgSa_GM1 level.

The initial DSs induced by GM1 have a considerable effect on the expected LRs resulting from GM2. Accordingly, the following observations are summarized:

• When the SB frame has initially DS1 due to GM1, a LR of 1.00 is expected to occur at avgSa_GM2 of 2.11 g (Fig. 12a), whereas the LR is equal to 0.83 only at the same level of IM when the structure is undamaged;

• Similar observations are made for the other frames since the SI, PB and PI frames reach an LR of 1.00 at avgSa_GM2 of 3.22 g, 1.05 g and 2.18 g respectively when they are in DS1 due GM1 (Fig. 12b–d). The same levels of IM correspond to LRs of only 0.76, 0.87 and 0.84 respectively when the structures are undamaged;

• When the pre-code frames are in DS2 due to GM1, very high LRs are attained even at low avgSa_GM2 (Fig. 12c, d). For instance, the PB and PI frames reach a LR of 1.00 at an avgSa_GM2 of 0.58 g and 0.87 g respectively, whilst for the undamaged case, the LRs at the same IMs are only 0.37 and 0.16, which are significantly lower.

The previous results demonstrate that when the case-study frames undergo a slight level of damage (i.e., DS1) due to GM1, the increase in the LR ranges between 15–31% at high levels avgSa_GM2 compared to the undamaged state, while the increase in LR is very slight at low IM levels (i.e.., less than 0.5 g) as shown in Fig. 12. On the contrary, when the initial DS reaches the DS2, the increase in LR is very significant at all levels of avgSa_GM2 as shown by Fig. 12c, d.

State-dependent vulnerability relationships can be also represented by vulnerability surfaces (García de Quevedo Iñarritu 2018), as illustrated in Fig. 13. The IMs of the GM1 and GM2 are plotted on the horizontal axes, while the vertical axis shows the resulting LR. Each vulnerability surface is derived based on the relevant state-dependent vulnerability functions (i.e., the known slices of the surface), which are fixed at their intersection points with the undamaged-state vulnerability function. This intersection represents the avgSa_GM1 required to produce a LR equal to the initial LR in the state-dependent vulnerability functions. The intermediate points are then estimated using linear interpolation of 3D scatter of points. The vulnerability surfaces for all the case-study frames, as reported in Fig. 13, allow the estimation of the final expected LRs due to ground-motion sequences. The proposed surfaces can be particularly useful when implemented in a simulation-based time-dependent risk assessment framework of building portfolios. It should be noted that the vulnerability surfaces in this study assume that the ground-motion sequences are composed of two motions only, i.e. GM1 and GM2 and this represents a limitation of the obtained results. Further research is needed investigating the effects of multiple earthquakes within each sequence (e.g., Panchireddi and Ghosh 2019).

Fig. 13

Vulnerability surfaces for a the special-code bare (SB) frame; b special-code infilled (SI) frame; c the pre-code bare (PB) frame; and d the pre-code infilled (PI) frame

5 Conclusions

This study investigated the vulnerability of case-study RC frames against seismic sequences. Four case-study frames were selected including special-code frames designed according to EC8 and pre-code frames designed to resist gravity loads only. Both bare and uniformly-infilled frame configurations were considered. Artificial sequences composed of two ground motions (i.e., GM1 and GM2) were assembled via a randomized approach and then implemented in the NLTHAs via a sequential cloud-based approach. Fragility and vulnerability relationships were developed for the case-study structures in two states: (1) the undamaged state; (2) the damaged state in which the frames underwent an initial DS due to GM1.

Analysis of the undamaged frames showed that the special-code frames remained undamaged or experienced low DSs in the majority of analysis cases, with very few cases of near-collapse, unlike the pre-code frames that developed high DSs (including near-collapse) in a considerably-larger portion of cases. Since the infill-frame interaction did not cause any local plastic mechanism for the frames, the positive impact of infills in both the pre-code and special-code frames was obvious as it lowered the number of cases with high DSs. The resulting fragility curves confirmed these outcomes: the special-code frames were significantly less fragile than the pre-code frames while the presence of infills lowered the fragility of both categories of frames.

The performance of the damaged frames was then investigated by analyzing the full sequences, which showed that both MIDR values and DSs increased in many analysis cases compared to those obtained from the analysis of the undamaged frames. Such increase is much more prominent in the pre-code frames rather than special-code ones. It is noticed also that the sequences in which the frames experience an increase in the DS during GM2 compared to that acquired in GM1, are mostly the ones with avgSa_GM2 values that are at least comparable or greater than avgSa_GM1.

State-dependent fragility curves were developed for the frames initially damaged by GM1. These curves showed that the presence of initial DSs due to GM1 increased the fragility of the structures when they were subjected to GM2. However, such increase was not major when the initial DS was low (i.e., DS1) as, for example, the μ_DS values were reduced by less than 21.6% compared to those for the undamaged case. Conversely, the fragility increase when the initial DS was large (i.e., DS2) was very significant and, in fact, the μ_DS showed a reduction up to 42% in the pre-code frames.

Finally, the study derived vulnerability relationships to evaluate the expected LRs over a range of IM values. The resulting functions were in line with the fragility curves as the pre-code frames experienced higher LR values than the special-code frames, with a clear positive contribution of the infills in reducing the expected losses in both types of frames. It was also shown that the presence of initial DSs induced by GM1 led to increasing the expected LR values due to GM2. For example, when the frames were in DS1 due to GM1, the increase in the LR ranged from 15–31% at high IM levels compared to the undamaged state. This increase in the LR values was much more significant when the initial DS due to GM1 was larger (i.e., DS2 in the pre-code frames). Vector-valued vulnerability relationships were also generated for the considered frames in the form of a three-dimensional vulnerability surfaces that can be adopted to quantify the final expected seismic losses due to a ground-motion sequence that can be easily adopted in portfolio risk assessment.

Overall, this study provides insights about the importance of considering the effects of ground-motion sequences and damage accumulation on the fragility and vulnerability of the considered case-study frames. Additional research is needed to account for different frame layouts, heights, design levels in terms of material properties and details, and to study other EDPs that are better-correlated with damage accumulation due to sequential ground shaking compared to the MIDR implemented in this study.