Impact of Long-Duration Earthquakes on Successive Earthquake-Tsunami Fragilities for Reinforced Concrete Frame Archetypes

Offshore earthquakes have the potential to trigger tsunamis which can affect coastal communities near and far. While less frequent than the source earthquake itself, these events can cause widespread damage, loss of life, and subsequent socioeconomic impacts. For example, the Tohoku earthquake and subsequent tsunami that struck Japan in 2011 resulted in over 19,000 deaths, 6,000 injuries, and left 342,000 people houseless. In addition, according to Kajitani et al. (2013), about 119,000 residential buildings were completely destroyed, while 685,000 other buildings were partially damaged, resulting in direct losses of $211 billion (USD). As a result of these long-term impacts, many communities are exploring community resilience planning to ensure their community can recover from hazard events. Ideally, this process requires modeling using a comprehensive numerical platform encompassing component-to-community modeling procedures (Koliou et al. 2020). This resilience-based analysis can help to mitigate the impact of natural hazards on coastal communities both before and after an adverse event.

Presidential Policy Directive 21 (PPD 21, 2014) defines resilience as the ability to prepare for, withstand, and recover rapidly from adverse events. This means anticipating potential challenges, mitigating their impact, and promptly recovering from any resulting setbacks (Wang et al. 2021). Community resilience frameworks often use fragility functions as a first step to estimate the initial damage that a target community may experience after a single or compound hazard scenario, such as an earthquake, a tsunami, or both. Based on these initial risk-based damage assessments, scientists and engineers can evaluate the current functionality of a community after an event as well as the recovery path of complex socioeconomic institutions within that community (Watson et al. 2020). In this way, the fragility functions developed for existing structural typologies in a target community can use the output of a hazard analysis for a selected scenario event to assess the initial damage in a community.

Earthquake fragility curves have been developed for different structural systems in the past using both deterministic and fully probabilistic approaches (Bandyopadhyay et al. 1990; Chen et al. 2022; Li et al. 2014). However, tsunami fragility functions are not as prevalent in the literature, and most of them have been derived through empirical and statistical methods (Chua et al. 2021; Suppasri et al. 2011). Furthermore, tsunami fragility functions developed through physics-based simulations are typically for far-field tsunamis, which are based on the loadings from the tsunami hazard (Alam et al. 2018; Attary et al. 2017b; Petrone et al. 2017). These past studies employed a nonlinear static procedure (NSP) on a nonlinear structural model (deterministic or probabilistic), with a constant or adjustable lateral loading profile created based on an estimated equivalent tsunami force.

In contrast to the previously mentioned studies focused on far-field tsunamis, several investigations (Alam et al. 2019; Attary et al. 2021; Carey et al. 2019; Park et al. 2012; Petrone et al. 2020; Rossetto et al. 2019; Xu et al. 2021) have developed near-field tsunami fragility functions that can capture the cascading damage effects between earthquake and tsunami by taking into account loadings from both hazards. These studies typically utilized a two-phase successive analysis procedure consisting of a nonlinear time history analysis (NLTHA) for the earthquake loading, followed by a force-based NSP analysis for the tsunami hazard. However, these studies did not incorporate long-duration earthquakes in the NLTHA part of the successive simulation. Evidence from previous studies (Ando and Nakamura 2013; Carvajal et al. 2017) revealed that earthquakes with a large moment magnitude (M) can cause tsunami events, which were found to have a high positive correlation with the significant ground motion (GM) duration (Kempton and Stewart 2006; Mashayekhi et al. 2020; Salmon et al. 1992).

Previous research on the effect of earthquake strong motion duration across different structural systems, including reinforced concrete (RC) structures (Belejo et al. 2017; Hancock and Bommer 2007; Li et al. 2022), masonry (Bommer et al. 2004), dams (Xu et al. 2018; Zhang et al. 2013), steel frames (Barbosa et al. 2017; Chandramohan et al. 2016; Wang and Li 2022), and bridges (Su et al. 2023; Todorov and Billah 2021), has identified a marked positive correlation between earthquake duration and structural damage or collapse probability (Chandramohan et al. 2016). Furthermore, structures exhibiting deteriorative responses to seismic shaking, notably RC and masonry structures, have been shown to be more adversely affected by strong duration of earthquake motions, which aligns with the highest correlations observed in structural damage (Mashayekhi et al. 2019). The degree of this correlation is significantly influenced by the chosen engineering demand parameter (EDP) in these studies, encompassing both local and global damage indices. Local indices covered cumulative damage metrics such as the Park-Ang damage index (Park and Ang 1985), patterns of element cracks (Zhang et al. 2013), and hysteresis strength reduction (Bommer et al. 2004). Global indices incorporated measures like horizontal dam displacement (Xu et al. 2018), peak inter-story drift ratio (Belejo et al. 2017; Chandramohan et al. 2016), and floor acceleration (Li et al. 2022). While local indices, especially cumulative damage metrics, have shown a more pronounced correlation with the duration of strong earthquake motions (Hancock and Bommer 2007; Mashayekhi et al. 2019), peak interstory drift ratio (Belejo et al. 2017; Chandramohan et al. 2016) has been effectively employed in collapse and fragility assessments for higher earthquake intensity levels as an EDP, adeptly quantifying the impact of earthquake strong motion duration. It is important to note that these studies employed spectral matching to isolate the effect of earthquake strong motion duration from other amplitude-based earthquake parameters, aiming to solely quantify the impact of earthquake duration. This process involves adjusting or manipulating the frequency content of the original acceleration time series using spectral matching methods to match their response spectra to a target spectrum across all or a range of considered natural periods of vibration (Harati et al. 2021).

