Analytical fragility curves for masonry school building portfolios in Nepal

Schools represent a reference point for communities in any part of the world. Therefore, their safety and resilience against natural catastrophes is of paramount importance. The recent 2015 Gorkha earthquake has unfortunately shown that Nepalese school buildings are highly vulnerable to seismic actions. Most of them are indeed constituted by low-quality unreinforced masonry (URM). The quantification of URM vulnerability is fundamental to estimate the risk associated to school building portfolios at territorial scale. This work discusses statistics available for Nepalese schools and then presents analytical fragility curves for three recurrent URM typologies covering more than 50% of the school building stock. The methodology adopted to derive fragilities is spectral-based and accounts for out-of-plane and in-plane damage potential in a single easy-to-use analytical framework. Inter-building, intra-building and record-to-record variabilities are directly considered in the analysis. The obtained fragilities integrate the studies available for the region and can be used for pre-/post-earthquake risk assessment and prioritization of interventions at country level.


Introduction
Nepal is among the poorest economies in the world (Muzzini and Aparicio 2013) and one of the less urbanized nations in South Asia with 80.3% of the population living in rural areas (United Nations Department of Economic and Social Affairs Population Division 2018). Non-engineered and unsafe constructions (mostly unreinforced masonry, URM) still constitute most of the building stock and, analogously, school facilities are characterized by a high level of seismic vulnerability. According to a 2014 report by the Asian Development Bank (ADB), 18% of the nation's schools need to be reconstructed, while 43% require major retrofitting interventions (Asian Development Bank 2014).
Unfortunately, the last catastrophic 2015 seismic event worsened this scenario. According to different post-earthquake reconnaissance reports (e.g., Aon Benfield 2015;Build Change 2015;Government of Nepal 2015;Paci-Green et al. 2015;EERI 2016), approximately 6000-8200 schools were destroyed by the 2015 sequence of events. Post-earthquake surveys carried out adopting the inspection form from the National Society for Earthquake Technology (NSET) resulted in 6000 school buildings tagged with a damage grade (DG) between 4 or 5 (i.e., very heavy damage or destruction) and 11,000 tagged with DG2 or DG3 (moderate or heavy damage). In the Sindhupalchowk district, 99% of the classrooms exhibited a degree of seismic damage while, in Gorkha district, 85% of the classrooms collapsed. The Government of Nepal estimated the cost of damage and wider loss in the education sector to be of the order of $300-$400 million (Government of Nepal 2015).
In the face of this situation, governmental authorities and international non-governmental organizations (NGOs) have increased their effort in promoting coherent and coordinated actions on school safety (United Nations Office for Disaster Risk Reduction 2017). To achieve this goal, a systematic plan for assessing school facilities and prioritizing investments should be implemented by adopting state-of-the-art risk evaluation techniques (World Bank 2017). Therefore, the derivation of specific fragility curves for recurrent school building typologies is required (Gentile and Galasso 2019).
From a mathematical point of view, fragility curves express building physical damage in a probabilistic manner by describing the relationship between a relevant seismic intensity measure IM (e.g., peak ground acceleration, PGA) and the probability that a specific level of damage (DS i ) is exceeded (D'Ayala et al. 2015): Fragility curves are most commonly derived: (i) from expert opinions (expert/judgmental-based fragility curves); (ii) from statistical processing of post-quake reconnaissance data (empirical/observational fragility curves); or (iii) through analytical/numerical simulations (analytical fragility curves). Excluding the first approach, that could lack of scientific robustness (Aspinall 2010), the latter two techniques have been largely used in the last decades (e.g. Lagomarsino and Giovinazzi 2006;Jaiswal et al. 2011). In principle, observational fragilities should be considered the best option since they are based on damage evidence rather than simulations. However, their empirical nature is an intrinsic weakness (D'Ayala et al. 2015). Low quality of damage data and inaccuracy of earthquake shake-maps play an important role on the reliability of the results especially in data-scarce regions (McGowan et al. 2017).
Analytical fragilities are an alternative to overcome the limitations of empirical methods (Silva et al. 2019). Ideally, they should include all the sources of uncertainty e.g., geometry, material properties, static loads, ground motion, etc. (Calvi et al. 2006). Analytical fragilities based on numerical analyses, such as Finite Element Method (FEM) analyses, are generally time consuming as multiple FEM models need to be generated and analysed to include aleatory uncertainty. Therefore, several analytical fragility studies are based on simplified methods that require less computational effort.
In the Nepali context, some studies were carried out in recent years with the aim to develop observational and numerical fragility relationships. Guragain (2015) conducted a numerical study on four masonry types by adopting the Applied Element Method (AEM). Given the high computational cost of the AEM and the complexity of the calibration process (Malomo et al. 2018), fragility curves were derived considering eight models where random variation of geometry (e.g. walls thickness, interstory height, floor span, etc.) and (1) P(DS ≥ DS i |IM) mechanical properties (e.g., Young modulus, shear modulus, tensile strength, etc.) was not considered. Five scaled Ground Motions (GMs) were used to account for seismic input uncertainty. Chaulagain et al. (2016) derived fragility curves for general reinforced concrete buildings by adopting a numerical FEM approach. The same study includes some fragility curves for URM buildings whose source and typology is not reported. Gautam et al. (2018) developed observational fragility curves for three classes of Nepali residential buildings, namely, reinforced concrete (RC), brick masonry and stone masonry. This last work has the typical shortcomings of observational fragilities that use historical damage data (namely from 1934Bihar-Nepal earthquake, 1980Chainpur earthquake, 1988 Eastern Nepal earthquake, 2011 Eastern Nepal earthquake and 2015 Gorkha earthquake). More recently, the World Bank, within the Global Program for Safer Schools (ARUP 2015;World Bank 2017), has developed a global library of school infrastructure that includes information relevant to the Nepali context (World Bank 2019b). In particular, the database includes a set of index buildings for which fragility and vulnerability curves are reported. Additionally, thanks to a collaboration between the University of Bristol and the World Bank, a new empirical fragility model based on 18,000 school building damage data of the 2015 Gorkha earthquake has been recently produced (Giordano et al. 2020c). These works represent an important step towards the implementation of comprehensive seismic risk studies in Nepal (Robinson et al. 2018). However, our understanding of the vulnerability of the building stock in Nepal is still limited with respect to other regional contexts (Yepes-Estrada et al. 2016). A thorough risk/loss estimation should be performed with different fragility models to display a confidence interval on the results and to account for the epistemic uncertainties of the model (Hancilar et al. 2020). Additionally, referring to schools, a lack of tailored fragility curves could undermine the reliability of loss assessment and hinders the prioritization of building portfolios strengthening.
Given the widespread presence of masonry schools and their high vulnerability (Brando et al. 2017), this paper focuses on URM typologies and provides new fragility curves accounting for in-plane (IP) and out-of-plane (OOP) damage potential in a single easy-to-use analytical framework. Particularly, the OOP response is quantified with a recent model developed by the authors (Giordano et al 2020a). In details, Sect. 2 presents the characteristics of the URM school building stock by referring to data and information on the main structural typologies. Section 3 briefly reviews the existing methodologies for URM fragility assessment and presents the methodology adopted in this study. Section 4 presents the analysis of three recurrent masonry school typologies i.e. rubble stone in mud-mortar, brick in mud-mortar and brick in cement-mortar. Probabilistic input variables for the derivation of fragility curves are also discussed in the section. Lastly, new fragility curves for the considered URM typologies are presented and compared with existing studies relevant to the Nepali context (Sect. 5).

