Development of a Software Tool for Seismic Probabilistic Safety Assessment Quantification with a Sufficiently Large Number of Bins for Enhanced Accuracy

: Quantification in seismic probabilistic safety assessment (PSA) includes the convolution of hazard and fragility curves. Because it is difficult to find the closed-form integration for the convolution of two curves, numerical methods are widely used in practice. In practical applications, the number of ground motion level bins in numerical methods is limited, and it is not clear whether the limited number of bins leads to conservative or optimistic results. In this study, the effect of the number of bins on the quantification results with simplified assumptions is investigated. It is found that the quantification results mostly decrease as the number of bins increases. To enhance accuracy in the quantification results of seismic PSA, a method and a software tool that enable a sufficiently large number of bins to be used for the quantification of seismic PSA models are developed. The application of the developed software tool to an example seismic PSA model demonstrates how the quantification results approach accurate results as the number of bins increases. The software tool developed in this study is expected to enhance the accuracy of seismic PSA quantification results.


Introduction
A seismic event is an extremely interesting event in the evaluation of the risk of a nuclear power plant.This is because a seismic event triggers initiating events and simultaneously damages systems, structures, and components (SSCs).Seismic events also destroy many redundant SSCs.Because of the characteristics of a seismic event, it is necessary to perform hazard and fragility analyses and reflect their curves in seismic probabilistic safety assessment (PSA).The curve derived from the hazard analysis is generally obtained from available seismic hazard data and ground motion prediction equations and presented as the exceedance frequency of a seismic event depending on the ground motion level, such as peak ground acceleration and spectral acceleration.Seismic hazard seismic analysis experts, seismic PSA is performed by PSA experts to quantify the seismic risk of a nuclear power plant.
Several methods can be used to estimate the seismic risk by combining the hazard and fragility curves.Because it is difficult to find the closed-form integration for the convolution of two curves, alternative seismic PSA quantification methods such as the Monte Carlo simulation method, Latin hypercube sampling method, and discrete probability distribution method [4] have been widely used.PRASSE [5] is a quantification code that implements these methods.The alternative seismic PSA quantification methods may provide accurate results in the estimation of the branch frequencies of the primary seismic event tree; however, they have limitations in preserving the logical links between the primary seismic event tree and secondary event trees.The quantification results may be inaccurate when such logical links between the primary and secondary event trees need to be properly considered.In addition, limited information would be included in minimal cutsets owing to the loss of logical links.
Recently, the approach of developing discrete hazard interval models has been gradually adopted in the PSA community and has given guidance for seismic PSA [6,7].The approach is basically a numerical method that discretizes the ground motion level into several bins and develops seismic PSA models for the bins.The approach can be implemented with conventional PSA modeling and quantification codes, augmented by quantification codes based on a binary decision diagram (BDD).
One of the most important advantages of the approach that uses discrete hazard interval models to estimate the seismic risk is the preservation of the logical links between the primary seismic event tree and secondary event trees.By preserving the logical links, the minimal cutsets can include necessary information on the combination of failures that leads to core damage.
By contrast, in manually developed discrete hazard interval models, a limited number of ground motion level bins has been utilized in practice.For example, three and eight ground motion level bins were used in the Risk Assessment of Operational Events Handbook [7] and Surry seismic PSA [8], respectively.The number of ground motion level bins has often been determined by analysts with insufficient backgrounds.In addition, it is not clear whether the quantified risks are conservative or optimistic with a limited number of ground motion level bins.
Section 2 investigates the change in the quantification results as the number of ground motion level bins increases.Section 3 describes the method and a software tool that enables a sufficiently large number of ground motion level bins to be used in the evaluation of the seismic risk and enables the preservation of the logical links between the primary seismic event tree and secondary event trees to enhance the accuracy of the quantification results.In Section 4, the software tool is applied to an example seismic PSA model to demonstrate the level of accuracy in the quantification results with an increase in the number of ground motion level bins.Section 5 presents the conclusions of this study.