This paper employs a methodology proposed by Harati and van de Lindt (2024) to quantify the impact of long-duration earthquakes on the fragility surfaces of successive earthquake and tsunami hazards. The approach considers uncertainty characterization associated with the structural modeling procedure, input ground motions (GMs) for the NLTHA procedure, and NSP analysis for the equivalent tsunami force and its distribution for the lateral loading profiles. In this paper, we investigate and quantify the impact of long-duration earthquakes on the topological behavior of earthquake-tsunami fragility surfaces by applying both long-duration and regular-duration earthquake motions to a portfolio of RC frames, including both ductile and nonductile structures. The RC frame archetypes used in this study represent an RC building portfolio that could be used to populate a community-level model for a coastal community in areas exposed to tsunami and long-duration source earthquake risk.

There are various definitions in the literature for the strong motion duration of earthquake records. This study has selected significant duration for the definition of motion duration. Significant duration is based on a well-known integration-based relationship called Arias Intensity (AI), which provides a good representation of the earthquake input energy imparted to a structure. For a given earthquake motion, the AI can be calculated aswhere the AI of a motion is calculated using the acceleration time history of an earthquake record, denoted as $a(t)$, and a time window for the earthquake signal $[0\text{-}t]$, where $t$ is less than or equal to the total duration of motion. The AI of an earthquake record is a function of time, which can be calculated based on the $a(t)$ for the entire motion duration. Fig. 1 demonstrates the normalized AI of the 1989 Loma Prieta earthquake that has been calculated and plotted against motion duration.

As illustrated in Fig. 1, significant duration is typically denoted as DX-Y, where X and Y represent the lower and upper bounds defined on the AI function of an earthquake motion, respectively (e.g., 5% and 95% in this figure). In this case, the duration between 5% and 95% of the AI plot of an earthquake record is considered as the significant duration (referred to as D5-95 hereafter). While other forms of significant duration (e.g., D5-75) exist in the literature (Bommer and Martínez-Pereira 1999), the D5-95 is more commonly used in the GM selection procedure (e.g., Afshari and Stewart 2016), so it has been selected as the duration-related metric for this study.

Tsunamis and earthquakes are related natural hazards, but the damage they cause is not necessarily directly correlated. For example, a building exposed to near-field tsunami risk would first endure damage from the earthquake before being affected by the tsunami itself. Thus, it is crucial to consider earthquake and tsunami hazards in succession for communities susceptible to these types of multihazard events.

The buildings analyzed in this study focus on a suite of RC buildings, which include both ductile and nonductile RC MRFs, different heights, configurations, and design details. The MRF systems were chosen to represent both older and newer RC frame buildings located in an arbitrary US west coast coastal community vulnerable to seismic hazards.

This study specifically examines six RC ductile frames and four RC nonductile frames with space frame lateral resisting systems that had previously been designed and analyzed in other studies (Haselton et al. 2011; Liel et al. 2011). It is important to note that all 10 structural systems were designed as office buildings with a flat-surface floor system of 20.32 cm (8 inches) in thickness. RC ductile frames range from one-story to twenty-story levels, whereas RC nonductile frames range from two-story to twelve-story levels. Table 1 summarizes the geometric configurations of the lateral systems, as well as the corresponding design methods employed for each frame. The study, consistent with Liel and Deierlein (2008) assumes that nonductile RC frames were constructed before 1960 and designed according to the highest seismic zone in California at the time, which is equivalent to Zone 3 of the 1967 UBC code. Therefore, the design of these frames followed essential building code requirements such as maximum and minimum reinforcement ratios, stirrup designs, and maximum spacings. For more information on the nonductile MRF systems employed in this study, readers can refer to a scientific report by Liel and Deierlein (2008). On the other hand, the ductile RC frames were consistent with Haselton (2006) and designed according to more recent codes like International Building Code, ASCE07 (ASCE 2002), and ACI 318 (ACI 2002) to satisfy the applicable code requirements for strength, stiffness, design capacity, and dimensional details for structural members. Essentially, more ductility for structural members is provided in this case by implementing provisions such as strong column-weak beam ratios, enhanced joint shear capacities, and regulations regarding transverse confinement in plastic hinge zones. Despite adjustments to the period-based equation for designing base shear, the base shear coefficients for nonductile and ductile RC frames of equivalent height are relatively similar (Liel et al. 2011). This similarity holds except in the case of the shortest structures, where notable differences in the base shear coefficient may be observed. For more details about the design procedure for the ductile RC frames, Haselton (2006) provided a comprehensive report.