URM school building stock
According to previous studies (NSET 2000;ARUP 2015) there are mainly six school structural typologies in Nepal (Fig. 1). URM bearing wall buildings are generally made of: (i) field stones in mud-mortar [URM-SM] (Fig. 1a); (ii) fired bricks in mud-mortar [URM-BM] (Fig. 1b) Despite the large number of projects carried out on school facilities by governmental agencies and NGOs, a clear figure of the distribution of these structural typologies over the country is still lacking. ADB gives a rough indication of the school building types in different areas of the nation (Asian Development Bank 2014). In the plain area (terai), the RC framed structures are just 10% in favour of a more consistent presence of brick masonry (85%). These percentage distributions consistently change in the hilly and mountain regions where most of the buildings (respectively 87% and 97%) are recognized as "other type" i.e., stone masonry, adobe, wooden or mixed-system structures. Similar data are reported in post-2015-earthquake survey reports where most of the school buildings in Nepal (89%) are constituted by URM material (Paci-Green et al. 2015). Moreover, in mountain areas, because of the lack of construction materials such as fired bricks, cement and steel, 50% of these constructions are made of dry/mud-mortar rubble-stone masonry.
On the contrary, thanks to the vicinity to industrial activities, the distribution of school structural typologies in the Kathmandu Valley (KV) is different from the rest of the country. As reported by ADB (Asian Development Bank 2014), in the KV 30% of the schools have a RC frame bearing system while 65% have brick masonry bearing walls. A more detailed estimation of the school structural taxonomy and distribution in the KV dates back to 2000 (NSET 2000). This work gives an indication of the construction techniques present in the region starting from a survey-based catalogue of 909 school buildings. According to the study, at the time of the survey there were 643 public schools in Bhaktapur, Lalitpur and Kathmandu districts with approximately 1 to 9 buildings for each institute. Figure 2a reports the pie chart of the school dataset subdivided by structural typologies. It can be observed that about 65% of the school buildings are constituted by masonry material, split in 13% [URM-SM] and 52% [URM-BM/BC]. The second most recurrent school building type is the one-story steel structure. This typology represents 22% of the school building dataset. Lastly, reinforced concrete (RC) buildings account for 8% of the stock. For comparison purposes Fig. 2 reports additional building percentage distributions from other studies. In Fig. 2b, data from Chaulagain et al. (2016) show the general building stock distribution in the Kathmandu Metropolitan area. While the URM constructions still represent the 59% of the total stock, RC constructions account for 39% of the total. At national scale the comparison between school buildings and the general building portfolio is less accurate. ADB (Asian Development Bank 2014) provides a rough figure of school types (Fig. 2c) where 34% of the total is brick masonry, 5% is RC and the remaining part (61%) refers to other typologies i.e., adobe, stone masonry, timber, mixed structures, etc. For residential buildings (Fig. 2d), Gautam et al. (2018) provides an indication from the last census data where buildings are categorized according to their foundation material: 46% of the building stock have brick-mud or stone-mud foundation, 18% brick-cement or stone-cement foundation, 10% have RC foundation and the remaining part (26%) have timber pillars.
According to these statistics, URM buildings represent more than half the school/general building stock both at national scale and within the Kathmandu Valley. General characteristics of the three principal masonry school typologies (namely URM-SM, URM-BM and URM-BC) and recurrent seismic damage are discussed in the following subsections. 1 3

Rubble stone-mud masonry buildings [URM-SM]
Stone-mud URM schools are mostly present in hills and mountain areas (Robinson et al. 2018). Generally, they are one-story buildings with story height between 1.8 and 2.4 m. Walls are usually built with field stones whose shape and mechanical characteristics varies upon local availability. Their thickness is around 45-60 cm. The wall has usually a multiple-leaf morphology and lacks of through-stones (Wang et al 2018). Adequate wall-to-wall and wall-to-floor connections are also absent in most of the cases with the exceptions of state-of-the-art constructions with adequately built rectangular corner stone (e.g., World Bank 2019b). Floors are constituted by a 10 cm mud-layer supported by wooden planks and timber/bamboo joists. Roof is commonly made of Corrugated Galvanized Iron (CGI) sheets supported by small timber beams. As pointed out in numerous studies URM-SM are highly vulnerable to seismic actions (Gautam and Chaulagain 2016). Firstly, the absence of wall-to-wall and wall-to-floor connections results in frequent out-of-plane (OOP) damage (Sharma et al. 2016) as reported in Fig. 3a. Secondly, multi-leaf walls are affected by delamination under seismic shaking and, for this reason, the risk of OOP failure increases. Lastly, box-type response is essentially absent due to the high in-plane flexibility of the horizontal structures, so the walls behave almost independently being more vulnerable to OOP loads (Brando et al. 2017). These phenomena are obviously aggravated by the low and uncertain mechanical response of the mud-mortar.