Effect of Number of Bins on the Quantification Results
The quantification of seismic risks involves the convolution of hazard and fragility curves, as shown in Equation (1) [9][10][11] where is the ground motion level, ( ) is the exceedance frequency at ground motion level , and ( ) is the seismic failure probability at ground motion level .
The approach that uses discrete hazard interval models approximates Equation (1) with Equation (2) by discretizing the hazard and fragility curves into several ground mo-tion level bins.When the ground motion level is discretized into , … , , the approximation can be obtained by summing the products of the difference of hazard between and , that is, ( ) − ( ), and the failure probability at the representative ground motion level, as follows: where is the number of ground motion level bins, and , is the representative ground motion level for the i-th bin.
The hazard curve is a decreasing curve and is assumed to be approximated by the power law, as shown in Equation (3) [9][10][11]: where is an appropriate constant, and is the slope parameter.The hazard curve is used to assign the frequency of the seismic initiating event in each bin.
The fragility curve is lognormally distributed and is given by Equation ( 4): where is the median capacity, and is the composite variability, which is defined by + .and are the two uncertainty variabilities for the inherent randomness and the modeling uncertainty, respectively.The composite curve is the mean curve [12].The fragility curve is used to assign the seismic failure probability of the basic event in a seismic PSA model.
To investigate the effect of an increase in the number of ground motion level bins on the quantification results, a case is examined in which the number of ground motion level bins is increased by two times.For any ground motion level interval from a to e, as shown in Figure 1, two cases are examined depending on the number of bins for the ground motion level interval.The first case has a single bin of the ground motion level from a to e.The failure frequency for the ground motion level interval is given by Equation (5).
The second case has two bins of the ground motion level, one from a to c and the other from c to e.The failure frequency is given by Equation (6).
In this analysis, the representative ground motion level for calculating the seismic failure probability of the SSC from its fragility curve is assumed to be the arithmetic mean between the minimum and maximum ground motion levels of the bin.To compare the results, Equation ( 7) is obtained by subtracting Equation (6) from Equation (5).
The sum of the hazard terms in Equation (6), that is, [ ( ) − ( )] + [ ( ) − ( )], is equal to the hazard term in Equation (5), that is, ( ) − ( ).The second fragility term in Equation (7), that is, ( ) − ( ), is negative.If Equation (7) provides a positive result, the quantification result of the single-bin case is greater than that of the two-bin case, thereby implying that the quantification result may decrease as the number of ground motion level bins increases.By contrast, if Equation (7) provides a negative result, the quantification result of the two-bin case is greater than that of the single-bin case, thereby implying that the quantification result may increase as the number of ground motion level bins increases.By comparing the values of the terms in Equation ( 7), the failure frequencies of the single-bin and two-bin cases can be compared.Because a hazard curve is a decreasing curve and its slope increases as the ground motion level increases, the first hazard term of Equation ( 7) is greater than the second hazard term of the equation, that is, For the fragility terms in Equation ( 7), whether the absolute value of the first fragility term, that is, ( ) − ( ), is greater than that of the second fragility term, that is, ( ) − ( ), depends on the change in the fragility curve's slope.For a single SSC, the fragility curve is given by Equation (4).However, for a system with multiple SSCs, the analytic form of the fragility curve cannot be readily obtained.The fragility curve of a system is divided into three ground motion level intervals.
When the change in the fragility curve's slope is negative, that is, when the fragility curve's slope is decreasing, the first fragility term is greater than the second fragility term, that is, In this case, it is determined that Equation ( 7) is always positive, that is, the failure frequency of the ground motion level interval calculated with a single bin is always greater than that with two bins.
When the change in the fragility curve's slope is close to zero, that is, when the fragility curve is almost linear, the first and second fragility terms do not differ significantly, that is, In this case, because the first hazard term is greater than the second hazard term, it is expected that Equation ( 7) becomes positive, that is, the failure frequency of the ground motion level interval calculated with a single bin is always greater than that with two bins.
When the change in the fragility curve's slope is positive, that is, when the fragility curve's slope is increasing, the first fragility term is less than the second fragility term, that is, Because the first hazard term is greater than the second hazard term, whereas the first fragility term is less than the second fragility term, as shown in Equations ( 8) and (11), whether the result obtained from Equation ( 7) is positive or negative depends on the magnitudes of the four terms.The hazard term decreases exponentially following the power law shown in Equation ( 3), and the result obtained from Equation ( 7) becomes negative only when the first fragility term is significantly less than the second fragility term, which usually occurs when the ground motion level is significantly less than the median capacity ( ).
Therefore, in most of the cases, it can be concluded that the failure frequency of the ground motion level interval calculated with a single bin is greater than that with two bins, except when the ground motion level is significantly less than the median capacity ( ).In a similar manner, the conclusion can be generalized when the ground motion level interval is divided into more than two bins.