In this study, RC frames were numerically modeled using OpenSeesPy (Zhu 2021). OpenSeesPy was initially written in the TCL programming language but has now been converted to a Python library, allowing for the integration of structural analysis scripts into Python codes for parallel computations and postprocessing tasks. To reduce the computational demands of the two-phase successive analysis mentioned earlier, the 3D structural systems were transformed into 2D finite element (FE) models. However, it is important to note that this 2D modeling approach inherently implies that earthquakes and tsunamis are applied in the same direction. While this may be considered a conservative assumption, it does present a limitation of the current study by potentially simplifying the multidirectional nature of these hazards, but nevertheless the discrepancy between regular- and long-duration earthquakes is still investigated.

Fig. 3 presents the 2D structural model used for simulating successive hazard events. It contrasts two methods: the fiber and concentrated plasticity approaches, for modeling structures under tsunami or earthquake impacts. The fiber method is crucial for depicting inelastic behaviors along a structure, especially beyond column ends, essential for handling distributed tsunami loads. However, a major concern is shear failure in columns, particularly under lower inundation levels where shear demands are high, potentially leading to failure before reaching the first floor. Beam and column shear capacities, dependent on axial loading, vary across earthquake-tsunami phases, necessitating explicit simulation or post-processing for shear failure analysis. Force- or displacement-based distributed plasticity elements are effective here, allowing for moment-axial-shear interactions and a comprehensive view of structural responses to multihazard events. This method, demonstrated in prior studies (e.g., Alam et al. 2019, 2018; Carey et al. 2019; Petrone et al. 2020, 2017; Rossetto et al. 2019), incorporates both flexural and shear failures, using adjusted shear links or shear failure checks in postprocessing, providing a robust assessment framework for this multihazard scenario.

Conversely, other studies (e.g., Attary et al. 2021, 2017a, b; Haselton 2006; Haselton et al. 2011; Liel et al. 2011) have utilized the concentrated plastic hinge method. This technique involves lumped plasticity elements and finite joint shear springs, selected for their superior performance in capturing strain softening due to rebar buckling and spalling, essential for simulating the collapse of RC frames (Haselton 2006; Haselton et al. 2011). Based on empirical data from extensive testing on both ductile and nonductile RC columns, this method ensures accurate modeling. However, it primarily calibrates model parameters to axial compression levels in columns due to gravity loads and does not fully account for axial-flexure-shear interaction during the analysis, a factor that may be important in taller structures (Liel et al. 2011). While the fiber-spring model offers an alternative, it demonstrates less numerical stability in dynamic analysis and struggles with simulating severe strain softening, a critical failure mechanism in collapse scenarios (Haselton 2006; Haselton et al. 2011). In contrast, the lumped-plasticity model effectively captures this aspect, exhibiting severe strain softening at a roof drift ratio of approximately 3% or higher (Haselton 2006), thus proving more fit for seismic load simulations.

It is important to recognize that the concentrated plastic hinge method may not adequately simulate tsunami responses due to the development of inelasticity away from the seaward column ends under distributed loading. However, in the simulations, structures first undergo earthquake-induced damage, serving as an initial damage condition for subsequent analysis. This necessitates accurate modeling of the collapse or near-collapse response of the system. To address the limitations in simulating tsunami impacts, a postprocessing procedure is employed to incorporate the shear failure of seaward columns, thereby mitigating the trade-off errors, and enhancing the model’s overall accuracy in representing this multi-hazard scenario. As shown in Fig. 3, the structural elements for beams and columns are created using a linear elastic element connected to two nonlinear rotational springs at the endpoints. The nonlinear rotational springs provide regions with concentrated plasticity for these elements, allowing the formation of plastic hinges at the endpoints throughout the analysis. The seismic mass mounted to linear beam elements is calculated according to the guidelines provided by ASCE07 (ASCE 2022) regulations.

TThe modified Ibarra-Medina-Krawinkler (IMK) model (Ibarra et al. 2005; Lignos and Krawinkler 2011) is employed in this study as the nonlinear rotational springs, which can provide a deterioration model for the behavior of plastic hinges. As shown in Fig. 4, the IMK model can effectively control the strength and stiffness cyclic degradation of plastic hinges in both ductile and non-ductile RC members. This is mainly due to its ability to adjust based on the expected ductility an element can demonstrate during a numerical earthquake simulation. For example, the governing parameters of the IMK model can be tuned in such a way that it represents the hysteresis response of a nonductile RC column with limited ductility capacity. In addition, as shown in Fig. 4, the plastic hinge model is capable of capturing the cyclic behavior of elements with asymmetric hysteresis responses, such as RC beams working with concrete slabs in a fully composite manner.