Brick-mud masonry buildings [URM-BM]
In terai/hilly rural areas and in older city districts, it is likely to find URM school buildings made of fired bricks and mud-mortar (Robinson et al. 2018). URM-BMs have usually twostories and inter-story height around 2.7 m. The thickness of the walls is around 35-45 cm (ARUP 2015). Fired brick mechanical properties can differ from region to region and are highly influenced by the production process, i.e., manual or industrial. Floors and roofs are usually realized with traditional techniques such as timber-mud and timber-CGI. In some cases, traditional floors have been replaced with RC slabs.
URM-BM buildings are affected by many structural weaknesses that influence their seismic response. Despite bricks are arranged in regular patterns (e.g. English bond), the construction quality is generally low (Brando et al. 2017). Lack of wall-to-wall and wallto-floor connections, lack of bonding between wall leaves, lack of seismic detailing (e.g. tie 1 3 rods, anchors, ring beams, etc.) are common deficits of URM-BMs ). Due to the horizontal structure configuration, monolithic box-behavior is generally absent. For these reasons the predominant damage mode of URM-BM typology is the OOP, sometimes preceded by delamination (Brando et al. 2017). It has been observed that non-loadbearing walls and walls located at the upper stories are more vulnerable to OOP damage (Fig. 3b). This has been largely documented in numerous post 2015 earthquake reconnaissance reports (Chiaro et al. 2015;EERI 2016).

Brick-cement masonry buildings [URM-BC]
School buildings made of fired bricks in cement-mortar have been increasingly built in the last 30 years, especially in urban areas (NSET 2000). They are generally single storied buildings but, in some cases, can go up to four stories. Their average inter-story height is 2.7 m and are characterized by large openings. According to NSET (2000) and ARUP (2015) walls are generally one brick thick (23-24 cm) for single story and multi-story buildings while World Bank (2019b) reports multi-story URM-BC with thicker walls. URM-BCs have usually cast-in-place RC floors and roof.
Thanks to the presence of RC horizontal structures, they generally exhibit a box-behaviour response under seismic actions and in-plane (IP) damage of piers and spandrels is more recurrent (Fig. 3c). However, sometimes they lack adequate wall-to-wall/wall-to-floor connections and, consequently, suffer from OOP damage as well.