Quantification with a Sufficiently Large Number of Bins
For each ground motion level bin, the seismic initiating event frequency and seismic failure probabilities of SSCs need to be calculated to develop the seismic PSA model for the bin.Therefore, the number of seismic PSA models and quantifications to be the same as the number of ground motion level bins is needed.It may require a considerable amount of time and effort to analyze the seismic risks as the number of ground motion level bins increases.Therefore, a software tool is developed to assist seismic PSA quantification with a sufficiently large number of ground motion level bins by automatically generating seismic PSA models for the bins and performing repeated quantifications.
The seismic initiating event frequencies of the ground motion level intervals are calculated from the hazard curve.A hazard curve is generally presented as the exceedance frequency of a ground motion level.When seismic hazard data are provided as point values to represent the hazard curve, linear interpolation is performed on a semi-log scale.The difference in the interpolated hazard curve is calculated to derive the seismic initiating event frequency for each ground motion level bin.
The basic events for the seismic failures of SSCs are added to the seismic PSA models in each ground motion level bin.The seismic failure probabilities of SSCs are dependent on the representative ground motion level of the bin and are calculated using Equation (4).As the ground motion level increases, the seismic failure probabilities of SSCs also increase.
Human error probabilities (HEPs) are also affected by the ground motion level.Several methods can be used for human reliability analysis under seismic events.The Swiss Federal Nuclear Safety Inspectorate (ENSI) [13] assumes that the human failure probabilities are not significantly affected when the ground motion levels are less than 0.2 g, increase linearly for ground motion levels from 0.2 g to 0.6 g, and are equal to one when the ground motion levels are more than 0.6 g.The seismic human reliability analysis method for a Surry nuclear power plant, described in EPRI TR 1020756 [8], utilizes multipliers depending on the ground motion level and the time after the seismic event.
During the development of the seismic PSA model for a nuclear power plant, separate models for the ground motion level bins are developed, and the quantification results for all the bins are added to calculate the quantification result for the entire ground motion level considered for the seismic PSA.
A software tool implementing the aforementioned method is developed to assist the seismic PSA quantification.Figure 2 compares the existing method and the novel method utilizing the software tool when the seismic PSA model is developed using AIMS-PSA.AIMS-PSA [14] is one of the most widely used PSA modeling and quantification packages in Korea.In the existing method, the ground motion level from 0.1 g to 1.0 g is divided into nine bins, and separate seismic models for the bins are developed, usually manually, to quantify the seismic risks for the bins.
The software tool utilizes a script language called SIMA (Script Interpreter for Mapping Algorithm) and the sensitivity analysis function of AIMS-PSA to perform quantification with a sufficiently large number of bins.The SIMA script enables analysts to add or change the logic and basic event probabilities in the single top fault tree model, which combines all the event tree and fault tree models.By using the SIMA script, the seismic PSA model for each bin can be generated from the base model by reflecting different values for the initiating event frequency, the seismic failure probabilities of SSCs, and HEPs under seismic events.Similar script languages or features can be found in other PSA modeling and quantification packages.
The software tool performs a number of quantifications for the generated seismic PSA models for the ground motion level bins using the sensitivity analysis function of AIMS-PSA.The sensitivity analysis function is originally provided for sensitivity analysis by changing the values of the variables.In Figure 2, the ground motion level is divided into 90 bins, and SIMA scripts for the 90 bins are provided.It is necessary to generate the SIMA script for each bin.Therefore, the software tool automatically creates the SIMA script file for each bin.Moreover, the software tool generates input that can be copied to the window shown in Figure 2, which facilitates automatic quantification for a sufficiently large number of bins.
Figures 3 and 4 show the user interface of the software tool.It is assumed that the hazard analysis and fragility analysis are finished and the hazard data points and seismic fragility data of SSCs are provided by seismic analysis experts.Because seismic PSA is performed based on given seismic hazard curve and fragility curves, the application of the proposed method and software tool is not affected by how the hazard and fragility curves are obtained.The software tool receives the hazard data points, user-defined ground motion level intervals, and seismic fragility data of SSCs and base HEPs.The seismic hazard curve is constructed with hazard data points by applying linear interpolation on a semi-log scale.Ground motion level bins are defined according to the user-defined ground motion level intervals.The seismic fragility data of the SSCs include the median capacity and the two uncertainty variabilities: one for the inherent randomness ( ) and the other for modeling uncertainty ( ).Seismic failure probabilities of the SSCs are calculated based on the seismic fragility data and the representative ground motion level of the bins.HEPs for the bins are calculated with the base HEPs and user-defined rules for adjusting HEPs under seismic events.Figure 4 shows the case where ENSI guidance is applied.
Based on the calculated hazard curve, seismic failure probabilities of SSCs, and HEPs under seismic events, the SIMA script for each bin is automatically prepared for use, as shown in Figure 2. The SIMA script is used to modify the base seismic PSA model to the seismic PSA model for the bin.Each seismic PSA model for a bin, reflecting the initiating event frequency and seismic failure probabilities of SSCs and HEPs under seismic events for the bin, is analyzed with fault tree solvers such as FTREX [15] to identify the minimal cutsets.When the single top fault tree model for each bin includes the logical links between the primary seismic event tree and secondary event trees (which include associated fault trees) for the bin, the logical links between the primary seismic event tree and secondary event trees are preserved, and the ground motion level is divided into a sufficiently large number of seismic bins.This is an important advantage over the widely used methods based on Monte Carlo simulation or Latin hypercube sampling, in which the logical links between the primary seismic event tree and secondary event trees are not properly preserved.
Because high-failure probability events are involved in the primary seismic event tree, it is important to preserve the success branches in a primary seismic event tree properly while generating the minimal cutsets for each bin.Kim and Kim [16] present the lessons learned with the effects of preserving the success branches in a primary seismic event tree during the quantification of plant-level seismic PSA models.
After generating minimal cutsets for all the bins, BDD-based quantification [17] is then applied to the minimal cutsets of each bin to provide the seismic risk quantification result for the bin.Currently, the application of BDD-based quantification to minimal cutsets needs to be performed outside the AIMS-PSA environment.A separate command line script is used for repeated application of BDD-based quantification to the minimal cutsets of all bins.
The software tool enables practical analysis of the seismic risk with a sufficiently large number of bins to enhance the accuracy of the quantification results; moreover, it significantly reduces the work involved in manual modeling while preserving the logical links between the primary event tree and secondary event trees.