In addition to defining plastic hinge regions in beams and columns, the modeling procedure also considers the shear behavior of beam-column joints. This is accomplished using a 2D joint element in OpenSeesPy as suggested by Haselton et al. (2011), which has four perimeter nodes that can be connected to the beam and column elements using IMK plastic hinges, as well as a central node that accounts for the shear deformation of the joint. For nonductile RC frames, a uniaxial material model recommended by Celik and Ellingwood (2008) has been assigned to the central node of the 2D joint elements to model the shear deformation behavior of the beam-column joints. On the other hand, in ductile systems, a uniaxial material with large elastic properties was used for the central node of the joints, as ATC (2017) reported that joints in RC ductile frames are compatible with the provisions of current seismic codes and exhibit a rigid behavior. In order to account for the interaction between the foundation and superstructure, linear spring links have been used in this study. The properties of these springs were determined based on the grade-beam design and soil stiffness properties (Liel et al. 2011). Furthermore, the structural level simulation incorporates viscous damping through the Rayleigh damping method, utilizing a random damping ratio with a mean value of 0.065, applied to the beam and column elements of the entire structural system.

Fig. 5 illustrates the pushover curve for the RC3 building’s structural model developed in this research, in comparison with the curve documented by Haselton et al. (2011). A detailed analysis of Fig. 5 reveals a notable consistency between these two pushover curves, indicating a robust correlation in their structural responses. The derivation of these pushover results in both studies stems from the application of median values of the input parameters used in the structural modeling process. For this study, the structural model is formulated on the foundation of the modified IMK model (Lignos and Krawinkler 2011). This contrasts with the approach of Haselton et al. (2011), who employed an earlier version of the IMK plastic hinge element (Ibarra et al. 2005) in their modeling.

Fig. 6 shows the scalar earthquake fragility function for the complete damage state, derived using the Multiple-Stripe Analysis (MSA) procedure and uncertainty modeling outlined later in this paper, depicted by a black line. The dashed blue curve illustrates the fragility findings of Haselton et al. (2011), considering the record-to-record (RTR) variability in earthquake motions. Furthermore, another fragility curve by Haselton et al. (2011) is shown with a solid blue line, integrating both the RTR variability and the structural uncertainties present in the finite element (FE) modeling process. Observations from Fig. 6 reveal a general alignment between these earthquake fragility curves. However, noticeable differences exist in certain segments, likely attributable to the varying assumptions and methodologies employed in their development. While the fragility curves in this study are derived from the MSA method and MC simulations, Haselton et al. (2011) based their earthquake fragility assessment on IDA and the concept of dynamic instability. These kinds of disparities in approach and results have been noted and discussed in other studies (Baker 2015; Pang and Wang 2021), highlighting the diverse methodologies and interpretations within this field of research.

Nonlinear verification procedures were uniformly applied across all structural models of the archetypes in this study, ensuring consistency across various typologies. Refer to Table 1 for details on structural geometry, design methodologies, and eigenvalues of archetypes.

To enable the formation of VH tsunami load profiles in the second phase of the successive analysis procedure, columns on the seaward side of the structure were divided into six segments at each story. This allows the tsunami loading profile to evolve vertically through the segments of each discretized element. These columns are expected to be subjected to relatively large shear demand from the tsunami flow. However, the presented modeling procedure explained herein fails to account for the direct shear failure of these columns during the analysis. To address this limitation, the shear demands of the seaward columns are checked against the computed shear capacity through a post-processing procedure, as discussed earlier in this paper. Specifically, if the peak drift ratio in one analysis realization is less than the predefined damage state, but shear demand in one or more columns exceeds the pertinent shear capacity, the entire analysis would be interpreted as a failed case. Table 2 provides detailed acceptance criteria for various building damage states, with PIDRs for DS1 to DS3 sourced from HAZUS. Importantly, the criterion for DS3 centers on the collapse potential of seaward columns, which is used to define the column shear failure (CSF) performance metric. It should also be mentioned that the modeling procedure explained herein does not consider the contribution of infill walls in the structural response, particularly under sequential seismic loadings. This omission represents a limitation of the study, as infill walls may influence the seismic performance of the structure, and their damage could also impact the forces induced by the tsunami. Therefore, the numerical simulation considers primarily the sideways collapse mode of the system, with an acknowledgment of the need for future research to comprehensively address the role of infills in this context.

As previously discussed, earthquakes with large moment magnitudes can generate tsunamis when they occur in seas or oceans. As illustrated in Fig. 7, produced through millions of GM simulations by using the relationship developed by Mashayekhi et al. (2020), earthquakes with longer motion durations can contribute to such multihazard events triggered by source earthquakes with larger moment magnitudes. Fig. 7 shows that the median value of strong motion durations of earthquakes varies exponentially with their moment magnitude in different soil types, indicating that larger significant durations (i.e., D5-75) should be expected from earthquakes with larger moment magnitudes. Therefore, selecting source earthquakes with long motion durations appears to be a more rational approach for numerical simulations of fragility functions in multi-hazard events, such as tsunamis. In this study, two different GM suits were selected: long-duration earthquakes and regular-duration earthquakes.