Seismic assessment methodology
The derivation of analytical fragility curves firstly requires the implementation of a proper seismic assessment methodology. This consist of four steps (Calvi et al. 2006): (i) Definition of a numerical/analytical model for the considered structural typology.
In this phase, the analyst selects the level of complexity of the structural simulation, ranging from analytical formulations to detailed numerical models. Several studies on URMs relies on simplified equations to define seismic capacity (D'Ayala and Speranza 2003; Lagomarsino and Giovinazzi 2006;Borzi et al. 2008; D'Ayala and Paganoni 2011). This seems a reasonable approach to be followed in the context of Nepal where advanced numerical techniques can be limited by the lack of building/ material data. Analytical formulations adopted in this study are reported in Sects. 3.1 and 3.2. (ii) Selection of a seismic analysis method. Looking at previous literature studies, three main methods have been adopted to derive fragilities for URMs: nonlinear time history analysis (e.g. Park et al. 2009), incremental dynamic analysis (e.g. Rota et al. 2010) and nonlinear static analysis (e.g. Borzi et al. 2008;Frankie et al. 2012)). The latter method involves simplified capacity spectrum techniques such as the CSM (Freeman 1978) or the N2 (Fajfar 2000). Spectral approaches are certainly less advanced with respect to the other methods, but computationally faster and still quite reliable (e.g. Lagomarsino 2015; Lagomarsino and Cattari 2015). For these reasons, they seem a valid option for fragility assessment of Nepali URMs. Details on the seismic assessment method adopted in this study are reported in Sects. 3.3 and 3.4.
(iii) Definition of criteria to identify damage states (DSs). In analytical fragility studies it is common practice to identify damage states as displacement thresholds on the capacity curve (e.g. Lagomarsino and Giovinazzi 2006;Simoes et al. 2015;Giordano et al. 2020b The assessment methodology adopted in this study considers OOP and/or IP damage potential depending on the horizontal structure typology. As observed in the aftermath of the 2015 Gorkha earthquake (e.g. Brando et al. 2017), masonry buildings with traditional flexible floors cannot respond in a box-like manner and, therefore, are mostly affected by OOP damage (Giordano et al. 2020a). On the contrary, rigid RC slabs guarantee (EERI 2016): (i) a better wall-to-floor restraint against OOP and (ii) a uniform distribution of the seismic forces among bearing walls that, consequently, suffer from IP damage. Three recurrent URM categories are analysed in this work: URM-SM with traditional flexible floors, URM-BM with traditional flexible floors, URM-BC with rigid RC slabs. It is worth mentioning that these categories do not cover the entire masonry school building stock but, according to the literature (NSET 2000;ARUP 2015), represent a consistent proportion of the total. Interested readers can consult the Global Library of School Infrastructure (World Bank 2019b) for an extended collection of index school buildings in Nepal. Table 1 summarizes the damage considered for the different cases. Subsequently, CSM is adopted to estimate IMs for different DSs, as for previous studies (Lagomarsino and Cattari 2015;Giordano et al. 2020a).

Out-of-plane capacity
The OOP capacity of the constituting walls of the building is estimated with a closed-form mechanical-based formulation recently proposed by Giordano et al. (2019bGiordano et al. ( , 2020a) that relies on the following assumptions: • The OOP response of a vertically spanning URM wall is purely governed by bending (Griffith et al. 2004;Shawa et al. 2012;Ferreira et al. 2015; Degli Abbati and Lagomarsino 2017). • Since the nonlinear flexural deformations localize in the areas with maximum bending moment (Ferreira et al. 2015; Degli Abbati and Lagomarsino 2017), the wall is modeled as a system of rigid bodies and nonlinear hinges. • The moment-rotation relationship of the nonlinear hinge is computed starting from the moment-curvature (M-χ) diagram of the critical cross-section. The M-χ is calculated under the assumption that axial strains behave linearly in bending i.e., URM sections remain plane (Cavaleri et al. 2005;Brencich and Felice 2009;Parisi et al. 2016;Giordano et al. 2017). • The rotation θ is calculated through a constant integration of the critical cross-section's curvature over the integration length L i .
The model is capable to represent three boundary configurations i.e., cantilever, pinned-pinned and clamped-clamped wall by introducing the quantity h LV , i.e., the shear length of the wall. The integration length L i was calibrated through analytical-experimental comparisons (Giordano et al. 2020a).
where F is the base shear of the wall, Δ is the maximum horizontal displacement (top displacement for cantilever configuration and mid-height displacement for pinned-pinned and clamped-clamped configurations), E m is the masonry elastic modulus, B is the width of the wall, t is the thickness, N is the vertical load at the top of the wall, W = B⋅t⋅h⋅γ m is the self-weight of the wall (h and γ m are height and unit weight of the wall respectively) and f mb is the compressive limit of the units. α h is a non-dimensional coefficient ranging from 0 to 1, which defines the position of the resulting seismic horizontal force over the height of the wall. For instance, it is equal to 2/3 when an inverse triangular distribution is assumed (Doherty et al. 2002).
To adopt the CSM, force-displacement curves defined by Eqs. 2-4 are reported in a spectral acceleration (S a ) versus spectral displacement (S d ) plane (capacity curves, CC). As for Doherty et al. (2002), the transformation is carried out with the following equations: where M e is the effective mass of the equivalent SDOF system estimated by discretizing the wall into a finite number of elements with mass m i and dimensionless modal displacement δ i ; Δ e is the effective dimensionless displacement of the SDOF. The CCs are independent from the width B of the wall.

In-plane capacity
The IP capacity is estimated through the approach proposed by Lagomarsino and Giovinazzi (2006). The method was formulated for non-engineered masonry buildings. The CC of the building is assumed as elastoplastic bilinear. In details, the maximum pseudo-spectral acceleration capacity S a,y of the structure is defined as follow: where N S is the number of stories, α is the ratio between the effective resisting area of the building and its footprint, ξ is a reduction parameter which considers the non-uniform contribution of piers and ranges from 0.8 to 1. Values of τ and σ 0 are defined with the following equations: where τ 0 is the masonry shear strength, p is the floor overload and h S is the interstory height. The elastic limit pseudo-displacement S d,y is defined starting from the Eurocode 8 (European Commitee for Standardization 2004) simplified estimation of the period of a masonry structure: where H is the total height of the building. Subsequently S d,y is estimated from S a,y with the following well-known equation: The ultimate pseudo-displacement S d,u is calculated as the minimum between the one related to soft-story collapse S d,us and the one referred to uniform collapse mode S d,uu .
where δ u is the masonry ultimate drift ratio and Λ is the modal participation factor of the first natural mode calculated as follow:

Seismic demand
According to the CSM, the seismic demand is represented in terms of a spectral-acceleration versus spectral-displacement response spectrum, so called ADRS. The spectral shape of the ADRS is generally provided by building codes (European Commitee for Standardization 2004; ASCE/SEI 2010). For Eurocode 8, the spectral shape remains unchanged regardless the demand intensity. Therefore, when code-compliant spectral shapes are used for fragility curves derivation, record-to-record variability cannot be considered. To overcome this limitation, some studies adopted ADRS from recorded GMs. For instance, Shinozuka et al. (2000) expressed the seismic demand by: (i) selecting 80 real GMs and grouping them in intervals of 0.05 g; (ii) calculating mean and mean ± standard deviation ADRS spectra for each interval; (iii) intersecting these ADRS with a RC bridge CC to obtain corresponding fragilities. More recently, Rossetto et al. (2016) have included record-to-record variability in a nonlinear static procedure to derive fragility curves of RC buildings using 150 GMs from the SIMBAD (Selected Input Motions for displacement-Based Assessment and Design) database (Smerzini et al. 2014). They have also successfully validated the results with nonlinear time history analyses. The use of real ADRS has some drawbacks. Referring to masonry structures in OOP loading, Lagomarsino (2015) pointed out that the jagged shape could lead to inconsistent IM values for subsequent DSs. This issue does not occur when smooth code-like spectra are employed (Giordano et al. 2020a). A possible way to account for different ADRS shapes is to process GM records with smoothing spectrum procedures. For instance, Malhotra (2006) proposed coefficients for ADRS estimation from PGA, PGV (peak ground velocity) and PGD (peak ground displacement) of the GM record. This methodology is herein adopted to process and smoothen 467 × 2 = 934 horizontal GM records from the SIMBAD database (Smerzini et al. 2014). Resulting spectra are then used to account for record-to-record variability as described in Sect. 5. Figure 4a shows the subdivision of the SIMBAD GMs with respect to their PGA while Fig. 4b reports an application of the smoothing procedure for a selected record (January 17, 1994 Northridge earthquake, Los Angeles-UCLA Grounds-station code 24,688). It is worth mentioning that when assessing the OOP damage of walls located at upper stories, the spectral demand should be properly modified to take into account the filtering effect of the supporting structure (Suarez and Singh 1987;Calvi and Sullivan 2014). In the present study the equations proposed by Lagomarsino (2015) are adopted. The amplified displacement response spectra S d,amp. (T,z) at level z of the building is defined as: where S d (T) is the GM displacement response spectrum, and S d,Z (T,z) is: where T b is the period of the supporting building as for Eq. (10), ψ b is the corresponding modal shape (linear shape assumption), γ b = 3N S /(2N S + 1) is the coefficient of participation, η(ξ) and η(ξ b ) are damping correction factors of the wall and of the supporting building respectively. The function η(ξ) is defined in EC8, Part 1 (European Commitee for Standardization 2004). For illustrative purposes, Fig. 4b reports the amplified ADRS estimated at the first floor of a three-story building (total height 7.5 m).
From the elastic response spectra, the inelastic demand is quantified by adopting the Capacity Spectrum Method (Freeman 1978). The CSM relies on the definition of an equivalent hysteretic damping of the structure that has been largely adopted to the case of masonry structures (e.g. Giordano et al. 2020a;Lagomarsino 2015;Borzi et al 2008). Additional details on the quantification of the overdamped spectra are reported in Sect. 3.4, point 4).

Estimation of intensity measures for different DSs
Once capacity and demand are defined, the estimation of IMs for different DSs is performed. According to previous fragility studies on masonry, PGA is adopted as IM (e.g. Borzi et al. 2008;Frankie et al. 2012;Simoes et al. 2015). This choice is motivated by several experimental campaigns (e.g. Benedetti et al. 1998) where PGA resulted a good damage predictor for stiff URM buildings subjected to IP damage. More controversial is the reliability of PGA in predicting OOP damage (Dimitrakopoulos and Paraskeva 2015). However, a recent research by Giresini et al. (2018) on the OOP fragility assessment of . 4 a Subdivision of SIMBAD GM records by PGA, b comparison between jagged ADRS from GM record and smooth ADRS derived according to (Malhotra 2006). The grey-dotted line indicates the amplified spectrum at the first floor of a three-story building masonry walls has shown that PGA remains one of the most relevant IMs to predict OOP damage. Figure 5 reports the conceptual flowchart of the OOP and IP assessment procedures. For the OOP the following steps are needed (Fig. 5a) (Giordano et al. 2020a): (1) Walls classification Identification of vulnerable walls with respect to their boundary conditions and overburden load. Referring to URMs with flexible floors/roof: • Loadbearing walls located at the top-story (purple) are assumed with cantilever boundary configuration since flexible CGI roofs do not provide transversal top restrain to the walls (Brando et al. 2017). The presence of the roof translates in a vertical force on top of the wall and in a horizontal inertial load since the overburden is assumed to be unrestrained (Vaculik and Griffith 2017). • Non-loadbearing walls located at the upper stories (green) are characterized by cantilever configuration that involves one or more stories. No vertical forces are transferred from the roof/floors to these walls. • Non-loadbearing full-height walls (red) have cantilever boundary condition that spans over the entire height of the building. No vertical forces are transferred from the roof/floors to these walls. • Loadbearing walls located at the lower stories (blue) are assumed with clampedclamped boundary configuration.
These boundary assumptions are consistent with previous studies on the topic (Doherty et al. 2002;Ceran and Erberik 2013). URMs with rigid floors/roof are characterized by clamped-clamped loadbearing configurations since RC bidirectional slabs generally insists on all the walls.
(2) Capacity curves calculation For any configuration, S a -S d CCs are extracted with Eqs. 2-6. (3) DSs estimation As underlined by Lagomarsino (2015) the definition of actual OOP physical damage of URM walls is quite complex when using analytical methods. Theoretically only two DSs should be considered: the rocking activation and the overturning. In between these two states, the wall should ideally rock without accumulating permanent damage. However, experimental evidence has shown that rocking generates cracks, crushes and sliding of units. In absence of specific physical-based definitions, in this work OOP DSs are calculated in terms of displacement thresholds on the CCs, as for Lagomarsino (2015): • DS1 (slight damage): displacement at 70% peak S a capacity; • DS2 (moderate damage) displacement at peak S a capacity; • DS3 (severe damage): 25% of S d,u i.e. displacement at null S a or displacement at cross-section failure (upper limit of Δ as for Eqs. 2-4); • DS4 (near collapse): 40% of S d,u .
(4) PGA evaluation for any DS The PGA estimation is carried out with the CSM (Freeman 1978) as suggested by Lagomarsino and Cattari (2015). Starting from the (S d,DSi ; S a,DSi ) points over the CC, corresponding periods T DSi are estimated: subsequently the equivalent damping coefficient ξ DSi for any DS is calculated: where ξ 0 is the initial damping of the masonry structure, ξ h,MAX is the asymptote of the hysteretic damping which depends on the structure type and μ ς DSi is the displacement ductility defined as S d,DSi /S d,DS1 (ς ranges between 1 and 2). Lastly, the PGA is estimated by calculating the normalized ADRS (PGA = 1.0 g, damping equal to ξ DSi ) and by homothetically scaling the spectrum to intersect the CC at the given DS (Lagomarsino and Cattari 2015;Giordano et al. 2020a). As discussed in Sect. 3.3, for walls located at the upper stories, the ADRS shape is modified to account for the filtering effect of the supporting building. Once DS PGAs are calculated for each wall configuration, at any DS the minimum is extracted as representative of the whole building.
The IP assessment is executed with the following procedure ( Fig. 5b): (1) Estimation of the building resisting area (Sect. 3.2).
(3) DSs estimation. Displacement thresholds are evaluated over the CC with the criteria by Lagomarsino and Giovinazzi (2006) where corresponding physical damage is defined according to EMS-98 (Grünthal 1998 (4) PGA evaluation for any DS. PGAs calculation is performed as for OOP assessment.

Probabilistic input variables for Monte Carlo simulations
OOP and IP assessment procedures described in Sect. 3.4 are used to derive fragility curves for Nepalese URM school buildings by performing Monte Carlo simulations where input parameters variability is directly taken into account. Monte Carlo is a standard statistical technique that is mainly adopted in combination with low-computational cost models (D'Ayala et al. 2015) such as the ones described in Sect. 3. Four main sources of uncertainty are considered: (i) intra-building variability, e.g., mechanical properties, material density, acting loads, etc.; (ii) inter-building variability e.g., number of stories, interstory height, walls thickness etc.; (iii) CSM parameters variability; (iv) record-to-record variability.
Depending on the floor typology, OOP or OOP + IP assessment is carried out (Table 1). When OOP + IP is executed, resulting DS PGAs are the minimum values between OOP and IP estimations for any Monte Carlo realization. In Sects. 4.1, 4.2 and 4.3 points (i), (ii) and (iii) are discussed for the three considered structural typologies while point (iv) is addressed in Sect. 5.

Rubble stone-mud masonry buildings [URM-SM]
As reported in Table 1, URM-SMs are assessed with respect to the OOP damage potential. For the geometrical parameters, indications provided by ARUP (2015) and NSET (2000) are considered. The story height h S is assumed normally distributed around a mean (μ) of 2.1 m with lower/upper truncation at 1.8 m and 2.4 m respectively. The coefficient of variation (CoV) is assumed equal to 0.3 as in Ahmad et al. (2018). The midspan of the flooring system L M , which is required to compute the walls acting load, is represented by a truncated normal distribution with μ = 2.0, CoV = 0.3 and max/min at 1.5 m and 2.5 m. Lastly, the walls thickness is uniformly distributed between 45 and 60 cm (ARUP 2015). An important aspect of the OOP analysis is to quantify the effect of inadequate transversal bond on the wall capacity (i.e., presence of multiple and poorly connected wythes). Lagomarsino (2015) suggests decreasing the wall flexural stiffness (and consequently the OOP capacity) by introducing a modified thickness t' = κ t, where κ is a non-dimensional reduction factor (κ < 1). t' does not modify the self-weight of the wall (still estimated with t) but solely applies when calculating the flexural stiffness and strength of the cross-section. An extensive numerical study by de Felice (2011) has shown that OOP strength and displacement capacity is highly affected by the bond between masonry layers. Particularly, a comparison between rigid-body analysis and discrete element analysis which accounted for the real cross-section geometry was carried out. It was observed that: (i) the average discrepancy of the two modelling approaches is around 30%; (ii) when the transversal bond is particularly weak, as for Nepalese stone walls (Brando et al. 2017), the overestimation of the rigid-body can go up to 70%. Recent experimental investigations by Giaretton et al. (2017) has also shown a systematic delamination of stone walls where approximately 2/3 of the cross section detaches and fails in OOP. Given these observations, for the URM-SM typology, a uniform distribution of κ is considered with lower/upper limit equal to 0.3 and 0.7 respectively. The modified thickness t' (and consequently t and κ) is a key input parameter for OOP assessment. This has been observed in a recent study by the authors on the sensitivity of out-of-plane capacity to input parameters of Nepali walls (Giordano et al. 2019a). Therefore, the probabilistic model of κ should be improved as soon as regional OOP experimental results become available.
Masonry properties required for the OOP assessment are taken from previous literature studies. As suggested in the Probabilistic Model Code (JCSS 2011), lognormal distributions are considered for material characteristics. The elastic modulus E m of rubble stone masonry in mud mortar provided in World Bank (2019a) is adopted. This results in μ = 240 MPa. In absence of specific information on the dispersion, CoV = 0.30 is considered. The compressive strength of the units f mb is derived from the stone properties reported in the Indian Standards (Bureau of Indian Standards 1974). This results in μ = 25 MPa and CoV = 0.28. Lastly the average specific weight of masonry γ m is assumed at 22 kN/m 3 (Bureau of Indian Standards 1987) while the CoV is assumed equal to 0.05 (JCSS 2001). Due to the lack of experimental data on Nepali masonry and the uncertainties around its response, the mechanical parameters have been assumed as non-correlated. An analogous simplification can be found in previous literature works (e.g. Ahmad et al 2018;Park et al 2009).
Variability of the CGI roof load p R is estimated from the Indian Standards (Bureau of Indian Standards 1987) by considering different type of CGI sheets and geometrical variation of the supporting timber/bamboo beams. This results in a lognormal distribution with μ = 0.15 kN/m 2 and CoV = 0.22. The floor load p F considers variability from: (i) geometry of the wooden supporting system (i.e. beam spacing, beam cross-section, planks thickness); (ii) timber specific weight (lognormal distribution with μ = 5 kN/m 3 and CoV = 0.10 (JCSS 2001)); (iii) mud-layer specific weight (uniform distribution between 17.25 and 18.85 kN/ m 3 (Bureau of Indian Standards 1987)); (iv) mud-layer thickness (uniform distribution between 8 and 12 cm); (v) prescribed live load (normal distribution with μ = 0.3 × 3 kN/ m 2 = 0.9 kN/m 2 as for Indian Standards (Bureau of Indian Standards 1987) and CoV = 0.25 (JCSS 2001)). The combination of these probabilities results in a lognormal distribution of p F with μ = 3.12 kN/m 2 and CoV = 0.10.

Brick-mud masonry buildings [URM-BM]
For the case of URM-SM, geometric properties are taken from ARUP (2015) and NSET (2000). Story height h S is assumed normally distributed with μ = 2.7 m, CoV = 0.30 (Ahmad et al. 2018) and lower/upper truncation at 2.4 m and 3.0 m. L M is assumed with the same distribution function of URM-SM. Walls thickness ranges between 35 and 45 cm while κ is assumed uniformly distributed between 0.7 and 1 to account for possible ineffective bond at cross-section level (Brando et al. 2017).
Lognormal distributions are assumed for material properties with mean values derived from the experimental tests on Nepali brick-mud masonry carried out by Rits-DMUCH   Table 3 summarizes the input probability functions for the OOP assessment. It is specified that probabilistic distributions and mean values adopted in this work are like the ones reported in a previous study by the authors (Giordano et al. 2019b) while dispersions are taken form JCSS (2001,2011).

Brick-cement masonry buildings [URM-BC]
Given the presence of rigid RC slabs, URM-BC buildings are assessed for IP + OOP damage potential. Geometric parameters are derived from ARUP (2015) and NSET (2000). Specifically, h S and L M are assumed with same distributions as URM-BM. Since brick cement walls are mostly one-brick thick, a narrow band for the thickness distribution is considered (23 cm-24 cm) and κ is assumed equal to 1. To execute the IP assessment, α (i.e. the ratio between the ground floor walls area and the footprint of the building) is also required. This quantity is generally correlated to the main geometric characteristics of the building (e.g., number of stories, horizontal structures span, etc.). However, in informal/ non-engineered URM Nepalese constructions this correlation can be very weak (for example it is quite common to see URMs where stories have been added during the years without making modification at the ground floor bearing walls). Therefore, in absence of specific data for Nepal, α is assumed uniformly distributed between 0.10 and 0.2 as suggested by Lagomarsino and Giovinazzi (2006).
As pointed out by Adhikari and D'Ayala (2019), material properties of Nepalese masonry can vary considerably due to the different state of deterioration of the material. As a consequence, in-situ tests generally provides lower values with respect to laboratory tests. To be consistent with the previous typologies, existing laboratory tests are adopted to characterize brick-cement URM. Material properties are considered lognormally distributed with CoV taken from JCSS (2011). In absence of specific data for Nepalese BC-URM, mean values are extrapolated from experimental works on typical brick URMs in the Himalayan region (Shahzada et al. 2012;Javed et al. 2015). Therefore: elastic modulus E m is assumed with μ = 1263 MPa and CoV = 0.25; shear strength τ 0 has μ = 0.14 MPa and CoV = 0.40. Compressive strength of the bricks f mb is assumed as for URM-BM. Lastly, the ultimate drift capacity ratio δ u needs to be defined for the IP analysis. Given the lack of experimental data on Nepalese URM-BC buildings, literature works on Pakistani URMs are used as a reference for δ u . Particularly: (i) building-scale test performed by Shahzada et al. (2012) resulted in δ u = 0.50%; (ii) 4 perforated wall-scale tests by Sajid et al. (2018) gave δ u values in the range of 0.30% and 0.50%; (iii) with an extended experimental campaign, Javed et al. (2015) estimated δ u at pier-scale obtaining μ = 0.54%, CoV = 0.37 and lower/upper limits of 0.28% and 1.05% respectively. Given the reasonable consistency of the three studies, the last statistical parameters are considered in this work to generate a truncated lognormal distribution for δ u .
Average specific weight of URM-BC is assumed equal to 18 kN/m 3 (Bureau of Indian Standards 1987) and CoV = 0.05 (JCSS 2001). Variability of the RC horizontal structure load is estimated considering: (i) nominal floor thickness between 10 and 20 cm (Bureau of Lastly, CSM parameters variability is treated differently for the OOP and IP assessment since box-behavior response lead to consistently larger ξ h,MAX values (Lagomarsino and Cattari 2015). Resulting probabilistic distribution parameters for URM-BC are reported in Table 4.

Monte Carlo analysis results
Monte Carlo simulations are performed by randomly generating N = 10,000 combinations of the input parameters in Tables 2, 3, 4. According to Byrne (2013), this allows to achieve a confidence level of 95% to obtain a mean result within 1% of the exact solution for CoV up to 0.5. For one structural typology, Monte Carlo tests have been performed with larger N showing that 10,000 is sufficient to achieve the prescribed accuracy. Record-to-record variability is accounted by adopting the following procedure: (i) For a given DSi displacement threshold on the CC of one Monte Carlo realization, a set of PGA DSi,k , k = 1: 934, is estimated with the procedure represented in Fig. 5a, b (step (4)) and described in Sect. 3.4. (point 4)). 934 represents the number of SIMBAD records (Smerzini et al. 2014) converted to smoothed ADRS spectra with the procedure by Malhotra (2006). In Fig. 6a an example of PGA DSi,k assessment is reported. The CC refers to the OOP of a cantilever wall of a single-story URM-SM building. It is observed that the 934 SIMBAD smoothed spectra intersect the CC at the given displacement damage threshold, in this case DS1. (ii) Each of these 934 spectra provides a PGA DS1,k and a corresponding scale factor SF DS1,k = PGA DS1,k /PGA GM,k , where PGA GM,k is the original PGA of the k-th GM record. Depending on the direction of the scaling, SF DS1,k can be larger than one (the k-th spectrum is scaled up) or lower than one (the k-th spectrum is scaled down). To quantify the absolute scaling, the quantity max{SF DS1,k ; 1/SF DS1,k } ≥ 1 is defined. Figure 6b reports the empirical cumulative distribution of this value for the considered example. The majority of the SIMBAD spectra intersect the CC with large and unrealistic absolute scale factors. It is also observed that the 10 th percentile of the distribution is 1.3. This means that by selecting a subset of 100 spectra (approximately 10% of the SIMBAD database) the absolute scale factor is maintained close to one and considerably lower than the maximum allowable values suggested in the literature (Bommer and Acevedo 2004). Therefore, a subset of 100 PGA DS1,k is selected. The dimension of the subset is comparable to the number of spectra used for record-to-record variability in previous studies (e.g., Rossetto et al. 2016). (iii) The last step of the procedure is to randomly select one PGA value of the subset to define the final PGA DSi . Referring to the example, Fig. 6c reports the empirical frequency distributions of PGA DS1,k . It is observed that the subset of 100-GM PGAs ranges in the same interval of the total set of 934 SIMBAD-GM PGAs. This means that the selection by minimum absolute scale factors doesn't affect consistently the record-to-record variability.  To ensure the consistency of the results, PGA DS# of each Monte Carlo realization should respect the consequential order of DS i.e. PGA DSi ≥ PGA DSi-1 , i = 2, 3, 4. This is always valid when adopting a single spectral shape for all the DS. Therefore, when different ADRS shapes are used, conflicting PGA DS could occur. To overcome this issue, the following check is performed: To quantify the different vulnerability of one-story and multistory buildings, the Monte Carlo analysis is repeated for six cases: (i) URM-SM1, single story; (ii) URM-SM23, twoto-three stories; (iii) URM-BM1, single story; (iv) URM-BM24, two-to-four stories; (v) URM-BC1, single story; (vi) URM-BC24, two-to-four stories. Among the various probability models suitable for fragility curves definition (De Risi et al. 2017), lognormal distribution is considered in this study. Figure 7 reports the cumulative distribution functions from the Monte Carlo analysis of URM-SM1 and the comparison with lognormal fits, where the method of maximum moments is adopted for the fitting. Figure 8 reports the fragility curves obtained from the Monte Carlo analyses while Tables 5, 6, 7 summarize the corresponding statistical parameters (median η and logarithmic standard deviation β). As expected, it can be observed that the presence of rigid  Table 8 reports the breakdown of the Monte Carlo results with respect to the governing damage mechanisms. Single-story buildings with flexible lightweight roof (URM-SM1 and URM-BM1) are governed by the OOP of loadbearing and non-loadbearing walls in equal parts. This result reflects the post-earthquake field observations. The overload of a CGI roof is almost irrelevant and does not provide any kinematic restraint against OOP displacement. Multi-story buildings with flexible floor/roof (URM-SM23 and URM-BM24) are ruled by three damage mechanisms. Particularly, the OOP of non-loadbearing walls at upper stories (green walls in Fig. 5a) is the governing mechanism in most of the analyses (41% for URM-SM23 and 48.7% for URM-BM24). Loadbearing walls at lower levels (blue walls in Fig. 5a) are less vulnerable to OOP loads and, consequently, do not affect the fragility results. OOP is the governing damage mode in 63.5% of the analyses of single-story URM-BC buildings while IP is more critical in 36.5% of the cases. For multistorey URM-BC buildings the percentage of IP and OOP is 56% and 43.7% respectively. As for the previous building categories, the OOP damage of loadbearing walls at lower levels is not a critical damage pattern. In general, the analytical results are in line with post-earthquake observations (e.g. Brando et al. 2017) that report consistent OOP damage of non-loadbearing walls and at the higher stories of URM-SM/BM buildings (Giordano et al. 2020a). On the contrary, IP damage resulted critical for URM-BC buildings with stiff floors/roof (EERI 2016).
To compare the results of this study with the existing fragility curves for Nepali URMs, an equivalence for different definitions of damage states is needed. As previously said, DS1, DS2, DS3 and DS4 of this study are consistent with the EMS-98 definition (Grünthal 1998). This agrees relatively well with the DS classification adopted by Guragain (2015) (i.e. slight, moderate, severe, collapse). On the contrary, Chaulagain et al. (2016) and Gautam et al. (2018) studies are based on three DSs, therefore the comparison is carried out for DS2, DS3 and DS4.
• Comparison with observational fragility curves for residential buildings by Gautam et al. (2018)-URM-SM1 and URM-SM-23 (Fig. 7a) results are conservative with respect to observational fragilities derived by (Gautam et al. 2018) (Fig. 7b). For any DS, median values are approximately 50% of the corresponding empirical estimations. Furthermore, lognormal standard deviations estimated analytically are lower  than observational values. The comparison between URM-BM fragilities and the corresponding observational ones shows instead a better agreement. These results are in line with what is observed in other regional contexts (e.g. Rota et al. 2008;Borzi et al. 2008) where analytical fragilities are generally more conservative and with less dispersion. Usually these differences depend on two aspects: (i) analytical fragilities neglect the epistemic uncertainty of the model, (ii) in the context of developing countries, empirical relationships are usually derived from incomplete damage datasets and rather poor shake maps (McGowan et al. 2017). • Comparison with fragility curves reported by Chaulagain et al. (2016)-Fragility curves reported by Chaulagain et al. (2016) are characterized by larger dispersion when compared with the analytical curves of this study. In terms of median PGAs there is a good agreement at DS2 (Moderate) and DS3 (Extensive) with one-story URM-SM/BM buildings. At DS4 (Collapse) η result larger for any building typology. Given the lack of information on the methodology used to derive these curves, it is not possible to draw any conclusion about the source of the discrepancies. • Comparison with analytical fragility curves by Guragain (2015)-For any structural type, median PGAs by Guragain (2015) appear in good agreement with the analytical results of this study. Dispersion vales are generally lower, probably due to the limited number of GMs adopted by Guragain (2015) to account for record-to-record variability.

Conclusions
In this work an analysis of the available data on the Nepalese school building stock has been produced by comparing statistics from different sources and geographical areas. Fragility curves for Nepali unreinforced masonry (URM) school typologies have been derived adopting a spectral-based methodology that accounts for in-plane/out-of-plane damage potential in a single easy-to-use analytical framework. Particularly, a recent closed-form solution for out-of-plane assessment has been integrated (Giordano et al. 2020a) with an existing approach for in-plane damage evaluation. Three URM structural types have been considered, namely, stone in mud (URM-SM), brick in mud (URM-BM), and brick in cement (URM-BC). The analyses have been carried out for singlestory and multi-story building classes to increase the granularity of the fragility model. Inter-building, intra-building and record-to-record variabilities have been accounted in the analysis. The analytical results have been compared with Nepal-specific fragility studies available in the literature. The new fragilities show a good agreement with respect to the analytical study by Guragain (2015). On the contrary, the comparison with observational fragilities (Gautam et al. 2018) shows larger discrepancies since: (i) the analytical results do not account for the epistemic uncertainty of the model and (ii) in the context of Nepal, empirical relationships are negatively affected by low-quality of damage data and poor shake maps. This study underlines that more experimental research on traditional Nepali masonry is needed to release some of the assumptions of the model. By adoptiong the same methodology, the present results can be easily updated when additional country-specific data become available. The new set of fragility curves, in combination with the existing studies, increases our understanding of the vulnerability of Nepali masonry schools and can be used for seismic loss assessment at building and portfolio scale.