Fault displacement hazard estimation at lifeline–fault crossings: A simplified approach for engineering applications

Lifelines, such as pipelines, roads, and tunnels, are critical infrastructure and when crossing active tectonic faults, a reliable estimation of the fault displacement in case of an earthquake is required. The first and simplest approach is to use empirical fault scaling relations to compute the design fault displacement, but this may result in an unknown level of safety. Thus, the probabilistic fault displacement hazard analysis (PFDHA) is the appropriate tool to assess the fault displacement hazard within a performance-based framework. Based upon an established PFDHA model, we present a simplified approach for engineering applications focusing on the lifeline–fault crossing along with appropriate simplifications and assumptions to extend its applicability to numerous faults. The aim is to provide a structure-independent approach of PFDHA that can be used when a site-specific study is not required, not possible (e.g., absence of recent sediments for dating past events), or too cumbersome, e.g., for lifeline route selection. Additionally, an in-depth investigation is presented on the key parameters, such as maximum earthquake magnitude, fault length, recurrence rate of all earthquakes above a minimum magnitude, and lifeline-fault crossing site, and how they affect the hazard level. This approach will be the basis for deriving hazard-consistent expressions to approximate fault displacement for use within the Eurocodes. The latter is intended to serve as a compromise between hazard-agnostic fault scaling relations and a comprehensive PFDHA, which requires detailed calculations and site-specific seismological data.