In this study, a set of synthetic GMs (referred to as the M9 Project records hereafter) generated for a variety of earthquake rupture scenarios on the Cascadia subduction zone (megathrust) approaching magnitude 9.0 have been chosen specifically for long-duration earthquakes (Frankel et al. 2018). For regular-duration earthquakes, the far-field ground motions recommended by FEMA P-695 (FEMA 2011) have also been selected as the reference earthquake motions in this investigation.

Spectral matching techniques (Hancock and Bommer 2007; Harati et al. 2021) were not implemented for GM records in this study to avoid altering earthquake time series to match a predefined response spectrum. It is important to recognize that response spectra are not the only critical intensity measure of earthquakes. Modifying earthquake records to conform to a single response spectrum can inadvertently affect other vital earthquake characteristics, such as frequency content and the displacement and velocity time series of ground motions (Harati et al. 2021). Other investigations (Chandramohan et al. 2016; Hancock and Bommer 2007; Harati et al. 2019) aimed to separate strong motion duration from other earthquake characteristics to determine its independent effect on structural damage. As noted by Mashayekhi et al. (2019), amplitude-based earthquake parameters (like PGA or SA) and strong motion duration (e.g., D5-75) are interrelated. Separating the duration component can alter the earthquake’s frequency content and modify the amplitude-based characteristics. Isolating the duration allows for assessing its sole impact on structural response metrics. However, this approach alters the final outcome, such as fragility functions. These modified vulnerability functions, having adjusted GMs in terms of frequency content and amplitude-based parameters, may not be reliable for use in risk and community resilience analysis due to potential inaccuracies in representing real earthquake characteristics. Fig. 8 presents the earthquake response spectra for the selected ground motions.

Two distinct categories for earthquake selection were utilized in this study: ground motions with long duration (M9 Project) and those with regular duration (FEMA P-695). In other words, the study investigates the impact of long-duration earthquakes (not the duration parameter of these motions) with their duration-based and amplitude-based characteristics against a set of regular-duration motions (and keeping all of its parameters) being employed as a reference GM set for developing fragility functions in the western United States. The normalized AI plots for the two selected GM suites are shown in Fig. 9, revealing two distinct patterns of input earthquake energy between the GMs from FEMA P-695 and the M9 Project.

Fig. 10(b) presents a scatter plot illustrating the relationship between significant durations (D5-95) and the PGA values of the records in both GM sets (long and non-long suits) selected for this study. It is clear from both Figs. 10(a and b) that these two record sets are distinctly separated in the vicinity of the D5-95 equal to 60 s, through inspection of the fitted lognormal distributions for their data in Fig. 10(b).

This study uses the methodology proposed by Harati and van de Lindt (2024) to develop earthquake-tsunami fragility surfaces for representative RC frame archetypes exposed to long-duration earthquakes and subsequent tsunami loading. A two-phase, probabilistic successive analysis approach, grounded in the Monte Carlo simulation (MCS), is used to determine joint failure probabilities for specific damage states, represented as discrete data points. Initially, MCS selects and scales a GM based on a joint IM, where the spectral acceleration (Sai) serves as the earthquake IM. The method involves generating uncertain parameters for both the FE nonlinear structural model and the calculation of equivalent tsunami forces. The FE integrates random variables, including those relevant to the development of plastic hinges, dead and live loads, damping ratios, and beam-column joint modeling, as indicated in Table 3. This uncertainty propagation ensures a realistic representation of structural behavior under varying conditions. Alam et al. (2018) and Attary et al. (2021) identified two primary variables for the equivalent tsunami force calculation, which are probabilistically modeled using uniform distribution as specified in Table 3.

Next, a nonlinear structural model is developed using the distribution functions from Table 3, and the equivalent tsunami force is calculated using Eq. (3). This equation considers random variables $C_{cx}$ and $Cd$, along with the momentum flux (MFi)—the tsunami IM in this study. The bay length (B) in Eq. (3) is constant, representing the structural system’s bay length perpendicular to the 2D building model. The fluid density parameter $\rho s$, accounting for water density with sediments in tsunami scenarios, is conservatively set at $1.2\mathrm{g}/{\mathrm{cm}}^3$ (Hashimoto et al. 2023). The uncertainty in parameters like inundation depth ($h$) and velocity ($u$) in Eq. (3) is incorporated in the second phase of the simulation. The algorithm calculates flow velocity from MFi and checks the Froude number (Fr) to ensure velocity is within acceptable limits ($\mathrm{Fr}\le 2$). If Fr is ≤ 2, the algorithm calculates the equivalent tsunami force for creating a lateral tsunami load profile. Successive analyses incrementally increase the tsunami flow depth ($h$), starting from zero or the water depth above the foundation, up to the building’s maximum height ($H_{\mathrm{m}\mathrm{a}\mathrm{x}}$). The flow depth’s initial and increment values are based on the columns’ discretization on the seaward side at each story, resulting in various VH profiles with different tsunami flow depths but identical MFi values.