Application to an Example Seismic PSA Model
To develop an example seismic PSA model, an example internal event PSA model is necessary.The pilot plant PSA model is a simplified PSA model based on the PSA model for optimized power reactor 1000, developed by the Korea Atomic Energy Research Institute for research and educational purposes.Thirteen initiating events are considered, including large, medium, and small loss of coolant accidents; loss of component cooling water; general transients (GTRNs); loss of offsite power (LOOP); and station blackout.Important front-line systems and supporting systems are included, such as a high-pressure safety injection system, low-pressure safety injection system, auxiliary feedwater system, containment spray system, component cooling water system, essential service water (ESW) system, and electrical systems.The one-top (single top) model includes 465 events, including initiating events, random failures of SSCs, human failure events, and common cause failures.With a cut-off value of 10 -14 per year, 5466 minimal cutsets are identified, and the core damage frequency is estimated as 8.432 × 10 -6 per year with rare event approximation.The pilot plant PSA model can be used as a common research platform to discuss new ideas in the PSA field.
The example seismic PSA model is developed based on the hazard and seismic fragility data provided by Park et al. [18] with generic fragility data for LOOP from EPRI TR 3002000709 [6].The seismic hazard data for the example seismic PSA model are provided in Table 1.
A primary seismic event tree is shown in Figure 5.The seismic fragility data for the primary seismic event tree are provided in Table 2.In Sequence 4, the occurrence of any of the three events leads to the failure branch.The sequences for GTRN and LOOP, that is, Sequences 1 and 2, are connected to corresponding secondary event trees that are adopted from the pilot plant PSA model for internal events with additional seismic failures of SSCs.
Secondary event trees for GTRN and LOOP include both random failures and seismic failures.In the example model, only the seismic failures of emergency diesel generators (EDGs) and ESW pumps (ESWPs) are included in the secondary event trees for the purpose of demonstration.Other seismic failures of SSCs may be added if necessary.The seismic fragility data for the secondary event trees are provided in Table 3. Owing to the existence of seismic failures in the secondary event trees, quantification of the example seismic PSA model needs to consider the logical links between the primary and secondary event trees.To understand the relation between the quantification results and the number of ground motion level bins, and to demonstrate the usefulness of the developed software tool, the quantification of the example seismic PSA model is performed while increasing the number of bins.
Table 4 illustrates how the number of bins is increased in this study.In Park et al. [18], the ground motion level interval that is greater than 0.1 g is divided into six bins.By further dividing each of the five ground motion level intervals between 0.1 g and 1.0 g into bins, the number of bins yields the following form: 5 + 1.For example, when = 2, that is, each of the five ground motion level intervals between 0.1 g and 1.0 g is divided into two bins, the number of bins becomes 11, as shown in Table 4.By applying = 1, 2, 3, 4, ⋯ , 20, the number of seismic bins becomes 6, 11, 16, 21, and up to 101.
For the quantification of seismic risks, minimal cutsets are generated with FTREX [15] while preserving the success branches in the primary seismic event tree shown in Figure 5.The truncation limit (cutoff value) of 10 -14 is applied to all bins.Table 5 lists the top five minimal cutsets for the bin with the ground motion level interval of 0.3 g-0.5 g, which has the largest failure frequency when the ground motion level is divided into six bins.The five minimal cutsets take up to 99% of the failure frequency of the bin.In Table 5, '/' in front of an event name denotes the NOT logic, that is, the success of the event, and it can observed that success events are well preserved while generating minimal cutsets.Whereas the first, second, third, and fifth minimal cutsets include only seismic failures in the primary seismic event tree, the fourth minimal cutset includes a seismic failure in the primary seismic event tree, a random failure, and a seismic failure in the secondary event trees.The failure information in the primary seismic event tree is logically linked to the failure information in the secondary event trees in this approach, whereas the conventional approaches based on Monte Carlo simulation or Latin hypercube sampling cannot provide such logical links.As the number of bins increases, the quantification results become more accurate.The sequence frequencies of Sequences 3 and 4 in Figure 5 calculated with the support of the software tool and the use of 101 seismic bins are 2.11 × 10 -6 and 5.65 × 10 -6 per year, respectively, where those calculated with an alternative seismic PSA quantification method are 2.11 × 10 -6 and 5.66 × 10 -6 per year, respectively.When 101 seismic bins are used for quantification, the core damage frequency is calculated to be 8.779 × 10 -6 per year.The error in Figure 6 is defined as the difference between the quantification results for the number of seismic bins and 101 seismic bins divided by the quantification result for 101 seismic bins.Figure 6 depicts how the errors in quantification results change as the number of bins increases.From the error of 10.7% when six bins are used, the error gradually decreases and approaches zero as the number of bins increases.Figure 6 quantitatively shows how the quantification results approach the accurate quantification results as the number of bins increases.If more seismic failures of SSCs are considered, the errors in case of small numbers of bins are expected to be greater than those shown in Figure 6, where the errors related to the seismic failures of only seven SSCs are considered.As mentioned in Section 1, three and eight seismic motion level bins were used in Risk Assessment of Operational Events Handbook [7] and Surry seismic PSA [8], respectively.Even though direct comparison would not be possible, Figure 6 provides implications on possible errors caused by small number of seismic motion level bins.By referring to the 10.7% error when six bins are used in Figure 6, it is expected that similar amount of errors can be reduced if sufficiently large number of bins is used by using the method and software tool described in this paper.
The newly developed software tool makes analysis of the change in the quantification result with the increase in the number of bins possible, as shown in Figure 6.The use of a sufficiently large number of bins is recommended for more accurate quantification of seismic PSA models, unless the use of fewer seismic bins can be justified for a specific seismic PSA model.

