Collapse fragility analysis of historical masonry buildings considering in-plane and out-of-plane response of masonry walls.

Engineering Structures 319 (2024) 118804 Available online 24 August 2024 0141-0296/© 2024 The Authors. Published by Elsevier Ltd. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/). Collapse fragility analysis of historical masonry buildings considering in-plane and out-of-plane response of masonry walls. Daniel Caicedo a,* , Igor Tomi´ c b , Shaghayegh Karimzadeh a , Vasco Bernardo a , Katrin Beyer b , Paulo B. Lourenço a a University of Minho, ISISE, ARISE, Department of Civil Engineering, Guimar˜ aes, Portugal b ´ Ecole Polytechnique F´ ed´ erale de Lausanne (EPFL), Earthquake Engineering and Structural Dynamics Laboratory (EESD), Lausanne, Switzerland ARTICLE INFO Keywords: Historical masonry buildings Incremental dynamic analysis Out-of-plane response Real ground motion records Record-to-record variability Fragility curves ABSTRACT This paper explores the seismic fragility assessment of historical masonry buildings in the context of performance-based earthquake engineering, using multiple-record incremental dynamic analyses (MR-IDA). The methodology is applied to two case studies, the first being a stiff monumental structure, and the second a tall and slender masonry building. Both case studies are modelled in the OpenSees software using three-dimensional macroelements that account for the in-plane (IP) and out-of-plane (OOP) response of masonry walls. Simultaneously, the numerical models account for the influence of non-linear connections in floor-to-wall interfaces and wall-to-wall interlocking. The record-to-record variability induced by the random nature of the ground motion is addressed through earthquake selection consistent with the variability of European hazard. Ground motion records are normalised at a uniform level of seismic intensity to be sequentially increased for the computation of IDA curves. Once IDA curves and fragility functions are derived, the failure results are thoroughly discussed. The fragility functions accounting for the uncertainties arising from the record-to-record variability are contrasted against the results of the seismic assessments affected by multiple sources of uncertainty in material and modelling parameters from a previous study. Finally, the fragility curves derived from the methodology presented herein are directly compared with fragility functions from probabilistic seismic demand analysis (PSDA). Overall, IDA-based fragility functions are found to provide more conservative predictions than the ones obtained from the cloud-based approach using unscaled records. 1. Introduction Throughout Europe, historic urban centres are mostly composed of masonry constructions, some of them being landmarks of cultural heritage [1]. Masonry structures are known to be highly vulnerable to earthquake damage. On top of that, a large portion of the European building stock is located in regions of significant seismicity [2,3]. Depending on the age of construction or the structural typology, masonry buildings might exhibit a very different behaviour during major earthquake events, as observed after the Faial (1998) earthquake in Portugal [4], the L’Aquila (2009) and Emilia (2012) earthquakes, both in Italy, the Lesvos (2017) earthquake in Greece [5–7], and more recently the Kahramanmaras¸ Earthquakes (M w 7.7 and M w 7.6) in Turkey [8]. It is paramount to adopt efficient procedures for assessing the seismic vulnerability of historical masonry constructions, which may be a cumbersome task when dealing with epistemic and aleatory uncertainties [9–11]. In this sense, that deterministic assessment of such buildings may lead to unrealistic results given the significant sources of uncertainties in their structural characterization. Conversely, the probabilistic framework of Performance-based earthquake engineering (PBEE) [12,13] emphasizes several issues at the system level, such as collapse risk, repair costs, fatalities, and post-earthquake loss of functionality. Thus, this methodology denotes a useful tool for determining important outcomes and findings in the assessment of historical masonry construction. Besides, expected annual loss and collapse risk have been integrated lately in line with modern PBEE requirements [14]. In the context of PBEE, several methodologies have been implemented to tackle the collapse safety investigation of engineering structures, including singleand multi-record incremental dynamic analysis (IDA) [15], endurance time analysis (ETA) [16], and multiple stripe analysis (MSA) as a particular case of IDA with dynamic analyses performed at discrete levels of intensity [17,18]. * Corresponding author. E-mail address: [email protected] (D. Caicedo). Contents lists available at ScienceDirect Engineering Structures journal homepage: www.elsevier.com/locate/engstruct https://doi.org/10.1016/j.engstruct.2024.118804 Received 9 May 2024; Received in revised form 1 August 2024; Accepted 13 August 2024 Engineering Structures 319 (2024) 118804 2 IDA stands as one of the most efficient approaches for evaluation of the seismic resilience of the structural systems [19,20]. The fundamentals of IDA were rigorously established in [15], and its first practical application was depicted by Vamvatsikos et al. [21] on a medium-height older reinforced concrete (RC) structure. These contributions focused on aspects of structural behaviour revealed by IDA, such as the large record-to-record variability, the hardening and softening, the flatline, and the structural resurrection (further discussion on these concepts is presented in Chapter 2). Later on, Dolsek [22] introduced a set of models of a four-storey RC frame to address modelling uncertainties using the Latin hypercube sampling (LHS) method. It was found that including modelling uncertainty in extended IDA increases the computational cost and may have significant implications on the collapse capacity. Likewise, IDA methodology has been applied to seismic risk assessment of bridges, leading successfully to the identification of damage states (DSs) and derivation of fragility curves [23,24]. For the case of bridges, the PBEE probabilistic framework has also been adopted to estimate repair costs in recovery investments [25]. More recently, Repapis and Zeris [26] contrasted the results of static pushover analyses and IDA procedures for the case of reinforced concrete buildings with clay brick masonry infill walls, while He and Lu [27] focused on seismic fragility assessment of a super tall building with hybrid control using IDA methodology. Kazemi et al. [28,29] employed IDA techniques to study the seismic vulnerability of steel moment-resisting frames (SMRFs) and showed the improvement of performance levels when infill masonry is adopted as a retrofitting strategy for the steel-damaged buildings. Other engineering structures have been examined in the context of PBEE through IDA such as concrete dams and isolated nuclear power plants [30,31]. The literature survey reveals that there is a scarcity of application of IDA techniques for the fragility assessment of masonry buildings, mainly attributed to the extensive computational demand in large-scale modelling of them [32,33]. For instance, Abo-El-Ezz et al. [34] approached the fragility assessment of low-rise stone masonry buildings through a displacement-based damage model accounting for the characteristics of existing stone masonry buildings in Eastern Canada. Vanin and Beyer [35] proposed a logic tree approach to incorporate the epistemic uncertainty in the seismic assessment of masonry building and pointed out the displacement capacity as a major source of bias compared to record-to-record variability. Chieffo et al. [36] estimated the vulnerability of an isolated masonry building in terms of expected damage using macroelements and non-linear static analyses. The variability of structural and geometrical parameters in the fragility assessment of masonry building aggregates, through non-linear static analyses, was accounted for in the works of Battaglia et al. [37] and Angiolilli et al. [38]. Thanks to the advances in the macroelement modelling strategy introduced in [39], Tomi´ c et al. [40] successfully accounted for the randomness in material properties of historical masonry buildings within the single-record IDA framework. Recently, Bernardo et al. [41] provided analytical fragility curves, supported by non-linear static analysis, for representative building archetypes accounting for the uncertainty in materials properties, geometry, and seismic demand for several levels of seismicity in Portugal. While the focus of these investigations has been on uncertainty scenarios arising mainly from geometric and material properties the effects of ground motion variability have been left aside, and there is a clear need for their investigation. This paper adapts the multi-record IDA procedure to evaluate the seismic fragility of historical masonry buildings in the context of PBEE. The proposed framework is applied to two historical masonry constructions: (1) the Holsteiner Hof building, as an example of stiff monumental heritage structures; and (2) the Lausanne Malley, representative of tall and slender residential masonry buildings. The OpenSees software [42] is adopted for the numerical modelling of the case studies, in which three-dimensional macroelements, proposed recently by Vanin et al. [39], are used to account for IP and OOP effects of masonry walls. Moreover, the effects of non-linear connections in floor-to-wall interfaces and wall-to-wall interlocking are considered by implementing simplified frictional and tension-damage one-dimensional constitutive laws, respectively [43]. A set of accelerograms consistent with the European hazard is selected to cover the record-to-record variability, and non-linear time history analyses are conducted subsequently. Within this process, ground motion records are normalised first at a uniform level of intensity and scaled up sequentially for the computation of IDA curves. The fragility curves derived from this process, which account exclusively for uncertainty stemming from the record-to-record variability, are contrasted against the results of the fragility analysis for the same buildings accounting for modelling uncertainties through single-record IDA methodology [40]. In addition, a final comparison between the fragility functions computed by means of the framework proposed herein and the ones obtained from probabilistic seismic demand analysis (PSDA) [44,45] is presented. 2. Theoretical background The main goal of this investigation is to implement multiple-record IDA for the collapse fragility assessment of historical masonry buildings in the context of PBEE. Hence, it is needed to predict the structural performance of the buildings under analysis for a full range of credible earthquake ground motion intensities [46]. Through probabilistic seismic response analysis, one can determine the mean annual rate of exceedance of a structural response parameter, regarding and engineering demand parameter (EDP) exceeding certain level of edp [47]. Mathematically, this process can be expressed by means of the total probability theorem as: v(EDP >edp) = ∫∞ 0 (1−P[EDP <edp|IM])dv(IM) dIM dIM (1) where P[EDP <edp|IM] denotes the probability that the EDP is smaller than a certain level of edp at the ground motion intensity measure (IM), and v(IM) is the mean annual rate of exceedance of the ground motion intensity IM. The following subchapters explain the concepts of fragility functions as a mean to quantify the probability of reaching a particular level of DS—or collapse ultimately—at a certain level of IM, as well as MR-IDA as the adopted procedure for its computation. 2.1. Analytical collapse fragility functions Conceptually, fragility functions were introduced in the field of earthquake engineering by Kennedy et al. [48]. Fragility curves state the conditional probability of a structural system reaching or exceeding a DS or a predefined level of EDP given an IM parameter [49]. Mainly, these functions can be grouped into three types, empirical, analytical, and heuristic [50]. Empirical fragility curves are reproduced after earthquake damage observations [51,52] while heuristic functions are based on experts’ opinions [53]. Conversely, the analytical methodology requires structural modelling and numerical simulations since they are constructed based on non-linear dynamic analysis [54–57]. The analytical collapse fragility function is defined in Eq. (2). P[C|IM =x] = Φ(ln(x/ η ) β)(2) where P [C |IM =x] denotes the probability of a structural system to collapse under a ground motion with an intensity level x, Φ is the standard normal cumulative distribution function, η is the median of the fragility function, and β is the logarithmic standard deviation, also referred to as dispersion. Eq. (2) assumes that the values of IM of ground motions producing the collapse of the structural system follow a log-normal distribution. Previous investigations have confirmed the validity of this assumption D. Caicedo et al. Engineering Structures 319 (2024) 118804 3 [49,58–61]; nonetheless, alternative distributions can be implemented. Furthermore, the implementation of Eq. (2) requires the calibration of η and β from the results of structural analysis. The fitting of a fragility function to the analytical set of data points leads then to the estimation of  η and β. In this regard, Baker [18] exposed common statistical approaches for estimating these parameters from analytical data. Among them: (i) the method of moments (MM) in which the parameters  η and β are found such that the resulting distribution has the same moments (median and standard deviation) as the data points; (ii) the sum square error (SSE) that minimizes the error between the observed fractions of collapse data and probabilities observed within the fragility function; and (iii) the maximum likelihood estimation (MLE) in which the  η and βparameters in the resulting distribution have the highest likelihood of having produced the observed data. Prior to fitting, percentile IDA curves should be used to derive fragility curves depicted by single dots (i.e., each dot representing a ratio of the number of records that caused structural collapse to the total number of records at an IM level). Former investigations have addressed the fragility assessment of unreinforced masonry buildings by measuring structural damage in terms of drift ratio [37,38,62,63]. The values of drift limits reported in these investigations will serve as a reference to establish the asymptotic values of EDPs for deriving analytical fragility functions at collapse DS. 2.2. Multiple-record incremental dynamic analysis MR-IDA IDA methodology implies the sequential increment of intensity for each ground motion in a suite until they cause numerical collapse in the non-linear dynamic analysis, which is interpreted as structural collapse [15]. Overall, the methodology can be summarised as: (i) Develop a robust numerical model for the structural system under analysis able to capture potential sources of non-linearity (e.g., concentrated/distributed plasticity, second-order effects, etc.); (ii) Select an adequate suite of records to capture the ground motion variability [47]; (iii) For a given IM parameter and EDP, plot each single-record IDA (capacity curve) including the median and the 16 % and 84 % quantiles of the distribution, to study the seismic demand from low seismic intensity levels prior to yielding up to ultimate performance levels; (iv) Use the IDA data to analyse deeply and better understand the behaviour of the structure (e.g., identification of limit-state capacities, collapse fragility analysis, etc.). It should be noted that because of the complex structural behaviour frequently observed in masonry structures, influenced by multiple sources of non-linearity (i.e., materials, second-order effects, floor-towall, and wall-to-wall connections, etc.), and large record-to-record variability, IDA curves generated by different earthquake motions often result in quite dissimilar responses that are difficult to predict a priori [15]. Thus, fractile curves including the median, 16 % and 84 % percentiles, are a suitable approach for summarising IDA data. These percentile curves are much smoother than individual IDA curves and can better represent the global behaviour of a structure [23]. Besides, limit states can be easily identified on these curves. Immediate occupancy (IO) is violated when the building exceeds a certain level of EDP still within the elastic range; collapse prevention (CP) is reached when the local slope on the IDA curve decays to 20 % of the initial elastic slope; and finally, the global dynamic collapse (GC) is clearly observed at the flatline, where the structure responds with practically infinite EDP to any IM increase [64]. Depending on the particular dynamic behaviour of the structural system, IDA curves may exhibit a softening response; minor hardening; severe hardening; and also, a weaving behaviour (see Fig. 1a). Additionally, Fig. 1b depicts the resurrection phenomenon which seems to be a premature failure at a given seismic intensity level, but a safe response at a higher one [30]. The hunt & fill algorithm [15] has been used to optimally select the scaling levels of IMs and minimize the number of runs by finding prematurely the level of GC and performing additional analyses at intermediate levels. While this strategy has been successful in different applications [65–67], the capacity curves of masonry buildings are expected to exhibit discontinuous behaviour marked by the incidence of structural resurrection, which turns the early identification of GC into a more cumbersome task. Thus, the hunt & fill algorithm is not used, and discrete levels of IMs are adopted to characterise the overall structural behaviour of the case studies from elasticity to yielding and finally collapse. 3. Numerical modelling of buildings and ground motion characterisation As mentioned, IDA implies performing a series of non-linear dynamic analyses of a structure, using multiple records scaled up to higher levels of intensity, to cover the overall model’s behaviour, from elastic to yielding and inelastic incursion, finally leading to global dynamic Fig. 1. IDA curves typologies Source Vamvatsikos and Cornell [15]. D. Caicedo et al. Engineering Structures 319 (2024) 118804 4 collapse. Therefore, the structural models for IDA should be capable of simulating inelastic deformation, force redistribution due to sequentially deterioration in the inelastic range, contribution of higher modes, and structural collapse ultimately. Below, the description of the case studies and numerical modelling approach are presented, as well as the record selection strategy. 3.1. Numerical modelling of buildings Two buildings in Switzerland, representative of typical stiff monumental buildings, and tall and slender masonry buildings, are selected as case studies. The first one corresponds to the Holsteiner Hof building located in the city centre of Basel and considered a cultural heritage landmark according to the Swiss Inventory of Cultural Property of National and Regional Significance. The second case study resembles a residential URM building located in the city of Lausanne. Table 1 provides a general description of the buildings’ topology regarding the number of storeys, dimensions in plan, storey height, wall thickness, floor and roof system. and retrofitting. A more detailed description of the building can be found in [40,68]. Both structures are modelled in OpenSees [42] through the equivalent frame model (EFM) approach using a recently developed three-dimensional macroelement [39] to account for IP and OOP effects simultaneously. The macroelement is formulated as a one-dimensional element with two nodes at the element ends and one additional node at the midspan through which the non-linear shear/flexural response is accounted for [70]. The element can capture the IP and OOP response through three sectional models applied at the element ends and at the central section (P-Δ formulation is considered to capture the non-linear geometrical effects within compatibility equations). Drift values can be calculated individually for flexural and shear deformations by considering the rotations and lumped shear deformations at the central node. Exceeding the limits in drift values (δ c,flexure or δ c,shear ) will lead to the loss of lateral strength of the element. The floor system is modelled using orthotropic elastic membranes with higher stiffness in the direction of the beam span and lower stiffness in the other direction (i.e., membrane definition is given by the two moduli of elasticity in the orthogonal directions, shear modulus, and thickness of the diaphragm). Nonetheless, floor-to-wall connections are modelled to account for non-linear behaviour and potential connection failure at the beam support, which can result in the OOP failure of a pier element. In this regard, zero-length elements are adopted for modelling the floor-to-wall connections and a simple material model for frictional interfaces is defined. The local direction “x” of this material is given by the direction in which the vertical load is acting. In the perpendicular plane, frictional sliding is allowed ruled by a friction coefficient, μ , and no cohesion. The material can fail and lead to the loss of load-bearing capacity if a maximum slip is exceeded in the positive local “y” direction, defined by the orientation of the beam. Furthermore, the material can model pounding of the beam when the slip in the negative direction “y” or any direction “z” exceeds a predefined “gap” value. Likewise, wall-to-wall interfaces are modelled through zero-length elements to which a proper uniaxial material model is attached to simulate the formation of vertical cracks and separation of the orthogonal walls, which might lead to potential OOP failure of macroelements. The material features a linear elastic behaviour in compression, with no crushing, and a finite tensile strength with exponential softening. Fig. 2 describes the simplified frictional and tension-damage one-dimensional constitutive laws adopted for the modelling of floor-to-wall interfaces and wall-to-wall interlocking, respectively. To perform the non-linear dynamic simulations, a 5 % proportional Rayleigh damping is assumed, and a secant stiffness damping model is adopted to avoid overdamping in the OOP failure mechanism, and to capture the full OOP rocking response. Additionally, it should be mentioned that the dynamic simulations conducted herein do not take into account the role of soil deformability or non-linear soil-structure interaction (SSI). Table 2 summarises the modelling parameters adopted as the mean or median values reported by Tomi´ c et al. [40]. The symbol (*) over the values in the third column of Table 2 denotes the median value taken from a lognormal distribution. The median values of 1 reported for k floor , f w , and μ f-w , were implemented in [40] to vary the values of the floor stiffness, the stiffness of the interlocking interface, and the friction coefficient that governs the frictional sliding (see Fig. 2a). Nonetheless, their influence is neglected in this research since only the effect of record-to-record variability is investigated as main source of uncertainty. Fig. 3 illustrates the numerical models of the Holsteiner Hof and Lausanne Malley buildings, respectively, as well as their first three vibration periods, T 1 , T 2 , and T 3 . 3.2. Ground motion characterisation One of the main aspects of IDA and collapse capacity estimation lies in the utilisation of a proper record selection strategy for the inputs of inelastic time history analysis [23]. For Eq. (1) to be valid, the seismic inputs should be selected in such a way that the main seismological features (e.g., Moment magnitude, M w ; distance metrics, R JB ; shear wave velocity, V s30 ; and site conditions; etc.) affecting the levels of IMs, are adequately represented in the selected suite of motions. To this end, a large database of accelerograms is collected covering the hazard of the most relevant seismic-prone areas in Europe (i.e., Italy, Greece, Turkey, Portugal, etc.) [79–82]. A methodology based on unconditional selection, not dependent on structural periods, by Jayaram et al. [83] is implemented for the definition of the suite of motions, composed of 21 recordings whose response spectra target the mean and variance pre-defined by the confidence interval of the European dataset (i.e., means and covariances are determined based on the predefined target spectrum). Within the selection process, the seismological parameters are set as: 4.5 ≤M w ≤7.8; 90 m/s ≤V s30 ≤1050 m/s; R JB ≤185 km; and no pulse-like records [84,85]. Table 3 shows the 21 selected motions and provides their description regarding the name of the earthquake, event ID, latitude and longitude of the epicentre (Lat. epi and Long. epi), depth of the focus defined by the focal depth (FD), M w , R JB , and V s30 . Lastly, Fig. 4 shows the 5 % damped geometric mean spectral acceleration (S a ) of the selected records, including the mean, median, and 95 % confidence interval. Table 3 shows different soil conditions, mainly soil types C and D according to Eurocode 8 classification [86]. It should be remarked that the methodology adopted herein for ground motion selection does not target a specific soil type but is based on the mean target defined from the collected dataset, where the predominance of recorded accelerograms in stations with soil types C and D is notorious. Thus, the effect of soil conditions on fragility analyses of historical masonry constructions remains a pending issue and requires further investigation. 4. Analysis of results This section presents the discussion of IDA results in regard to the Table 1 Description of the buildings’ topology. Building Holsteiner Hof Lausanne Malley General Stiff monumental heritage structure Tall and slender residential masonry building Nº storeys 2 6 Dimensions in plan 26.00 m ×14.00 m 14.00 m ×12.00 m Storey height 4.50 m 2.80 −3.20 m Wall thickness 60 −30 cm 60 −25 cm Floor system Timber beams and planks Timber beams and planks Roof system Secondary wooden truss structure Secondary wooden truss structure Retrofit Minor during 1976–1979. Soundproofing and seismic performance[69]. D. Caicedo et al. Engineering Structures 319 (2024) 118804 5 Holsteiner Hof and Lausanne Malley buildings. To accomplish this task, ground motion records are normalised at a uniform level of intensity, and non-linear dynamic analyses are conducted increasing gradually the level of intensity at each new analysis. The large number of analyses provides ample collapse data which is examined thoroughly to characterise the phenomena. Subsequently, capacity curves from single IDA are derived and summarised using smooth percentiles (16 %, 50 %, and 84 %, respectively) in which performance points are identified to study the behaviour of structures throughout different DSs. The results computed through the methodology proposed herein account exclusively for the effects of record-to-record aleatory randomness while material properties are assumed to be deterministic (see Table 2). However, the investigation by Tomi´ c et al. [40] dig into the effects of the uncertainty stemming from limited data in material and modelling properties by reproducing 400 combinations of 11 modelling parameters through the Latin hypercube sampling (LHS) [87]. Hence, the results of collapse fragility analyses accounting independently for both sources of uncertainty are compared. A final comparison between IDA and PSDA approaches for historical masonry buildings is presented Fig. 2. One-dimensional constitutive laws to model non-linear connections [43]. Table 2 Masonry and modelling parameters. Parameter Definition Mean value Masonry Parameters E m [Pa] Modulus of elasticity[71–73] 3.5 ×10 9 G m [Pa] Shear modulus[71–73] 1.5 ×10 9 * f’ cm [Pa] Compressive strength[71–73] 1.3 ×10 6 c m [Pa] Cohesion[71–73] 2.33 ×10 5 * μ m [-] Friction coefficient[71–73] 0.25 * Modelling Parameters k floor [-] Floor stiffness factor[74,75] 1 * f w [-] Wall-to-wall connection factor[76] 1 * μ f-w [-] Floor-to-wall friction coefficient[43,77] 1 * δ c,flexure [-] Drift capacity in flexure[78] 1.04 %* δ c,shear [-] Drift capacity in shear[78] 0.70 %* ζ [-] Damping ratio[43] 5 % Fig. 3. Numerical models and fundamental vibration periods. D. Caicedo et al. Engineering Structures 319 (2024) 118804 6 using the results reported in Caicedo et al. [68]. 4.1. Holsteiner Hof A total of 399 analyses were conducted with intensity levels distributed as [0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40 0.45 0.50 0.55 0.60 0.65 0.70 0.80 0.90 1.00 1.10 1.20] out of which 271 led the building to global collapse; all of them at intensity levels higher than 0.20 g of peak ground acceleration (PGA). Collapse is observed mainly in two scenarios: (i) when P −Δ effect triggers the loss of the global equilibrium after reaching excessive levels of OOP deformation, and (ii) induced by a series of IP failures which makes it impossible to reach convergence and local equilibrium. Moreover, when the macroelement reaches either δ c,flexure or δ c,shear limits, the lateral stiffness and strength are set to zero, but the pier retains the ability to transfer axial load. When 50 % of the piers in the same direction of one storey reach their drift limits, the IP failure is labelled either as flexure or shear. Moreover, mixed IP failure can be identified when approximately the same amount of piers exceed the δ c,flexure and δ c,shear limits. Analogously, mixed IP-OOP failure may be identified as well. On this occasion, 44 % of the failure observations correspond to OOP, while the remaining 56 % are IP failures. Within the IP observations, flexure IP is predominant with 58 % of occurrences, followed by 25 % of mixed IP failures, and merely 17 % by pure shear IP (88, 38, and 26 collapse cases, respectively). Approximately, 60 % of the IP failures are evenly distributed between the second floor 2 in the “x” and “y” directions (F2 x-dir and F2 y-dir). 14 % of IP are mixed located while the remaining 29 % is distributed among first floor “x” and “y” directions and the third one in “x” direction (F1 x-dir, F1 ydir, and F3 x-dir). 54 % of the OOP occurs at F1 y-dir while 24 % and 19 % occur at F1 x-dir and F2 x-dir, respectively. Only a small portion of the failures is observed at F2 y-dir which correspond exclusively to 3 analyses. In general, all OOP failures correspond to overturning of the main façade and central node overturning. Results of failure analyses for the Holsteiner Hof building are summarised by means of the statistics portrayed in Fig. 5. 4.1.1. IDA curves and fragility functions Capacity curves from IDA are plotted in Fig. 6 showing the maximum values of EDP, adopted as the maximum values of average roof displacement (max(Δr)). The equivalent levels of global drift ratio are shown as well at the secondary axis on top. Because of the high levels of non-linearity, severe hardening, and weaving behaviour are predominant within individual curves. Structural resurrection is also detected—sometimes more than once within the same curve—and it is attributed simply to problems with convergence in the non-linear dynamic analyses. At least 11 individual IDA curves are necessary to accurately estimate the quartile curves. 16th, 50th, and 84th percentiles are useful for interpretation of these capacity curves, and to better capture the full performance of the structure through different DSs, from linear elastic to inelastic behaviour, until collapse at last. For this purpose, performance points are identified. For the Holsteiner Hof building, the limit for IO is specified approximately at the end of the elastic branch at max(Δr) lower than 5 mm or global drift ratios lower or equal to 0.05 %. CP, on the other hand, is recognized when the slope of the IDA curve reduces to 20 % of the initial elastic slope; all CP performance points are identified at drift ratios lower than 0.15 %. Finally, the GC is revealed by the appearance of the flatline at max(Δr) =40 mm, or global drift ratios approaching 0.45 %. It should be noted that these observations are coherent with the findings of previous investigations on the global drift ratio limits of unreinforced masonry buildings [37,38,62, 63] for the derivation of fragility functions at different levels of DS, as stated in Section 2.1. Furthermore, recent investigations [88–90] identified, numerically and experimentally, drift-based DSs for historical masonry structures, showing consistency with the numerical results attained herein. Additionally, Fig. 6 illustrates how the dispersion induced by the record-to-record variability (β R ) influences MR-IDA curves when considering different IMs. The values of β R are reported Table 3 Selected ground motions. Record No. Name Event ID Lat. epi [km] Long. epi [km] FD [km] M w R JB [km] V s30 [m/s] 1 TK.3805. HNX.D.INT−20230206_0000008. ACC.MP.ASC INT−20230206_0000008 37.17 37.08 20.00 7.8 175.45 384.00 2 IT.MSCT.00. HGX.D.EMSC−20161026_0000095. ACC.MP.ASC EMSC−20161026_0000095 42.90 13.09 9.60 5.9 38.80 652.00 3 HI.NAX1. HNX.D.EMSC−20201030_0000082. ACC.MP.ASC EMSC−20201030_0000082 37.91 26.84 10.00 7.0 134.88 461.26 4 IT.SDM.00. HGX.D.EMSC−20160824_0000006. ACC.MP.ASC EMSC−20160824_0000006 42.70 13.23 8.10 6.0 45.23 752.00 5 IT.CPS.00. HGX.D.EMSC−20161030_0000029. ACC.MP.ASC EMSC−20161030_0000029 42.83 13.11 10.00 6.6 64.34 732.00 6 IT.ASG.00. HNX.D.IT−1976 −0002. ACC.MP.ASC IT−1976 −0002 46.26 13.30 5.70 6.4 127.85 960.00 7 IT.CPC.00. HNX.D.IT−2012 −0011. ACC.MP.ASC IT−2012 −0011 44.84 11.07 8.10 6.0 54.31 174.00 8 HL.MSLA.00. HNX.D.GR−1993 −0027. ACC.MP.ASC GR−1993 −0027 38.20 21.77 12.40 5.6 32.93 255.40 9 IT.MSC.00. HGX.D.EMSC−20160824_0000013. ACC.MP.ASC EMSC−20160824_0000013 42.79 13.15 8.00 5.5 31.10 652.00 10 IT.FMG.00. HNX.D.IT−2009 −0009. ACC.MP.ASC IT−2009 −0009 42.34 13.38 8.30 6.1 19.26 790.00 11 IT.SSU.00. HNX.D.IT−2012 −0008. ACC.MP.ASC IT−2012 −0008 44.90 11.26 9.50 6.1 51.02 489.00 12 TK.3146. HNX.D.INT−20230206_0000222. ACC.MP.ASC INT−20230206_0000222 38.11 37.24 10.00 7.5 180.33 641.30 13 IT.TRL.00. HGX.D.EMSC−20170118_0000034. ACC.MP.ASC EMSC−20170118_0000034 42.53 13.28 9.60 5.5 23.88 380.00 14 IT.RTI.00. HGX.D.EMSC−20170118_0000037. ACC.MP.ASC EMSC−20170118_0000037 42.50 13.28 9.40 5.4 35.18 170.00 15 IT.NOR.00. HGX.D.EMSC−20170118_0000119. ACC.MP.ASC EMSC−20170118_0000119 42.47 13.27 9.50 5.0 36.51 423.00 16 HI.VAS2. HNX.D.EMSC−20151120_0000014. ACC.MP.ASC EMSC−20151120_0000014 38.47 20.49 12.00 4.7 19.50 260.70 17 IT.CTL.00. HNX.D.INT−20221109_0000046. ACC.MP.ASC INT−20221109_0000046 43.98 13.32 5.00 5.5 43.88 208.00 18 IT.GBP.00. HNX.D.IT−1998 −0063. ACC.MP.ASC IT−1998 −0063 43.18 12.77 10.00 4.8 20.11 224.00 19 IT.CSC.00. HNX.D.IT−1997 −0004. ACC.MP.ASC IT−1997 −0004 43.02 12.89 5.70 5.7 27.34 698.00 20 IT.AQK.00. HNX.D.IT−2009 −0174. ACC.MP.ASC IT−2009 −0174 42.50 13.38 9.00 5.0 15.41 705.00 21 IT.CLF.00. HGX.D.EMSC−20161103_0000003. ACC.MP.ASC EMSC−20161103_0000003 43.03 13.05 8.10 4.7 9.25 145.00 Fig. 4. 5 % damped geometric mean S a of selected records. D. Caicedo et al. Engineering Structures 319 (2024) 118804 7 in Table A1 for a full list of IMs identified previously in the research by Caicedo et al. [68]. Besides PGA, which defines the level of intensity within the performed IDAs, the acceleration spectrum intensity (ASI) is depicted in Fig. 6b since it shows the lowest value of β R , being outperformed just by the composed metric I comp (See Fig. 6c) proposed specifically for the dynamic assessment of masonry buildings [68]. Interestingly, through the values of β R the same IMs as in [68] are identified as the optimal ones to study the seismic dynamic behaviour of the Holsteiner Hof building. The best ranked IMs (i.e., lowest β R ) include PGA; ASI; Effective Peak Acceleration, EPA; Improve Effective Peak Acceleration, IEPA; Cordova intensity, S a *; and Vamvatsikos intensity, either Sa or Sa. Now, results from IDA are employed to build collapse fragility curves to measure the probability of failure at different levels of seismic intensity. Collapse observations represent the ratio of the number of records that induced GC to the total number of records at an IM level, to which the best distribution is fit using the SSE method, which exhibits the highest computationally efficiency among all the three methods reported in Section 2.1 and leads to practically same results [30]. The values of  η and β and corresponding analytical functions fitting the collapse data points through a lognormal cumulative distribution (LogCDF) are shown in Fig. 7 using the same IMs as for IDA curves. Alternatively, collapse fragility could be defined as the likelihood of exceeding the global drift ratio limit of 0.45 % defined in Fig. 6 for the GC. A deeper discussion of these functions is presented in the following subsections when contrasting the fragility analysis accounting for uncertainty in modelling parameters vs. record-to-record variability, and the IDA approach vs. PSDA with unscaled records. 4.1.2. Material and modelling uncertainties vs. record-to-record variability While the uncertainty induced by the record-to-record variability (i. e., random nature of the ground motion) is quantified by the β R parameter in this research, modelling parameters are assumed to be deterministic, thus, neglecting them as a source of uncertainty. Nonetheless, the investigation by Tomi´ c et al. [40] accounted for the uncertainty of material and modelling parameters by performing LHS to generate 400 variations of the 11 parameters reported in Table 2. As in the current research, collapse fragility curves were generated, assuming PGA as IM. In this regard, the fitting parameters for the Holsteiner Holf building accounting for modelling uncertainties are  η U =0.3592 and β U =0.0873. Such parameters showed consistency with the set of fitting parameters exposed in Fig. 7a for PGA. Ideally, both sources of uncertainty should be considered simultaneously. Nevertheless, accounting for both in the seismic assessment of historical masonry buildings can be extremely expensive from a computational standpoint, and it is the main reason for its consideration separately. Alternatively, Cornell et al. [91] Fig. 5. Failure statistics — Holsteiner Hof. D. Caicedo et al. Engineering Structures 319 (2024) 118804 8 proposed a simplified method to combine analytically both sources of uncertainty by means of Eq. (3): βRU = βR2+βU2 √(3) Hence, Fig. 8 illustrates the impact of modelling uncertainties and record-to-record variability on the collapse fragility of the Holsteiner Hof building, plus an additional curve accounting for both effects through Cornell et al. approach [91]. The median for the mixed effect was fixed as  η =0.3259 (shown in Fig. 7a), however, adopting  η U leads to practically indistinguishable results because of the proximity in the median values. Although similar, the effect of record-to-record variability leads to a slightly more brittle response, reaching, for instance, a collapse probability near 80 % at PGA =0.4 g, while at the same level of intensity the probability is lower than 70 % when accounting solely for modelling uncertainties. On the other hand, for the combined effects it is shown how, in contrast to the other curves, higher collapse probabilities are obtained up to an intensity level close to 0.4 g. Contrastingly, lower collapse probabilities are obtained after this threshold. 4.1.3. MR-IDA vs. PSDA Cloud-based fragility curves can be obtained after the aggregation of the results from probabilistic seismic demand models (PSDMs). Within this methodology, the structural system is examined under the action of unscaled ground motion records, in which the EDP vs. IM cloud response delivers the conditional probability of an EDP reaching or exceeding a certain level of edp, at a given level of seismic intensity, i.e., P(EDP ≥ edp|IM) [92]. Fig. 9 compares the fragility functions at collapse DS, obtained after MR-IDA in the current study, and the cloud-based approach using the results reported in the research by Caicedo et al. [68]. Clearly, the MR-IDA approach provides much more conservative results that, as previously shown, match entirely the actual distribution of numerical collapse observations reproduced herein. On top of that, the suite of non-linear dynamic analyses conducted in [68] led the Holsteiner Hof building to global collapse only in 9 cases which, evidently, are not enough to characterize the collapse distribution in the whole range of seismic intensity. Moreover, within the cloud-based approach, the moments of the distribution (i.e., median and standard deviation) are estimated using the sample moments from a set of data, and further considerations should be accounted for including numerical failures, leading to inaccurate results when dealing with lack of data [30, 93,94]. 4.2. Lausanne Malley For the second case study the bins of intensity are distributed as [0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40 0.50 0.60]. Different than for the first case study, at PGA =0.40 g mostly collapse is observed, thus, bins of 0.5 and 0.6 were utilised merely to discard any further structural resurrection. In this regard, 172 out of 210 analyses led the building to global instability. Fig. 10 presents the statistics of the failure analysis for the Lausanne Malley building. These failure observations show similarity with the results of the Holsteiner Hof building being distributed as 55 % OOP and 45 % IP. Also, similar to the first case study, flexure failure is predominant within IP collapse with 68 % of occurrences, followed by 29 % of mixed failures, and only 3 % (i.e., 2 collapse observations) by pure shear IP. Moreover, IP failures are predominantly observed at upper stories, actually, more than 70 % are concentrated at the last storey, either F5 xFig. 6. IDA curves displaying maximum values of EDPs — Holsteiner Hof. D. Caicedo et al. Engineering Structures 319 (2024) 118804 9 dir or F5 y-dir, but mostly in the “y” direction with almost 50 % of the counts. On the other hand, OOP failures induced by the overturning of the main façade and central node overturning are, as expected [95], observed at upper stories, with almost 50 % at F5 x-dir. Approximately, 30 % of OOP collapses are evenly distributed among F4 x-dir, F4 y-dir, and F5 y-dir; while less than 20 % are sparsely distributed among F1, F2, and F3 in both directions “x” and “y”. 4.2.1. IDA curves and fragility functions For this case study, the IO limit is established in the fractile curves at max(Δr) lower than 10 mm which is around 0.05 % of the equivalent global drift ratio. According to all three curves (i.e., 16th, 50th, and 84th percentiles) this limit is expected to be reached at a 0.05 g fraction of PGA. On the other hand, CP points are identified for max(Δr) values in the range of 20 mm to 40 mm, or global drift ratios ranging from 0.15 to 0.25 %, respectively. In terms of PGA, CP is reached between 0.10 to 0.20 g as appreciated in Fig. 11a. GC is observed when the structure approaches values of max(Δr) in the order of 90 mm or 0.55 % global drift. In this case, stable quartiles are estimated using 14 individual capacity curves. The observations of global drift limits to be used subsequently for deriving the analytical collapse fragility are once again consistent with the findings reported in the literature for building structures with similar features [88–90]. Moreover, Fig. 11b and c show the effect of dispersion, β R , on the capacity curves in terms of alternative Fig. 7. Fitting of analytical collapse fragility functions considering different IMs — Holsteiner Hof. Fig. 8. Impact of modelling uncertainties and record-to-record variability on collapse fragility curves — Holsteiner Hof. Fig. 9. Comparison of fragility curves from IDA and PSDA methodologies — Holsteiner Hof. D. Caicedo et al. Engineering Structures 319 (2024) 118804 16 [65] Vamvatsikos D, Cornell CA. Applied incremental dynamic analysis. Earthq Spectra 2004;20:523–53. [66] Vamvatsikos D. Performing incremental dynamic analysis in parallel. Comput Struct 2011;89. https://doi.org/10.1016/j.compstruc.2010.08.014. [67] Kazemi F, Asgarkhani N, Jankowski R. Machine learning-based seismic response and performance assessment of reinforced concrete buildings. Arch Civ Mech Eng 2023;23. https://doi.org/10.1007/s43452-023-00631-9. [68] Caicedo D, Tomi´ c I, Karimzadeh S, Bernardo V, Beyer K, Lourenço PB. Optimal intensity measure and probabilistic seismic demand model for the assessment of historical masonry buildings considering in-plane and out-of-plane response. Manuscr Submitt Publ Reliab Eng Syst Saf 2024. [69] Michel C, Karbassi A, Lestuzzi P. Evaluation of the seismic retrofitting of an unreinforced masonry building using numerical modeling and ambient vibration measurements. Eng Struct 2018;158:124–35. [70] Magenes G, Calvi GM. In-plane seismic response of brick masonry walls. Earthq Eng Struct Dyn 1997:26. https://doi.org/10.1002/(SICI)1096-9845(199711)26: 11<1091::AID-EQE693>3.0.CO;2-6. [71] Guerrini G, Senaldi I, Scherini S, Morganti S, Magenes G. Material characterization for the shaking-table test of the scaled prototype of a stone masonry building aggregate. Mater Charact Shaking-Table Test Scale Prototype a Stone Mason Build Aggreg 2017:105–15. [72] Senaldi I, Guerrini G, Scherini S, Morganti S, Magenes G, Beyer K, et al. Natural stone masonry characterization for the shaking-table test of a scaled building specimen. Proc Int Mason Soc Conf 2018;vol. 0. [73] Guerrini G, Senaldi I, Graziotti F, Magenes G, Beyer K, Penna A. Shake-Table Test of a Strengthened Stone Masonry Building Aggregate with Flexible Diaphragms. Int J Archit Herit 2019;13. https://doi.org/10.1080/15583058.2019.1635661. [74] Brignola A., Podest` a S., Pampanin S. In-plane stiffness of wooden floor. 2008 NZSEE Conference, Paper 49 2008. [75] Brignola A, Pampanin S, Podest` a S. Experimental evaluation of the in-plane stiffness of timber diaphragms. Earthq Spectra 2012;28:1687–709. [76] POLIMI. Critical Review of Methodologies and Tools for Assessment of Failure Mechanisms and Interventions, Deliverable 3.3. Workpackage 3: Damage Based Selection Of Technologies,; 2010. [77] Almeida JP, Beyer K, Brunner R, Wenk T. Characterization of mortar–timber and timber–timber cyclic friction in timber floor connections of masonry buildings. Mater Struct/Mater Et Constr 2020;53. https://doi.org/10.1617/s11527-02001483-y. [78] Vanin F, Zaganelli D, Penna A, Beyer K. Estimates for the stiffness, strength and drift capacity of stone masonry walls based on 123 quasi-static cyclic tests reported in the literature. Bull Earthq Eng 2017;15. https://doi.org/10.1007/ s10518-017-0188-5. [79] Akkar S, Sandıkkaya MA, S ¸enyurt M, Sisi AA, Ay B, Traversa P, et al. Reference database for seismic ground-motion in Europe (RESORCE). Bull Earthq Eng 2014; 12. https://doi.org/10.1007/s10518-013-9506-8. [80] Lanzano G, Sgobba S, Luzi L, Puglia R, Pacor F, Felicetta C, et al. The panEuropean Engineering Strong Motion (ESM) flatfile: compilation criteria and data statistics. Bull Earthq Eng 2019;17:561–82. [81] Disaster, Authority E.M. Turkish National Strong Motion Network 1973. 〈htt ps://doi.org/10.7914/SN/TK〉. [82] Luzi L, Hailemikael S, Bindi D, Pacor F, Mele F, Sabetta F. ITACA (ITalian ACcelerometric Archive): A web portal for the dissemination of the Italian strong motion data. Seismol Res Lett 2008. [83] Jayaram N, Lin T, Baker JW. A Computationally efficient ground-motion selection algorithm for matching a target response spectrum mean and variance. Earthq Spectra 2011;27. https://doi.org/10.1193/1.3608002. [84] Woessner J, Laurentiu D, Giardini D, Crowley H, Cotton F, Grünthal G, et al. The 2013 European seismic hazard model: key components and results. Bull Earthq Eng 2015;13:3553–96. [85] Danciu L., Nandan S., Reyes C., Basili R., Weatherill G., Beauval C., et al. The 2020 update of the European Seismic Hazard Model: Model Overview. EFEHR Technical Report 001, v1. 0.0 2021. [86] Code P. Eurocode 8: Design of structures for earthquake resistance-part 1: general rules, seismic actions and rules for buildings. Brussels: European Committee for Standardization,; 2005. [87] Loh WL. On latin hypercube sampling. Ann Stat 1996;24. https://doi.org/ 10.1214/aos/1069362310. [88] As ¸ıko˘ glu A, Vasconcelos G, Lourenço PB, Pant` o B. Pushover analysis of unreinforced irregular masonry buildings: Lessons from different modeling approaches. Eng Struct 2020;218. https://doi.org/10.1016/j. engstruct.2020.110830. [89] As ¸ıko˘ glu A, Avs ¸ar ¨ O. Investigation of Drift-based Damage Limit States for Historical Masonry Structures. Int J Archit Herit 2023;17. https://doi.org/ 10.1080/15583058.2022.2053612. [90] Kocaman ˙ I, Kazaz ˙ I. Global drift ratio limits for historical masonry mosques. Bull Earthq Eng 2023:21. https://doi.org/10.1007/s10518-023-01613-1. [91] Cornell CA, Jalayer F, Hamburger RO, Foutch DA. Probabilistic Basis for 2000 SAC Federal Emergency Management Agency Steel Moment Frame Guidelines. J Struct Eng 2002;128. https://doi.org/10.1061/(asce)0733-9445(2002)128:4 (526). [92] Jalayer F. Direct Probabilistic Seismic Analysis: Implementing Non-linear Dynamic Assessments. Stanford University,; 2003. [93] Zhang J, Huo Y. Evaluating effectiveness and optimum design of isolation devices for highway bridges using the fragility function method. Eng Struct 2009;31. https://doi.org/10.1016/j.engstruct.2009.02.017. [94] Jalayer F, De Risi R, Manfredi G. Bayesian Cloud Analysis: Efficient structural fragility assessment using linear regression. Bull Earthq Eng 2015;13. https://doi. org/10.1007/s10518-014-9692-z. [95] Costa AA, Penna A, Arˆ ede A, Costa A. Simulation of masonry out-of-plane failure modes by multi-body dynamics. Earthq Eng Struct Dyn 2015;44. https://doi.org/ 10.1002/eqe.2596. [96] Housner GW. Measures of severity of earthquake ground shaking. Proc US Natl Conf Earthq Eng 1975;1975:6. [97] Arias A. A measure of earthquake intensity. Seism Des Nucl Power Plants 1970. [98] Reed JW, Kassawara RP. A criterion for determining exceedance of the operating basis earthquake. Nucl Eng Des 1990;123. https://doi.org/10.1016/0029-5493 (90)90259-Z. [99] Riddell R, Garcia JE. Hysteretic energy spectrum and damage control. Earthq Eng Struct Dyn 2001;30:1791–816. [100] Fajfar P, Vidic T, Fischinger M. A measure of earthquake motion capacity to damage medium-period structures. Soil Dyn Earthq Eng 1990;9. https://doi.org/ 10.1016/S0267-7261(05)80002-8. [101] Cosenza E., Manfredi G. A seismic design method including damage effect. 11th European Conference on Earthquake Engineering, 1998, p. 6–11. [102] Park Y, Ang AH -S, Wen YK. Seismic Damage Analysis of Reinforced Concrete Buildings. J Struct Eng 1985:111. https://doi.org/10.1061/(asce)0733-9445 (1985)111:4(740). [103] Dashti S, Bray JD, Pestana JM, Riemer M, Wilson D. Centrifuge testing to evaluate and mitigate liquefaction-induced building settlement mechanisms. J Geotech Geoenviron Eng 2010;136:918–29. [104] Council AT, California SEA. of. Tentative Provisions for the Development of Seismic Regulations for Buildings: A Cooperative Effort with the Design Professions, Building Code Interests, and the Research Community. Department of Commerce, National Bureau of Standards,; 1978. [105] Yang D, Pan J, Li G. Non-structure-specific intensity measure parameters and characteristic period of near-fault ground motions. Earthq Eng Struct Dyn 2009; 38:1257–80. [106] Housner G.W. Intensity of ground motion during strong earthquakes 1952. [107] von Thun J.L., Roehm L.H., Scott G.A., Wilson J.A. Earthquake ground motions for design and analysis of dams. Earthquake Engineering and Soil Dynamics II - Recent Advances in Ground-Motion Evaluation: Proceedings of the Specialty Conference, 1988. [108] Yakut A, Yılmaz H. Correlation of deformation demands with ground motion intensity. J Struct Eng 2008;134:1818–28. [109] Cordova PP, Deierlein GG, Mehanny SSF, Cornell CA. Development of a twoparameter seismic intensity measure and probabilistic assessment procedure. Second US-Jpn Workshop Perform-Based Earthq Eng Methodol Reinf Concr Build Struct 2000;vol. 20:0. [110] Vamvatsikos D, Cornell CA. Developing efficient scalar and vector intensity measures for IDA capacity estimation by incorporating elastic spectral shape information. Earthq Eng Struct Dyn 2005;34:1573–600. [111] Hariri-Ardebili MA, Saouma VE. Probabilistic seismic demand model and optimal intensity measure for concrete dams. Struct Saf 2016;59:67–85. [112] Zhou Y., Su N., Lu X. An elastic spectral value-based intensity measure for the incremental dynamic analysis of tall buildings. Proceedings of the 5th Kwang-Hua Forum on innovations and implementations in earthquake engineering research, Shanghai, China, 2012. D. Caicedo et al.