After the successive analyses, the EDPs from both phases are compared against predefined damage states. Failure is declared for any MCS realization where the structure does not meet the damage state under the hazard scenario. Acceptance criteria for different building damage states are outlined in Table 2. To obtain a near-continuous function that governs the earthquake-tsunami fragility surface, a surface fitting procedure is performed on discrete data points of the joint failure probabilities. In this framework, surface fitting is carried out using an equation suggested by Kern et al. (2004)where the parameter $E$, along with the tuning parameters $C1$, $C2$, $\alpha$, and $n$ for the fragility surface over each desired damage state, were optimized using a nonlinear least squares (NLS) method. While Kern et al. (2004) found the value of $E$ to be 1.01, in this study, $E$ is also treated as a target variable for optimization, along with $C1$, $C2$, $\alpha$, and $n$. This approach results in five parameters being optimized to define the fragility surface. The NLS method was selected for its effectiveness in aligning model parameters with observed data, aiming to minimize the sum of squared differences between observed fragility data and the predictions of the model. This involved iterative adjustment of the parameter values, starting from initial estimates, to ensure the most accurate fit and minimize discrepancies.

Furthermore, a comprehensive sensitivity analysis was conducted to evaluate the impact of various optimization algorithms on the NLS problem. This included a systematic comparison of the Levenberg-Marquardt, Nelder-Mead, BFGS, and Conjugate Gradient algorithms, under consistent initial conditions and data sets. The convergence behavior, iteration count, stability, and sensitivity to initial parameter estimates of each algorithm were closely monitored. The robustness and efficiency of these algorithms were assessed in terms of solution reliability and computational resource demands, with the process repeated under different initial parameter estimates to test each algorithm’s performance under varied conditions. The findings indicated that the Levenberg-Marquardt algorithm was most effective for the specific NLS problem at hand. Combining aspects of the Gauss-Newton algorithm and gradient descent method, it showed superior convergence speed and accuracy. Its ability to adapt between the Gauss-Newton approach and gradient descent based on the problem’s nature and the residual sum of squares was particularly beneficial, given the complexity of the fragility function. The results of the sensitivity analysis guided the selection of the Levenberg-Marquardt algorithm, ensuring precision and reliability in predicting these multidimensional fragility functions.

Meanwhile, the employed methodology can indirectly generate two independent scalar fragility functions at the boundaries of the surface, one fragility curve for earthquake hazard and the other for tsunami hazard. These scalar curves appear on the edges of the horizontal plane composed of the relevant IM axes. A repetitive optimization-based procedure is used to extract the controlling parameters of the lognormal CDF, as follows:where $IM$ = variable for the discrete data points of failure probabilities for a particular hazard (earthquake or tsunami) under a specific damage state; and $\mu$ and $\sigma$ = mean and standard deviation of the vector IM considered. To achieve this, a function representing the lognormal CDF is defined, and the Nelder-Mead algorithm is used to identify the optimal parameters (mean and standard deviation) that provide the best fit. This process is applied across different data sets, each representing varying earthquake-only or tsunami-only hazard scenarios, ensuring a nuanced fitting of the lognormal CDF to the extracted edge data and accurately reflecting the statistical properties. The Nelder-Mead simplex algorithm is particularly well-suited for handling nonlinear optimization problems where derivatives of the objective function are not available or difficult to compute. Recognized as a well-established technique in the realm of unconstrained optimization, the Nelder-Mead method is especially known for its application in multidimensional optimization problems. Its primary advantage in this context is its ability to navigate complex, nonlinear spaces without necessitating the gradient of the objective function. The goal of the optimization is to minimize the discrepancy between empirical data and the theoretical lognormal CDF.

The optimized parameters $\mu$ and $\sigma$, obtained through the Nelder-Mead method, can be readily utilized to reproduce the scalar fragility curves for earthquake-only or tsunami-only hazards. This precision in parameter optimization, achieved through the effective use of the Nelder-Mead algorithm, ensures a more accurate and reliable representation of the scalar fragility curves, critical for assessing structural vulnerability in single-hazard scenarios.

The methodology outlined in this paper has been utilized to assess the earthquake-tsunami fragility surface of a 4-story ductile RC frame, which is one of the RC buildings included in this study. This is achieved by generating earthquake-tsunami fragility surfaces for the complete damage state using two different sets of GMs: the M9 Project for long-duration earthquakes, and the records from the FEMA P-695 for non-long-duration earthquakes. Fig. 11 displays the discrete data points of failure probabilities, along with the corresponding fitted surfaces for both sets of GMs.