Introduction
Lifelines, such as fuel, water, and sewage pipelines, bridges, and tunnels, are critical infrastructures for the operation and prosperity of communities. These structures extend over a long distance and therefore multiple crossings of minor and major faults in seismic areas might be inevitable (Rondoyanni et al. 2013). In such a case, the fault activation will force the structure to develop excessive deformation with potentially devastating consequences, as identified in previous earthquake events for buried pipelines (O'Rourke and Liu 2012; Girgin and Krausmann 2016), above-ground pipelines (Honegger et al. 2004), and underground tunnels (Anastasopoulos et al. 2008;Ghadimi Chermahini and Tahghighi 2019). These facts render the implementation of the performance-based earthquake engineering framework (Cornell and Krawinkler 2000) necessary for the quantitative assessment of lifeline safety against earthquake hazard. A risk reduction strategy for the design, construction, operation, and maintenance of critical lifelines is, also, now dictated by pertinent guidelines under the Sendai Framework for Disaster Risk Reduction 2015-2030(United Nations 2015. The design fault displacement for lifelines founded on the ground surface or buried in the ground is typically estimated via a set of empirical fault scaling relations (e.g., Wells and Coppersmith 1994;Leonard 2014). This deterministic approach inherently leads to an unknown level of conservatism and safety (Bommer 2002) because the fault displacement is computed solely based on the fault dimensions and/or a characteristic earthquake magnitude. Even though this is a straightforward and easy-to-apply approach, the actual seismicity of the fault, the location of the crossing site on the fault trace, and the rupture length, among others, are disregarded. Given the importance of lifelines and in order to overcome the aforementioned drawbacks, one could employ a case-and site-specific Probabilistic Fault Displacement Hazard Analysis (PFDHA), which is the appropriate tool to quantify the potential of coseismic surface fault displacement. A comparison between these two alternative paths for estimating the design fault displacement is given by Melissianos (2022) through the investigation of three characteristic examples of faults in Europe.
The necessary data for PFDHA have been reported by Coppersmith and Youngs (2000) and its basis was established by Youngs et al. (2003) as the outcome of an extensive seismic hazard investigation that was carried out for the planned nuclear repository at Yucca Mountain, Nevada, USA (Stepp et al. 2001). Treiman (2010) presented a discussion from the geological point-of-view on the problems arising from the estimation of future fault displacement, focusing on rupture, deformation, and slip anticipation. Petersen et al. (2011) extended the PFDHA for strike-slip fault mechanism, providing more detailed calculations for off-fault displacements and focusing on mapping accuracy, fault complexity, as well as the plan area of the structure (site of interest). Chen and Petersen (2011) employed the approach of Petersen et al. (2011) to develop an empirical slip model for the southern San Andreas fault in California, USA, to reduce the uncertainties confronted in the design. Moss and Ross (2011) analyzed data collected from reverse fault ruptures and provided appropriate fault displacement prediction equations and conditional probability of slip functions. Recently, Chen and Petersen (2019) proposed an improvement regarding the uncertainty of rupture location for PFDHA by increasing the correlation between the structure's footprint, i.e., site of interest, and the fault displacement hazard.
Extensive research on the assessment of fault displacement hazards for engineering applications has been conducted mainly for the nuclear industry (Gürpinar et al. 2017;International Atomic Energy Agency 2010;Valentini et al. 2021). Contrarily, the pertinent studies for lifeline-fault crossings are scarce. Indicatively, for the case of a buried pipeline-fault crossing, Angell et al. (2003) presented an initial attempt to estimate the fault displacement hazard for offshore pipelines using the newly established methodology of Youngs et al. (2003). Melissianos et al. (2017a) provided a performance-based framework for evaluating the seismic risk of buried pipelines at fault crossings. Then, Melissianos et al. (2017b) used this framework to evaluate on a risk basis the effectiveness of alternative seismic countermeasures for buried pipes subjected to faulting. At the same time, Cheng and Akkar (2017) also addressed the probabilistic risk assessment of buried pipelines by estimating the fault displacement hazard using Monte Carlo simulations.
The application of PFDHA requires detailed calculations and a large set of data, which may not be available for every single fault crossing along the route of a lifeline. At the same time, design engineers are typically unfamiliar with these calculations and data, and so they tend to use the simpler empirical fault scaling relations. In that sense, our intension is to find the middle ground between the simpler approach (i.e., deterministic via the scaling relations) and the more "accurate" probabilistic approach (i.e., PFDHA). Thus, the aim is to formulate a methodology for a hazard-consistent estimation of the fault displacement using limited data, while being also simple enough to be applied by engineers. Such a methodology can be incorporated in a design code, such as Eurocode 8 (EN 1998) for the earthquake-resistant design of structures. At this point, it is important emphasize that the current European and international seismic design codes lack such provisions for estimating the design displacement for lifelines crossing active tectonic faults.
Owing to the above, the present study provides the context towards this direction by presenting an application of PFDHA for lifeline-fault crossings, with an emphasis on the caclulation decisions and the assumptions made. This simplified approach for engineering applications is an extension of the previous work on the seismic performance-based assessment of buried pipelines at fault crossings presented by Melissianos et al. (2017a). The methodology presented is austere by modern seismic design recommendations and it is tailored for the typical practitioner. Thus, the focus is only on the main properties of an active fault, namely, dimensions (length, width, and area), mechanism, and seismicity, as they can be found in available fault databases. Additional fault complexity such as secondary effects, fault segmentation, and linkage are not considered in this analysis. Furthermore, a parametric investigation is offered on the seismological parameters that contribute to the hazard level and shape the hazard curve, followed by a discussion on their influence. The uncertainty associated with the selection of alternative fault scaling relations that are required in the hazard calculations is, also, discussed. The outcomes are applicable to thousands of faults given fault characteristics present in a typical seismic source model used in Probabilistic Seismic Hazard Analysis [PSHA (Cornell 1968)]. The proposed approach is generic and applicable at a regional scale, however, for site-specific applications, in-depth characterization of the active faults via, for example, geophysical prospecting or paleoseismological investigations, is required. As a final remark, the proposed approach is structure-independent. Consequently, the obtained fault displacement should be converted to imposed ground displacement on the lifeline by taking into account the fault dip angle and the lifeline-fault crossing angle. Youngs et al. (2003) introduced two options for the PFDHA, namely the "earthquake approach" that follows the PSHA methodology and the "displacement approach" that requires paleoseismic data. Herein, the former is adopted because acquiring paleoseismic data for every single fault crossing of a lifeline extending over hundreds of kilometers is (1) practically not viable and (2) incompatible with seismic design codes. The "principal faulting" is also adopted, neglecting any distributed faulting issues and displacements "outside" the fault trace. Principal faulting stands for the slip along the main plane of the earth's crust discontinuity that leads to the release of energy during an earthquake event. The seismic productivity of each fault is assumed to be constrained by a rupture zone within the fault geometry. This assumption can be overly conservative for a building located in this rupture zone. On the other hand, it is a fair approximation for the total displacement imposed on a continuous lifeline that crosses the entirety of this zone, encountering both primary and secondary faulting. Moreover, one could consider the coseismic deformations persisting several kilometers around the fault trace, as well as that the rake can vary along the rupture trace. These considerations can improve the calculations, yet require site-specific data that is foreign to the seismic code and to most engineers. The key parameters used (variables and abbreviations) throughout the study are listed in Table 1.

Fault displacement hazard
The PFDHA yields the Mean Annual Frequency (MAF) of exceeding a predefined fault displacement value : where Δ ( ) is the MAF of fault displacement exceeding value , v is the recurrence rate of all earthquakes on fault above a minimum magnitude of engineering significance ( M min ) per year assuming a Poisson model, f M (m) is the probability density function of earthquake magnitude m being between a minimum M min and a maximum M max magnitude, P(Δ > |m, r) is the conditional probability that given an earthquake of magnitude m at a distance r from the site of interest has occurred, the fault displacement will exceed value , and f (r|m) is the conditional probability density function for distance r from the site to an earthquake of magnitude m occurring at the source under examination.
In the case of a lifeline-fault crossing, where the site of interest is located on the fault trace ( Fig. 1), the fault displacement hazard depends on the following: • Earthquake (moment) magnitude ( M ) is the key factor describing the size of an earthquake, ranging between values M min and M max . Earthquake magnitude is the conditioning variable of the calculations. • Rupture length ( RL ) is the segment of the fault that ruptures, acknowledging that different earthquakes may rupture fault segments of different lengths. Earthquake magnitude and rupture length are correlated (Wells and Coppersmith 1994;Leonard 2010). Rupture length ranges within RL min ≤ RL ≤ LF , where RL min is the minimum rupture length derived from M min via fault scaling relations and LF is the fault subsurface length that is typically reported by seismic hazard models and thus adopted herein. It should be noted that the rupture length at depth can be up to 25% higher than the surface rupture length (Wells and Coppersmith 1994). The rupture length is correlated with earthquake magnitude [see Eq. (12)]. • The position ( Pos ) of the rupture length ( RL ) on the fault trace is reflecting the uncertainty of the rupture location on the fault trace. Herein, by the term fault trace, we specifically refer to the projection of the (subsurface) top boundary of the fault plane onto the surface. This is typically the one reported in seismological maps and it is also the fault length employed in seismic source models. The actual manifestation of a surface rupture may or may not fall on this trace, which is beyond the scope of our investigation. • The location of the crossing point on the fault trace stands as the variable to account for whether the lifeline is intercepted by the ruptured segment, which is represented by XL that is the ratio of the distance (along the fault trace) of the lifeline-fault crossing point ( Z ) to the closest fault-end over the rupture length with 0 < XL = Z∕LF ≤ 0.50 . Symmetry around XL = 0.50 is considered, resulting in a larger variability (Youngs et al. 2003).
Based on the "earthquake approach" for principal faulting, as developed by Youngs et al. (2003), we focus on the fault displacement hazard calculations on the lifeline-fault crossing site by accounting for the aforementioned four elements. Using the total probability theorem and conditioning the hazard calculations primarily on the earthquake magnitude (Melissianos et al. 2017a), the MAF of exceeding fault displacement values is: where P Δ > |m i is the probability that the fault displacement Δ exceeds a defined value at the crossing site and P M m i is the probability of magnitude value ( m i ) falling within the i-th bin, following an appropriate discretization of magnitude values between M min and M max . Τhe lower bound M min represents the value below which no impact on the structure is assumed (minimum magnitude of engineering significance), while the upper bound M max is typically constrained by the fault dimensions. This probability of earthquake magnitude M being lower than a specific value m is estimated using the Gutenberg-Richter Bounded Recurrence Law (Gutenberg and Richter 1944): The b-value of Eq. (3) defines the slope of the curve that provides the "expected" future earthquake magnitudes. The effect of b-value is illustrated in Fig. 2 for a magnitude range between 5.50 and 8.50 and b-value between 0.70 and 1.10 as given in the 2013 European Seismic Hazard Model [ESHM13 (Woessner et al. 2015)]. The double truncation magnitude frequency distribution is the default assumption used in this study (Parsons and Geist 2009;de Santis et al. 2011;Page and Felzer 2015). The probability term P Δ > |m i in Eq.
(2) is a complementary cumulative distribution function (CCDF) and is estimated via the expression: where subscripts denote the integration over variable values: i for earthquake magnitude ( M ), j for rupture length ( RL ), t for average surface displacement ( ADS ), and k for RL j positions on the fault trace ( Pos ). Please note that variables RL and M , as well as, ADS and M are correlated via empirical fault scaling relations. Also, RL and fault displacement are correlated but the correlation of the residuals of ADD and RL conditioned on the magnitude is not provided in any pertinent study.
The first right-hand term, P Δ > |m i , RL j , ADS t , Pos j,k , of Eq. (4) is a CCDF and is estimated after Youngs et al. (2003) as: The term P Δ > |RL j , ADS t , Pos j,k , Slip of Eq. (5) is estimated as: where F(y) is the Conditional Probability of Exceedance depending on the fault mechanism. F(y) is a cumulative distribution function (CDF), estimated via empirical data models like the common approach for the ground motion prediction equations. Using historic data, fault displacement values ( ) were normalized by the average surface displacement values ( ADS ). The distribution of the ratio ∕ADS t is a function of the crossing point XL . The fault displacement is higher in the middle of the LF compared to the ends. In nature, the distribution of slips along the rupture could be skewed toward one of the two fault-ends (Ward 1997). Which end is favored is often unknown; therefore, the assumption of symmetry will tend to increase the variability (Moss and Ross 2011). The distribution of ∕ADS t is fitted with an appropriate model depending on the fault mechanism. The average surface displacement for the normalization of fault displacement data is represented by the variable y in the following Eqs. (7) and (8). The gamma distribution proposed by Youngs et al. (2003) is used for normal and strike-slip faults: while the Weibull distribution proposed by Moss and Ross (2011) is used for reverse faults: where XL j,k = x∕RL j is the crossing point relevant to the position ( Pos j,k ) of the RL ( Fig. 1 and 3). It is noted that in case the crossing site is not intercepted by the RL j at position Pos j,k under examination, then by definition the crossing point XL j,k is not computed and A graphical example of considering alternative RLs and positions of RL concerning the crossing point is illustrated in Fig. 3.
At this point, it should be pointed out that Youngs et al. (2003) and Moss and Ross (2011) provided the aforementioned data normalization also by the maximum fault displacement (MD). Wells and Kulkarni (2014) reported that the MD is expected to occur over a limited reach of the fault, and thus, adopting MD is a conservative approach. On the other hand, the AD is more representative of the expected displacement along a larger part of the fault because it is calculated using multiple measurements of displacement along the rupture zone (Wells and Coppersmith 1994). Therefore, hereinafter only the AD is considered. It is noted that both AD and MD approaches were examined by Melissianos et al. (2017a) within a logic tree formulation.

3
The second right-hand term, P RL j , ADS t |m i , of Eq. (4) is the joint probability of RL j and ADS t conditioned on the magnitude, being a probability mass function (PMF) due to the discretization of the variables. Estimation of this term is based on fault scaling relations as per Sect. 3. If the corresponding relations are provided at depth rather than at the surface, then a transformation is required. In the examined case, a simple linear transformation of ADD t = 1.32ADS t is adopted, based on the findings of Wells and Coppersmith (1994). In turn, translating P RL j , ADS t |m i to P RL j , ADD t |m i become as simple as keeping the individual probability masses constant and changing the AD coordinates that these refer to. Given the simplicity of this change, in the following lines, only the derivation concerning subsurface properties will be offered. Fault scaling relations allow the distribution of RL and ADD to be expressed as a joint lognormal function with positive correlation f (RL, ADD|M) that deals with the fault rupture at the subsurface level and consequently both quantities, namely rupture length and average displacement, refer to the subsurface level. A sufficiently fine discretization of the variables into RL j and ADD t is adopted, allowing the probability to be approximated via a single probability density function (PDF) value at its center via the expression: Note that the sum of P i,j,t over j and t for all RL j and ADD t should equal one and therefore a renormalization is performed for the derived PMF to account for any minor numerical issues: The function f RL j , ADD t |m i of Eq. (9) is estimated as a multivariate joint lognormal PDF using the mean value and the covariance matrix of the corresponding normal. The mean value [ ] is: where log 10 RL and log 10 ADD are the mean log values of RL and ADD given an earthquake magnitude m i and are computed using empirical fault scaling relations, such as those of Wells and Coppersmith (1994), Leonard (2014), and Thingbaijam et al. (2017), in the general form of: whereby the tilde operator ( ~) denotes that the functions g(⋅) and h(⋅) provide the distribution conditioned on their argument, for example, as defined by regression. Then, the covariance matrix [ ] for Eq. (9) is estimated as: where log 10 RL and log 10 ADD are the log standard deviations derived from the regression that produced the scaling relations of Eqs. (12) and (13), respectively, while log 10 RL,log 10 ADD is the correlation coefficient. ADD values based on RL are sampled using a fault scaling relation in the general form of: Still, the correlation of the residuals of ADD and RL conditioned on the magnitude is not provided in any pertinent study. A discussion on the effect of this correlation is offered in Sect. 4.
Finally, the third right-hand term, P Pos j,k , of Eq. (4) is the probability of a fault segment of length RL j rupturing at the k-th position along the fault trace. Lacking any better information, equal-length RL j are assumed to be located on all possible positions on the trace with the same probability (Youngs et al. 2003). Therefore, this term is a PMF and is estimated as: where N j is the total number of potential positions for the RL j estimated as: where the floor(⋅) function provides the largest integer lower than or equal to its argument and the minimum rupture length of interest ( RL min ) is estimated using a pertinent fault scaling relation and considering the minimum earthquake magnitude that is of significance for the case at hand.
The second term of Eq. (5) is the Conditional Probability of Slip (CPS), P Slip|m i , that accounts for the probability of the rupture reaching the ground surface conditioned only on earthquake magnitude. Petersen et al. (2011) reported that, for example, after the 1989 Loma Prieta, California ( M = 6.9 ) and the 2002 Nenana Mountain, Alaska ( M = 6.7 ) earthquakes, no surface rupture was observed. The CPS is estimated after Wells and Coppersmith (1993) for normal (Youngs et al. 2003) and strike-slip  faults: and after Moss and Ross (2011) for reverse faults: A comparison of Eqs. (18) and (19) is illustrated in Fig. 4 for earthquake magnitude values ranging between 5.50 and 8.50. A significant difference is observed, because, for example, to have over 80% probability [ P(Slip|m) > 0.80 ] for the rupture to reach the surface, a magnitude M > 8.50 event is required in the case of a reverse fault, contrary to a magnitude M > 6.77 event in the case of normal or strike-slip fault.

(14)
[ ] = 2 log 10 RL log 10 RL,log 10 ADD log 10 RL log 10 ADD log 10 RL,log 10 ADD log 10 RL log 10 ADD 2 log 10 ADD The fault displacement hazard curve at the crossing site is computed using Eq. (2) for a range of input fault displacement values ( ). A discretization scheme is adopted for the variables under consideration. At first, the earthquake magnitude range is discretized into several bins from a minimum value ( M min ) to a maximum one ( M max ) that is typically estimated using a fault scaling relation of magnitude given fault length: Then, the implementation of Eq. (4) requires the discretization of LF to multiple RL . Practically, a minimum RL is considered after Eq. (12) for M min , while all larger rupture lengths are integer multiples of this minimum one. Then, every RL is considered at all possible positions on the fault trace, which is represented by the variable Pos , keeping track of those that intercept the lifeline (Fig. 3). Only the latter contribute to the fault displacement hazard, contrary to what earthquake engineers are familiar with for the ground shaking in PSHA, where the contribution of all nearby seismic sources is accounted for at the site under investigation.
A case study to illustrate the outcome of the hazard calculations is provided subsequently. The baseline parameters considered for the lifeline-fault crossing are an interplate normal fault with length LF = 100km , recurrence rate v = 0.0066 years −1 , b = 1.00 , M min = 5.50 and M max = 7.57 . The lifeline-fault crossing point is located at a distance Z = 30km from the fault closest fault-end (crossing point XL = 0.30 ). The estimated fault displacement hazard curve at the crossing site is presented in Fig. 5. It is worth noting that contrary to a hazard curve for ground shaking obtained from PSHA, which in the vicinity of a dominant fault with a rate of v would tend to saturate at v for low enough intensities, the PFDHA curve will encounter a plateau at lower rates; this is due to the CPS, as low-magnitude earthquakes will rarely lead to surface rupture (Youngs et al. 2003;Moss and Ross 2011;Valentini et al. 2021). Finally, it is stated that the recurrence rate v is in fact an external multiplier of the seismic hazard calculations as per Eq. (2) and consequently the increase/decrease of v essentially leads to a proportional upward/downward parallel "shift" of the hazard curve to higher/lower MAFs for the same displacement.
(20) M max ∼ w(LF) Fig. 4 Empirical conditional probability of slip (CPS), i.e., probability of the rupture reaching the ground surface conditioned on earthquake magnitude ( M ), with respect to the fault mechanism (Youngs et al. 2003;Moss and Ross 2011;Petersen et al. 2011) 1 3

Fault scaling relations
Fault scaling relations are developed using data from historic earthquakes and relate fault characteristics and properties, namely earthquake magnitude, fault length, fault width, rupture area, subsurface rupture length, average fault displacement, maximum fault displacement. These regression equations are typically used within a deterministic seismic hazard assessment framework to compute the fault displacement offset for the fault crossing at hand. A comprehensive review of the available scaling relations can be found in Wang (2018), while a new set of relations especially for strike-slip faults has been very recently published by Anderson et al. (2021).
Hereinafter, the fault scaling relations of Wells and Coppersmith (1994), Leonard (2014), and Thingbaijam et al. (2017), abbreviated as WC1994, L2014, and TMG2017, respectively, are examined. In brief: • Wells and Coppersmith (1994) systematically examined a database of 421 historical earthquakes with a magnitude greater than 4.50, from which they selected 244 events with known parameters of the earthquake rupture. The data were obtained from continental and intraplate tectonic environments, while subduction zones and oceanic slabs were excluded. The authors employed the ordinary least-square regression method and provided a set of independent fault scaling relations for normal, reverse, and strike-slip fault mechanisms. • Leonard (2010Leonard ( , 2012Leonard ( , 2014 proposed two scaling relations; one between rupture width and rupture length, and the other between rupture displacement and rupture area. By substituting these scaling relations in the definition of seismic moment, the author developed a series of self-consistent scaling relations between seismic moment (or moment magnitude) and different rupture characteristics by applying ordinary least-square regression. The author has combined the data of normal and reverse fault mechanisms to a single category, that of dip-slip. Additionally, a distinction was made regarding the tectonic environment (TECT), namely Interplate (INT) and Stable Continental Region (SCR). It is noted that the Leonard (2014) relations are currently also used in the development of the active faults sources as part of the 2020 update of the European Seismic Hazard Model [ESHM20 (Danciu et al. 2021)]. • Thingbaijam et al. (2017) published a set of relations following the typical fault mechanism characterization as normal, reverse, and strike-slip like (Wells and Coppersmith 1994) and differently from Leonard (2014) who analyzed normal and reverse together. The database did not include SCR events but included subduction interface events in a separate set of relations. The authors employed general orthogonal regression, which allows for inverting the variables, and developed relations that are purely empirical, i.e., no prior scaling model was assumed to fit the relation coefficients, similarly to Wells and Coppersmith (1994) and differently from Leonard (2010Leonard ( , 2012Leonard ( , 2014.
Fault scaling relations are used explicitly within the PFDHA calculations at the following steps: • Estimation of the minimum rupture length ( RL min ) under consideration via Eq. (12) considering the minimum earthquake magnitude ( M min ). • Computation of the joint probability f (RL, ADD|M) [Eqs. (9) and (10) Thus, it is deemed necessary to examine the effect of selecting one over another set of relations. The (subsurface) rupture length may be estimated as a function of earthquake magnitude using the relations listed in Table 2 for Eq. (12). The comparison is presented in Fig. 6 for normal, reverse, and strike-slip fault mechanisms within a range of earthquake  magnitude between 5.50 and 8.50, where non-negligible variability is observed, especially for lower magnitude earthquakes occurring at normal faults and for higher magnitudes at strike-slip faults. As a remark, WC1994 and TMG2017 directly provide the standard deviation of the output variable, while L2014 provides it for the intercept a 1 of the relation M = a 1 + 1 log 10 (RL) . With parameter 1 being a constant, this translates to (M) = a 1 .
In this application, the relation has been inverted to obtain Eq. (12) with the corresponding standard deviation is log 10 (RL) = a 1 ∕ 1 . It is noted that Leonard (2014) intended the equation ( Y = a + X ) to be invertible since it is not fitted in the classical sense. The multiplicative constant is actually fixed, coming from a physical model, while optimization is performed only on the intercept a . Thus, the linear part of the model is essentially treated as deterministic, allowing the inversion to be mathematically correct. The average subsurface fault displacement ( ADD ) may be computed using the earthquake magnitude ( M ) [Eq. (13)] using the scaling relations listed in Table 3 for all fault mechanisms and tectonic environments. Please note that the relation ADS ∼ h(M) provided by Wells and Coppersmith (1994) for reverse fault mechanisms is not significant at a 95% probability level, and thus, it is replaced by the relation provided by Moss and Ross (2011), abbreviated as MR2011.
The obtained ADS values by Moss and Ross (2011), are also transformed to ADD values. Moreover, regarding the standard deviations, WC1994 and MR2011 provide the relation ADS ∼ p(M) with standard deviation log 10 (ADS) and thus the required standard deviation for the subsurface average fault displacement is obtained as log 10 (ADD) = log 10 (1.32ADS) = log 10 (ADS) + log 10 (1.32) = log 10 (ADS) . Then, L2014 provides the of parameter a 2 for the relation M = a 2 + 2 log 10 (ADD) . With parameter 2 being a constant, this translates to (M) = a 2 . The relation has been inverted to obtain Eq. (13) and thus the corresponding standard deviation is log 10 (ADD) = a 2 ∕ 2 . Finally, TMG2017 directly provide the output standard deviation for the relation ADD ∼ h(M) . The comparison of the different median (mean log) estimates of ADD based on the earthquake magnitude are illustrated in Fig. 7 for normal, reverse, and strike-slip fault mechanisms within a range of earthquake magnitude values between 5.50 and 8.50. It is observed that roughly in all cases, the SCR case of Leonard (2014)   The required sampling of ADD given rupture length RL is carried out using the scaling relations listed in Table 4 after Eq. (15). The pertinent comparison is depicted in Fig. 8, indicating significant variability in the results. Moreover, regarding the standard deviations, L2014 and TMG2017 provide the for the relation ADD ∼ r(RL) . Conversely, WC1994 provide the relation log 10 (ADS) = a 3 + 3 log 10 (SRL) , where SRL is the surface rupture length. In such case, both sides of the relation should be transformed to the subsurface level. Based on the mode of the distribution of ratios of average subsurface displacement to average surface displacement calculated by Wells and Coppersmith (1994), the transformation ADD = 1.32ADS is used for the average displacement. Then, the transformation SRL = 0.75RL is used for the rupture length, as Wells and Coppersmith (1994) found that,  on average, the surface rupture length equals 75% of the subsurface rupture. Owing to the above, the median log 10 (ADD) value may be estimated as: The corresponding standard deviation can be estimated as: as the variance of a random variable plus a constant is equal to the variance of the former. The maximum earthquake magnitude ( M max ) considered in the ] is computed using the scaling relations listed in Table 5 after Eq. (20).
(21) log 10 (ADS) = a 3 + 3 log 10 (SRL) → → log 10 (ADD∕1.32) = a 3 + 3 log 10 (0.75RL) → → log 10 (ADD) = a * 3 + 3 log 10 (RL) with a * 3 = − log 10 (0.76) + a 3 + 3 log 10 (0.75) (22) log 10 (ADD) = log 10 (1.32ADS) = log 10 (ADS) + log 10 (1.32) → → log 10 (ADD) = log 10 (ADS) At this point, it is worth recalling that alternative fault scaling relations are examined to capture the pertinent uncertainty on selecting alternative sets. To quantify this effect, the baseline parameters of the case study of Sect. 2 are adopted. The obtained fault displacement hazard curves are illustrated in Fig. 9 for all fault mechanisms and considering alternative sets of empirical fault scaling relations. It is observed that adopting the set of Thingbaijam et al. (2017) leads to lower MAF values, while adopting the Leonard (2014) fault scaling relations results in higher MAF values, especially for the interplate tectonic environment. The hazard curve when the Wells and Coppersmith (1994) relations are adopted are roughly in the middle of the other cases. In any case, the differences in terms of MAF values increase as the fault displacement increases. Concluding, the variability of the hazard curves highlights the fact that the selection of a set of fault scaling relations should be treated as an important epistemic uncertainty.

Sensitivity analysis
The proposed approach of PFDHA is intended for engineering applications, thus being applicable to a large number of essentially different faults at a regional level. Engineers must therefore become familiar with the effects of various parameters and aspects of the hazard calculations. To that effect, the impact of the main model parameters, namely b-value (Sect. 4.1.1) and maximum earthquake magnitude (Sect. 4.1.2), which shape the hazard curve is first investigated. Then, the effects of fault length (Sect. 4.1.3), fault mechanism (Sect. 4.1.4), and location of the crossing site (Sect. 4.1.5) on the hazard level are studied. Finally, the critical aspects of the calculation process, namely the correlation of RL and ADD (Sect. 4.2.1), the conditional probability of the rupture reaching the surface (Sect. 4.2.2), and the crossing point (Sect. 4.2.3) are discussed to provide insight into their significance in comparison to ground shaking hazard calculations that engineers are more familiar with. The sensitivity analysis presented here is not intended to assess the impact of all potential sources of uncertainty, but rather to demonstrate the effect of different parameters. When a case-specific detailed fault displacement hazard analysis is required, additional parameters would be considered by the experienced seismologist conducting the analysis. It is noted that the sensitivity analysis is also structure-independent and consequently there is no judgement on the relevant importance of each parameter on the lifeline structural response. Ultimately, this is a designated objective of a case-and site-specific study.
The sensitivity analysis was conducted for the lifeline-fault crossing case study shown in Sect. 2. The baseline parameters are summarized in Table 6 for the sake of clarity. Fig. 9 Fault displacement hazard curves on lifeline crossing site for normal, reverse, and strike-slip fault mechanism using alternative fault scaling relations Throughout the sensitivity analysis, selected parameters are changed to demonstrate their effect on the hazard calculations. Due to the fact that a single fault is being examined, the slip rate must be constant. Therefore, a constant value of the fault slip rate along strike of S = 0.5mm∕year is adopted, which falls within the range of the most frequently encountered values in the seismically active faults of EFSM20 (Basili et al. 2022). Noteworthily, at 0.5 mm/year an active fault would have accumulated ~ 6 m total offset over the Holocene (~ 12,000 years), or ~ 60 m over the Late Pleistocene (125,000 years), thereby creating geomorphic features that could be possibly recognized in the field. The recurrence rate has to be rederived after any adjustment to the parameters in the sensitivity analysis to ensure consistency with the constant slip rate along the fault, thus offering a common basis for comparing the different changes in parameters. The recurrence rate v is obtained after Youngs and Coppersmith (1986): where = 3 ⋅ 10 11 dyne∕cm 2 is the fault rigidity after Bilek and Lay (1999), A = LF ⋅ W is the fault area with the fault width W = 20km assumed to be constant, S is the fault slip rate, b is the value of the Gutenberg-Richter Law, parameter is defined as = bln10 , M 0 u is the maximum moment corresponding to the maximum magnitude M max of the fault.
The maximum moment is assumed to be deterministic and estimated as M 0 u = 10 d+cM max with parameters c = 1.5 and d = 16.1 (Hanks and Kanamori 1979). The first part of the sensitivity analysis (Sect. 4.1) deals with the input parameters, namely the b-value of the Gutenberg-Richter law, the maximum earthquake magnitude, the fault length, the fault mechanism, and the location of the crossing site, while in the second part (Sect. 4.2) critical aspects of the hazard calculations are discussed.

Slope of the Gutenberg-Richter Law (the b-value)
The effect of the b-value to the fault displacement hazard estimates is illustrated in Fig. 10, where it can be observed that for low fault displacement values the b-value has a relatively low effect, while for higher displacement values, the effect tends to be insignificant. The increase of the b-value increases the ratio of low versus high magnitude events (Fig. 2).

Maximum earthquake magnitude
Alternative values of maximum earthquake magnitude ( M max ) are typically considered in any seismic hazard study. The M max is estimated by empirical fault geometry-magnitude scaling equations. The median value, a higher or a lower value can be adopted, using the range of parameter a 4 of Table 5 after Leonard (2014) as per Table 7. The corresponding hazard curves are illustrated in Fig. 11

Fault length
The fault length is a parameter directly related to other important fault parameters such as fault width, fault area, and maximum earthquake magnitude. For example, strikeslip faults may extend up to hundreds of kilometers in length, while being shallow in terms of width, contrary to reverse faults. Earthquake magnitude has also been found to  be a linear function of the log 10 LF . Herein, the LF under examination is ranging from 50 to 150 km and the effect of fault length on the displacement hazard curve is shown in Fig. 12 for all fault mechanisms. Generally, increased fault lengths lead to higher MAFs. This mild trend appears more clearly in Fig. 13 for the normal fault, where small MAF increments for given fault displacements are observed even when the fault length doubles or triples.

Fault mechanism
A comparison of the hazard curves for different fault mechanisms, namely normal, reverse, and strike-slip is depicted in Fig. 14. A notable lower hazard level is observed for reverse mechanism compared to normal and strike-slip, which is mainly attributed to the conditional probability of slip [see Eqs. (18) and (19)] that mainly differentiates this mechanism from the normal and strike-slip, as explicitly outlined in Fig. 4.

Location of the crossing site
The crossing site location on the fault trace, as thoroughly discussed in previous sections, is critical to the fault displacement hazard level on the crossing site. Moss andRoss (2011) andYoungs et al. (2003) have normalized fault displacement data with the fault average surface displacement ( ADS ) to estimate the conditional probability that the offset at a specific point will exceed a predefined value. The distribution of ratio ∕ADS is expressed as a function of the crossing point XL with function f ( ∕ADS|XL) being symmetric about a maximum value at x∕RL = 0.50 . Accordingly, as shown in Fig. 15 and regardless of the fault mechanism, the hazard curves for alternative crossing points ( XL ) are parallel and the highest MAF values are obtained for XL = 0.50 because in this case the crossing site is located at the middle of the fault trace length. The practical outcome for consideration within a route selection process of a lifeline (Seel et al. 2014;Hamid-Mosaku et al. 2020) is that the engineer should expect higher fault displacements if the crossing site is decided to be closer to the middle of the mapped fault trace length.

Discussion on input parameters
To offer an overview of the sensitivity analysis, focusing on the input parameter to the hazard calculations, the effects of b-value, maximum earthquake magnitude ( M max ), fault length ( LF ), and crossing point ( XL ) are examined in terms of the obtained MAF for predefined values of fault displacement with respect to the baseline case (see Table 6). Figure 16 shows the variation of MAF with the baseline case being reported in the vertical line.
The primary observation is that the location of the crossing site on the fault trace affects significantly the resulting MAFs, regardless of the fault displacement value, as concluded in Sect. 4.1.5. This is of primary importance for engineers during the route selection  Table 6) for predefined values of fault displacement regarding the input parameters b-value, M max , LF , and XL . Along the vertical axis, values of the parameter investigated in each chart are increased procedure, and requires some fault-specific knowledge to be accurately accounted for. Then, regarding the two main seismological input parameters, namely the b-value and M max , the influence of the latter on the MAF is significantly higher, highlighting the need for appropriate treatment within a logic tree formulation to account for its epistemic uncertainty. Finally, as observed in both Figs. 12 and 13, the higher the fault displacement value examined, the lower the effect of LF is.

Correlation of RL and ADD
The term P RL j , ADD t |m i in Eq. (10) is introduced to account for the joint probability of RL and ADD conditioned the magnitude. In particular, the estimation of the covariance matrix [ ] in Eq. (14) requires the correlation coefficient log 10 RL,log 10 ADD of the residuals of RL and ADD conditioned on the earthquake magnitude. Unfortunately, this value is not provided in pertinent studies (e.g., Leonard 2014;Thingbaijam et al. 2017;Wells and Coppersmith 1994). Therefore, the effect of this factor is examined in terms of the resulting hazard curve by considering no correlation ( log 10 RL,log 10 ADD = 0 ), partial positive and negative correlation (± 50%) ( log 10 RL,log 10 ADD = ± 0.50 ), and full negative and positive (± 100%) correlation ( log 10 RL,log 10 ADD = ±1.00 ). At the same time, though, the value of this correlation coefficient cannot be accurately determined without recourse to the data employed for fitting each model. As the sensitivity analysis presented in Fig. 17 shows, the increase in the correlation leads to higher MAF values. Therefore, in the absence of any data, the 0% correlation option is adopted. This choice can be seen as reasonable in light of a conservative design approach; it can be empirically suggested that a negative correlation may be more appropriate, as a certain amount of energy (or given M ) may be released over a large length with a small slip or over a shorter length with a larger slip.

Fig. 17
Fault displacement hazard curves on lifeline crossing site by considering alternative correlations for the residuals of RL and ADD given the magnitude (correlation coefficient log 10 RL,log 10 ADD ) for the term P RL j , ADD t |m i of Eq. (10)

Conditional probability of the rupture reaching the surface
The Conditional Probability of Slip (CPS) is calculated using the empirical relation of Eqs. (18) and (19) for each fault mechanism. These relations are logistic regression models produced by the analysis of worldwide data (Wells and Coppersmith 1993;Moss and Ross 2011) to consider whether the rupture reaches the surface. The effect of the CPS on the fault displacement hazard curve is depicted in Fig. 18, where the hazard curve is obtained with and without considering the probability of the rupture reaching the surface. As expected, considering the CPS yields lower MAF values for the entire range of fault displacement values. It should be reminded that the CPS is an integral part of PFDHA and stands as a substantial difference from typical PSHA calculations because in the latter all nearby seismic sources are considered to contribute to the hazard at the site under examination.

Effect of crossing point
The effect of considering the crossing point is examined because contrary to what is anticipated in a PSHA where all potential ruptures contribute to the hazard integration at a given site, this is not the case for lifeline crossings. This effect is illustrated in Fig. 19 where the hazard curves are plotted with and without considering the location of the crossing point relative to each RL, in other words, by considering intercepting ruptures versus all ruptures. In the first case, if the lifeline is not intercepted by the ruptured segment, then Δ ( ) = 0 for that particular RL. In the second case, the hazard on the crossing site is computed assuming all RLs contain the crossing site regardless of the actual location of the RL on the fault trace, i.e., in any case Δ ( ) > 0.

Conclusions
Safeguarding the integrity and ensuring the seismic resilience of lifelines that cross active tectonic faults requires the reliable estimation of the fault displacement. To achieve this, one could employ an existing dataset of mapped faults and a set of empirical fault scaling relations. Even though it is a relatively easy computation, the resulting displacement comes with an unknown level of safety. On the other hand, the Probabilistic Fault Displacement Hazard Analysis (PFDHA) could be used to achieve a balance between safety and cost-effectiveness while adhering to the performance-based engineering framework. However, this requires an extensive effort in terms of input data and calculations, as well as significant experience. In order to find a middle solution and towards formulating a methodology for computing the fault displacement that can be incorporated in a design code, we present in this study a structure-independent simplified approach of PFDHA for engineering applications. Based on the established PFDHA of Youngs et al. (2003) we focus on the lifeline-fault crossing and we introduce appropriate simplifications and assumptions, allowing our approach to be applicable to different faults at a regional level. In any case, the approach is not intended to replace a case-and site-specific PFDHA accounting for all relevant uncertainties and specialized data if such an analysis is required by the lifeline owner/operator and/or the regulatory authorities. We use the presented approach in a next study to derive a set of simplified hazard-consistent and code-compatible expressions, suitable for use in the lifeline route selection, for the earthquake-resistant preliminary design of lifelines, and in cases where carrying out a detailed PFDHA may not be possible, for example, when there are no recent sediments to date past earthquakes.
At first, the calculation steps are presented in detail, based on the work of Melissianos et al. (2017a) and explanations are provided on the hazard calculation aspects. A thorough discussion is then offered on the available alternative empirical fault scaling relations that are required at certain critical steps of the hazard computations and how they affect the results. Moreover, a sensitivity analysis on the influence of input parameters, such as the Gutenberg-Richter b-value, maximum earthquake magnitude, fault length, and the location of the crossing site provides engineers with insight into the effect of parameters and aspects on the hazard level calculations. Generally, higher mean annual frequency values of exceeding low fault displacement values were estimated for short faults, while vice versa higher mean annual frequencies were obtained for high fault displacement values in case of longer faults. The location of the crossing site on the fault trace was identified as a critical parameter because the closer to the middle of the fault the site is, the higher the displacement hazard. The M max (maximum earthquake magnitude) is strongly correlated with the fault dimensions, and thus its effect on the fault displacement hazard curve is rather minimum. Additionally, the empirical fault-magnitude scaling equations that are incorporated in the hazard calculations were evaluated, indicating a strong variability in results due to the different background of each set of relations. Finally, the b-value plays a moderate role in terms that the higher the value, the lower the hazard because small magnitudes have a low impact on the fault displacement hazard.
In conclusion, the study offers a comprehensive exercise on how to calculate the fault displacement hazard for lifelines at fault crossings, also considering the uncertainties associated with seismological parameters. It provides practitioners with insights by demonstrating the impact of key parameters on both the hazard level and the shape of the hazard curve. In addition, it serves as a background document for future reference, detailing the effects of the simplifications and assumptions made, which are pertinent and relevant for the engineering community.