Conclusions
In seismic PSA modeling and quantification, the approach with discrete hazard interval models has been gradually adopted over the widely used methods based on Monte Carlo simulation or Latin hypercube sampling.As only a limited number of ground motion level bins has been utilized in practice, it needs to be determined whether the quantified risks are conservative or optimistic with a limited number of bins.
In this study, it is investigated whether the failure frequency increases or decreases as the number of ground motion level bins increases.For two cases of a single bin and two bins for the same ground motion level interval, the failure frequency calculated with a single bin was generally greater than that calculated with two bins, except when the ground motion level was significantly less than the median capacity.This conclusion implied that the calculated core damage frequency mostly decreases as the number of bins increases, except that the median capacities of SSCs are significantly greater than the ground motion level interval of interest.
For more accurate quantification of seismic risks, a software tool is developed that supports the use of a sufficiently large number of bins in seismic PSA modeling and quantification by automatically generating the seismic PSA models for the bins, identifying minimal cutsets, applying a BDD-based quantification, and aggregating the quantification results.Through application to an example seismic PSA model, the software tool demonstrated how the quantification results approached the exact result as the number of bins increased.Because the error in the quantification results was larger than 10% when six bins were used for the simple example seismic PSA model with only seven seismic failures of SSCs, higher quantification errors are expected with the use of a small number of bins as the seismic PSA models for nuclear power plants include more and more seismic failures of SSCs.The proposed method and software tool are also applicable to different kind of earthquakes and installations as far as the hazard curve and fragility curves for SSCs are provided.The proposed method and software tool have the strength of enhancing the accuracy of seismic risk quantification while linearly increased time will be required for the quantification as the number of seismic bins increases.Future research may include how to reduce the required time for quantification as a sufficiently large number of seismic bins are used.
The relation between the number of ground motion level bins and the seismic risk quantification results is expected to expand our understanding of seismic risk quantification.The software tool for supporting seismic PSA modeling and quantification by allowing a sufficiently large number of bins is expected to contribute to a more accurate estimation of the seismic risks of nuclear power plants.

Figure 2 .
Figure 2. Comparison between the existing model and the new model with automatic generation and quantification.

Figure 3 .
Figure 3. Input to the software tool.

Figure 4 .
Figure 4. Calculation of the hazard curve, seismic failure probabilities of systems, structures, and components (SSCs) and human error probabilities (HEPs).

Figure 5 .
Figure 5.Primary seismic event tree in the example, General Transients (GTRN), Loss of Offsite Power (LOOP).

Figure 6 .
Figure 6.Errors in quantification results as a function of the number of seismic bins.

Table 1 .
Seismic hazard data for the example model.

Table 2 .
Fragility data for the primary seismic event tree.

Table 3 .
Fragility data for the secondary seismic event trees.