As evident from Fig. 11, the fragility surfaces obtained from the two sets of GMs for 4-story RC buildings (ductile and non-ductile frames) exhibit a noticeable difference, indicating a shift in the fragility surface toward the vertical or z-axis when derived based on the M9 Project ground motions. This shift implies a higher probability of complete damage to the structure for a given joint IM of an earthquake-tsunami event when long-duration earthquakes are selected for the simulation procedure. For example, following the vertical black line with the arrow downward in Fig. 11(b) for an $Sa$ of 2g in this structure, the intensity of the tsunami event would have to be very large to have an exceedance probability approaching unity when exposed to regular-duration earthquakes. However, when subjected to long-duration earthquakes, the RC frame building reaches this failure probability with any value of the tsunami IM, indicating that it collapses before facing the incoming tsunami waves simply from the long-duration earthquake. This pattern of failure is also observed for the nonductile frame (RC9) depicted in Fig. 11(a), which exhibits a higher failure probability for lower earthquake IMs compared to its ductile counterpart, highlighting the critical influence of structural ductility on seismic resilience.

The earthquake-tsunami fragility surfaces for the 4-story RC ductile system at different considered damage states—slight, moderate, and complete—are derived for both sets of GMs selected in this study, as shown in Figs. 12(a and b). These fragility surfaces are accompanied by governing equations that enable analysts to calculate exceedance failure probabilities for a desired joint IM under each considered damage state enabling researchers to populate a database for community-level modeling with a portfolio of RC buildings. The surfaces in Fig. 12 include damage state exceedance probability calculations for the joint IM of (1.95, 0) as a demonstration. As shown on the fragility surfaces, the structures exposed to long-duration earthquakes have significantly higher failure probabilities compared to their regular-duration earthquake.

The earthquake-tsunami fragility surfaces for both sets of GMs for the complete damage state (DS3) are presented side-by-side in Fig. 13. The scalar fragility functions for tsunami-only and earthquake-only hazards are extracted from the edges of these fragility surfaces and displayed in Figs. 13(b and c), respectively.

As shown in Fig. 13(b), the scalar fragility curves for the tsunami-only hazard for both GM sets are closely aligned and share overlapping data points, indicating that the structures have a similar vulnerability to different levels of tsunami IM when prior seismic damage from the source earthquake is not taken into account. Since no earthquake is run for the tsunami-only analysis, one can see that the two surfaces converge in the three-dimensional surface plot for both earthquake sets to form the tsunami-only fragility. But as seen in Fig. 13(c), the scalar earthquake fragility curve of the long-duration set is shifted to the left side, indicating that more structural damage can be expected from long-duration earthquakes.

The degree of the fragility curve shift (${\mathrm{X}}_{\mathrm{sh}}$) has been quantified as a percentage, using the 50% cumulative probability of exceedance reference point for $Sa({\mathrm{T}}_1)$ for regular-duration motions. A positive shift is defined if the fragility curve generated using the long-duration motions is shifted to the left side compared to the reference point, demonstrating more initial seismic damage from the earthquake phase in this successive simulation. For example, in the four-story ductile RC frame illustrated in Fig. 13, the positive shift is approximately 35%. This shift could be attributed to several factors, such as longer shaking duration leading to more cumulative damage to the structure or different spectral contents and their interaction with higher vibration modes in the structure of interest. This initial earthquake damage would affect the structural response for the tsunami in the second phase of the simulation, as demonstrated by the earthquake-tsunami fragility surface derived from long-duration motions.

Fig. 14 displays the shift in fragility (${\mathrm{X}}_{\mathrm{sh}}$) between long-duration and regular-duration earthquakes for both ductile and nonductile systems. As previously noted, an $\mathrm{X}sh$ of zero suggests a negligible shift in response due to long-duration earthquakes, implying almost identical responses for both sets of GM records. The impact on low- to mid-rise RC frames, irrespective of their design ductility levels, is more significant from long-duration earthquakes, showing a shift of more than 40%, compared to their response under regular-duration earthquakes. In contrast, high-rise RC structural systems show minimal shifting when subjected to long-duration motions, indicating that their fragility behavior is relatively similar when exposed to earthquake sets from both the FEMA P-695 and the M9 Project. This pattern suggests that in taller structures, long-duration earthquakes do not notably alter the structural response compared to regular-duration earthquakes, likely due to their specific dynamic characteristics or the seismic input’s nature at higher modes. Moreover, the ${\mathrm{X}}_{\mathrm{sh}}$ values are directly influenced by the initial damage each structure endures, as depicted in Fig. 11. This figure clearly demonstrates that nonductile frames are inherently more susceptible to earthquake damage, regardless of the type. Consequently, the variance observed in fragility shifts for ductile versus nonductile frames should be interpreted based on the initial and reference damage states (based on regular-duration earthquakes) and inherent vulnerabilities of the considered structural systems, rather than as an indication of the direct structural responses of the buildings being examined under different seismic events.

Fig. 15 shows scalar tsunami fragility functions for collapse or complete DS, conditioned on earthquake IMs for two distinct GM categories: FEMA P-695 and the M9 Project. These fragility functions are illustrated across various earthquake intensity levels, with specific focus on $Sa$ values of 0.5g, 1.0g, and 1.5g. Each section of the figure demonstrates the changes in fragility functions corresponding to different $Sa$ values, highlighting how these variations are key in determining the likelihood of varying levels of structural damage under tsunami loading post-earthquake. The comparative analysis of the two GM categories at each intensity level sheds light on the impact of different GM characteristics on structural vulnerability during tsunamis following earthquakes. Notably, the fragility curves for tsunamis, conditioned on each $Sa$ value, exhibit a vertical compression and initiate from a higher cumulative probability in scenarios where a near-field tsunami follows a long-duration earthquake.

This vertical shift is evident at the beginning of each tsunami fragility curve under long-duration GM sets. Therefore, for this specific building archetype (RC3), the probability of collapse following a near-field tsunami event is generally higher when the structure experiences a seismic event characterized by a long-duration earthquake. However, it is essential to acknowledge that these observed vertical shifts in the fragility curves occur when there is a positive horizontal shift (${\mathrm{X}}_{\mathrm{sh}}$) in earthquake fragility curves. In cases where ${\mathrm{X}}_{\mathrm{sh}}$ is zero or minimal, such vertical shifts in the scalar tsunami fragility (conditioned on $Sa$) would not be observed. Furthermore, the rate of change in this vertical shifting tends to accelerate with larger ${\mathrm{X}}_{\mathrm{sh}}$ values, indicating a direct correlation between the magnitude of horizontal shift and the rate of vertical shift in the fragility curves.

The key parameters of the long-duration earthquake-tsunami fragility surfaces for ductile RC structures, have been extracted and presented in Table 4. These parameters ($E$, $C1$, $C2$, $\alpha$, and $n$) govern Eq. (4) and can be utilized to replicate the fragility surfaces for buildings assigned these archetypes for community-level modeling. Additionally, Table 5 contains the controlling parameters of the governing equations of near-continuous surfaces of nonductile RC structures.

The parameters governing the scalar fragility curves for earthquake-only and tsunami-only hazards have been extracted and presented in Tables 6 and 7, respectively. These tables contain the parameters controlling the lognormal CDF function introduced in Eq. (5), which plays a critical role in determining the vulnerability of RC structures to earthquake-only and tsunami-only hazards. It is worth noting that the scalar earthquake fragilities for the RC frames when subjected to long-duration earthquakes can be constructed using the parameters presented in Table 6. The scalar earthquake fragility curves developed in this study can be utilized to model and assess the vulnerability of RC structures in a region anticipated to be exposed to long-duration earthquakes only, i.e., not near the coast, which can aid in community-level planning and decision-making for seismic risk reduction.

This study utilized a fully probabilistic simulation within a two-phase successive analysis procedure to investigate the potential impact of long-duration earthquakes on RC frame buildings. This was accomplished through the investigation in the changes in earthquake-tsunami fragility surfaces for a selected suite of typical RC frame buildings. The successive analysis includes a NLTHA procedure for earthquake loading and an NSP analysis for the tsunami hazard. Additionally, uncertainty is included in the simulation at different levels for both the structural resistance and modeling and the earthquake and tsunami loading. The selected RC frames in this study are representative of both ductile and nonductile systems found in US (coastal) communities and form an RC building portfolio, varying in design philosophy, bay geometry, story heights, and the number of stories.

The vector-valued fragility surfaces for the earthquake and tsunami for the suite of RC frame buildings, exposed to both long-duration and regular-duration earthquakes, were evaluated in this study. The information pertinent to the parameters of scalar fragility curves for earthquake-only (with long durations) or tsunami-only, can be extracted directly from the edges of the fragility surfaces in their respective planes. Parameters for both the surfaces and scalar fragilities are provided for community modeling purposes and can serve as a data set for this suite of 10 RC frame buildings.

Additionally, studies have shown that initial seismic damage resulting from the source earthquake in low- to mid-rise structures increases the probability of failure by up to 40% when subjected to long-duration earthquakes. Specifically, RC buildings up to eight stories in height, including both ductile and nonductile RC frames, are at a higher risk of structural damage and failure (defined with DS1, DS2, or DS3) in the earthquake before the tsunami strikes. Furthermore, equivalent tsunami forces in the second phase of the simulation exacerbate the condition of structures with more initial earthquake damage as one would expect. This study indicates that the conclusions can potentially be applicable to a wide range of buildings facing strength and stiffness degradation, although their current scope is confined to RC frame buildings. Further research is essential to determine their generalizability to other types of structures.

Some or all data, models, or code that support the findings of this study are available from the corresponding author upon reasonable request.

The Center for Risk-Based Community Resilience Planning is a NIST-funded Center of Excellence; the Center is funded through a cooperative agreement between the US National Institute of Standards and Technology and Colorado State University (NIST Financial Assistance Award Nos. 70NANB15H044 and 70NANB20H008). The views expressed are those of the presenter and may not represent the official position of the National Institute of Standards and Technology or the US Department of Commerce.