One dimensional seismic response analysis at the non-commercial nuclear reactor site, Serpong-Indonesia

One dimensional seismic response analysis on the ground surface of the Non-Commercial Power Reactor (RDNK) site based on the mean uniform hazard spectrum (UHS) and disaggregation analysis has been conducted. The study’s objective was to perform an analysis on site-specific response spectra on the ground surface based on existing mean UHS and disaggregation data of the site that correspond to a 1,000 and 10,000 year return period of earthquakes in compliance with the national nuclear regulatory body requirements of Indonesia. Detailed site characterization was defined based on secondary data of a geotechnical drill-hole, seismic cross-hole, downhole data, and microtremor array data. The dynamic site characteristic analysis was presented along with strong motion selection and processing using two types of strong motion datasets. An investigation of strong motion selection, spectral matching, and scaling has been presented as an essential step in ground motion processing. One-dimensional equivalent linear analysis simulation was performed by computing the processed ground motions. A seismic design spectrum and ground surface response spectra from the two datasets of strong motion, both corresponding to a 10,000 and 1,000 year return period, are presented at the end of this study. This study has shown that in order to establish the appropriate seismic response design spectrum, site-specific data and seismic hazard analysis must be immensely considered.


Introduction
In 2015, the National Nuclear Energy Agency (Badan Tenaga Nuklir Nasional, BATAN) had a plan to build a non-commercial nuclear power reactor (Reaktor Daya Non-Komersial, RDNK) based on a pebble-bed reactor high-temperature gas-cooled reactor (Setiadipura, T. et al., 2018). The selected site is located near (less than 1 km distance) the Serpong Nuclear Complex inside the Center for Research, Science, and Technology (Pusat Pengembangan Ilmu Pengetahuan dan Teknologi, PUS-PIPTEK) area. There are numerous and complex criteria for selecting a site for a nuclear reactor and other related facilities. For instance, to choose the location for the disposal of spent nuclear fuel (SNF), the site must comply with the geological, hydrogeological, hydrological, and other technical requirements (Veinović, Z. et al., 2019). The criteria list for the nuclear reactor itself is longer than the SNF repository. At least seven aspects should be considered during site selection: seismic, volcanology, geotechnical, meteorology, hydrological, dispersion, and human-induced events.
Seismic hazard at the site became the primary concern in investigating the planned reactor's natural external hazard of the planned reactor because it is located within a relatively moderate seismic activity region and surrounded by a relatively dense population. Based on the 1900-2016 earthquake compilation catalogs (Incorporated Research Institutions for Seismology-IRIS; Meteorological, Climatological, and Geophysical Agency of Indonesia (Badan Meteorologi, Klimatologi, dan Geofisika, BMKG), and United States Geological Survey-USGS, there were at least four moderate scale earthquakes with a magnitude of 4-6Mw recorded within a 25 km radius from the site. The depth of these earthquakes was within 60 km down to 300 km. Probabilistic Seismic Hazard Analysis (PSHA) was performed within a 300 km radius from the site, as shown in Figure 1. Data for PSHA at the site shows a mean Uniform Hazard Spectrum (UHS) at the bedrock level was acquired from a previous study (Badan Tenaga Nuklir Nasional, 2016). The study also presented the disaggregation analysis corresponding to a 10,000 year return period, which complies with the national nuclear regulatory body regulation (Badan Tenaga Nuklir Nasional, 2016).
Stability of a structural system relies upon the shapes, sizes and the performance of the materials used (Ja guljnjak Lazarević, A. et al., 2017). However, structural engineering design requires a comprehensive analysis of how earthquake characteristics on the ground surface affected the structure. It is well established that local soil conditions significantly affect the amplitude and response spectral characteristics of ground motions.

Suntoko, H. et al. (2019)
performed an analysis of the response design spectrum of RDNK reactor building. However, the study was done using only the mapped spectral response based on SNI 1726-2012 without considering the site-specific soil characteristics, PSHA and disaggregation analysis for the RDNK site. Indonesian National Standard (Standar Nasional Indonesia, SNI) 1726-2012 is the national standard that presents the design and construction criteria for general structures or buildings and could not be applied for specific or critical structures such as nuclear reactor installation. Thus, it was necessary to perform a site-specific investigation of design earthquake characteristics. This study's aim was to perform analyses on site-specific response spectra on the ground surface based on existing mean UHS and disaggregation data of the site that corresponds to a 1,000 and 10,000 year earthquake return period. The soil model for the equivalent linear analysis was developed using geological, geophysical, and geotechnical drill data. Previous study had performed preliminary analysis on spectral matching analysis at the site by using only Chi-Chi 1999 strong-motion dataset (Yuliastuti, et al., 2021). Current analysis uses two datasets of strong motion, later called datasets A and B, which have different selection criteria. Dataset A was selected only by magnitude and distance, while dataset B was selected based on fault mechanism, magnitude, and distance.

Methods
In the following, the analysis of dynamic site characteristics is presented along with strong motion selection and processing using two strong motion datasets. An investigation of recorded strong motion selection, spectral matching, and scaling has been essential in strong-motion processing. One-dimensional equivalent linear analysis simulation was performed by computing the processed strong-motions, or later called input motions.

Generic soil model development
The RDNK site was situated between Banten Blok and southwest Java Basin (Syaeful and Muhammad, 2017). showed that the lithological unit at the RDNK site includes limestone as a part of Bojongmanik Fm., tuff sandstone as part of Genteng Fm., conglomerate sedimentary rock as part of Serpong Fm., tuffs deposit, and alluvial sedimentary deposits (see Figure 2). In developing ground surface response spectra, dynamic site analysis was performed using geotechnical and geophysical data, as shown in Figure 3. Geotechnical data used in this study consists of Standard Penetration Tests (SPT) (27 metre depth), a seismic cross-hole (100 metre depth), and downhole (50 metre depth) test data. Seed et al. (1986), as well as Ohta and Goto (1978) utilized empirical correlations of shear-wave velocity (V s ) and standard penetration test blow count (SPT), to acquire V s from borehole data. Simultaneously, the geophysical data in the form of microtremor array was also used to cover a deeper soil profile. The shearwave velocity profile for the soil model was developed based on the averaging approach. The soil model was generated based on the averaged V s . Meanwhile, the lithology type for each layer was defined based on the geological coring and SPT results. A density log was used to determine the unit weight for each layer. From top to bottom, the soil model consists of a Serpong formation (topsoil, clay to silt, clayey to silty sand, and silty to clayey sand) and a Bojongmanik formation (Sandy claystone 3, Sandy claystone 2, Clayey sandstone, Sandy claystone 1, and Claystone). Topsoil in this model is the top layer of ground consisting of clay and organic materials. A generic one-dimensional soil model lithology column was generated consisting of nine soil layers, as shown in Figure 4.
A generic one-dimensional soil model, was then loaded into the DeepSoil code program (Hashash et al., 2016) as such to produce a frequency range of 30-79 Hz. At least three input parameters were required: thickness, unit weight, and shear wave velocity for each layer of the soil model. The highest frequency was identified as a very thin layer of claystone at a depth of about 45 metres. The effective stress of the soil model tends to increase with depth, which also increases the rate of soil compression. A relatively high effective vertical stress was identified below 45 metres of depth, reaching 2000 KPa at the bottom of the soil model. The high-level of effective stress induced a high level of maximum strain as well.
The highest frequency is the maximum frequency that the layer can propagate and is calculated as given in Equation 1 (Budhu, 2010). Meanwhile, total vertical stress (σ v ) at depth, z, is shown in Equation 2 (Budhu, 2010). Effective vertical stress is equal to the total stress minus pore water pressure. Effective vertical stress is displayed in Figure 5, where the value of effective vertical stress is at the soil layer's mid-depth.
(1) Where: V s -the shear wave velocity of the layer (m/s), d -the thickness of the layer (m).
Where: γ 1 , γ 2 , γ n -the unit weight for each soil layer. The input parameters required to perform site response analysis on the ground surface are damping and modulus reduction curves for each layer of the soil profile. The selection of the curves was based on the soil and rock characteristics of each layer. Seed et al. (1986), Vucetic and Dobry (1991), and EPRI 250-500 (Pyke and North, 1973) modulus and damping curves were applied to each of the soil layer based on the lithology class and secondary data of soil laboratory test. Figure 6 shows the material modulus and damping reference curves.

Input motion selection and scaling
Seismic hazard disaggregation analysis is a process to determine the dominant earthquake Magnitude (Mw) and rupture distance (R), which significantly contributes to the probabilistic seismic hazard analysis of the site at a specific return period. The disaggregation analysis of the RDNK site resulting in a dominant hazard coming from the subduction regime (Normal fault mechanism) with a mean magnitude of 7.88 Mw and a mean distance of 122.28 km (Badan Tenaga Nuklir Nasional, 2016). Based on this analysis, strong motions were selected.
A dataset of strong-motion data was acquired through several strong-motion databases such as NGA West (Ancheta et al., 2014), Center for Engineering Strong Motion Data (CESMD), and K-NET. All selected strong motions were located on the outcrop and consist of horizontal and vertical components. Input motion spectral matching and scaling were performed only for the horizontal component. In this study, a comparison of two types of strong motion selection were presented. The first type of selection, later called A, was conducted by considering all kinds of tectonic regimes, regardless of the fault mechanism, within the range of magnitude and distance based on disaggregation analysis, as mentioned before. Table 1 lists the selected strong motions based on this type of selection, consisting of active shallow crustal and megathrust/subduction interface tectonic regimes. The second type, later called B, was to select strong motions strictly based on fault mechanism (Normal/Reverse) and magnitude and distance following the disaggregation analysis result. Five strong motions have been selected based on this type of selection, as shown in Table 2.
Despite the out of range magnitude, the Tohoku earthquake (2011) was also utilized as the input motions for both of the ground motion selection types to consider the possibility of a 9 Mw earthquake at the Java megathrust as stated in the Indonesian National Standard 03-1726-2002. The 9.1 Mw Tohoku earthquake was one of the earthquake events that happened on the subduction interface.
Spectral matching and scaling processes were applied to each of the strong motions by preserving the non-stationary properties of the original accelerogram using Al Atik & Abrahamson algorithm available on Seismomatch 2020. The result of spectral matching and scaling for ground motions dataset A and B correspond to a 10,000 year and 1,000 year return period (see Figure 7).
As shown in Table 2, the ground motion dataset B consisted of the El Salvador (2001) earthquake due to the compatibility of its magnitude, distance, and source mechanism. Nevertheless, this particular earthquake generated a high rise to the spectral matching misfit reaching a 192% maximum misfit. Thus, the El Salvador (2001) earthquake was eliminated from the ground motion dataset B. Baseline correction was applied to all the selected input motions to avoid time series drifting. As stated in Boore (1999), a good baseline correction to maintain the actual real signal was required. Amplitude shifting at low and high frequency was noticeably due to the baseline correction process. Baseline correction was applied to all the selected input motions to avoid time series drifting. Figure 7 shows one example of a baseline correction application on Kepulauan Mentawai (2007) ground motion. The blue line shows the original motion data, while the purple line shows the baseline-corrected data. Fourier spectra comparison before and after baseline correction yielded 2834 correction points. Amplitude shifting at low and high frequency was noticeable due to the baseline correction process.

Results
One-dimensional total stress analysis using equivalent linear analysis in the frequency domain was con-ducted using DeepSoil codes (Hashash et al., 2016). Baseline corrected strong motions were utilized as the input motions for the equivalent linear analysis. Based on the input motion type used in this study, elastic halfspace was utilized to characterize the bedrock. The bedrock properties included at least three parameters, namely Vs, unit weight, and damping ratio. According to the geophysical data, Vs of 750 m/s, unit weight of 18.36 kN/m3, and 5% damping ratio was assigned as the bedrock properties.
Iteration for the analysis was performed fifteen times to get the convergence value. The input enabled the program to calculate initial shear modulus by correlating shear wave velocity and unit weight. Iterated shear modulus was calculated based on the frequency-dependent complex shear modulus formulation, which relates the shear modulus initial to shear modulus iterated by using the damping ratio (Hashash et al., 2016). The initial damping ratio was calculated based on selecting the damping curve, as mentioned in the previous chapter of dynamic site characteristics.
Spectral acceleration for each input motion has been plotted over the seismic design spectrum in Figure 9,  The lower table shows the average and maximum misfit resulted from the spectral matching process. Maximum acceleration is the maximum calculated acceleration from the selected strong motion. both corresponding to a 10,000 and 1,000 year return period. The design spectrum was developed following the Indonesian National Standard 03-1726-2002 on "Procedures on Earthquake Resistance Planning for Buildings" and ASCE 7-10. The seismic design spectrum parameters are shown in Table 3.    Figure 10 shows the mean spectra comparison of the top of rock layer and surface layer. The top of rock is defined at the bottom of soil model (Claystone with a unit weight of 16.08 kN/m 3 and V s of 698 m/s). The surface layer is defined as the topsoil layer (Clay and organic material with a unit weight of 15.69 kN/m 3 and V s of 229 m/s. The spectra were built based on pseudospectral acceleration data for each layer when subjected to each input motion corresponding to both, a 10,000 and 1,000 year return period. The soil amplification factor corresponds to a 10,000 year return period for dataset A, and B, and yields 1.65 and 1.67, respectively. Meanwhile, for a 1,000 year return period, the soil amplification factors were 1.69 and 1.68, respectively, for datasets A and B.

Discussion
For a 10,000 year return period, maximum pseudospectral acceleration of the input motions at a short period range (T = 0.2 -1 sec) was contributed by the Mentawai, Indonesia (2007) input motion, which yields 1.7678g. This particular earthquake is categorized as a shallow earthquake caused by a thrust faulting on the boundary plate of Indo-Australia and the Sunda plate (Hayes, G.P., et al., 2017). Meanwhile, at a longer-period (T > 1.5 sec), Tohoku (2011) input motion has a clear offset from the rest of the datasets. This phenomenon is because the large magnitude and distant earthquake(s) tend to have a moderate response spectra peak at a short period and a much higher peak at a longer period compare to the moderate magnitude and close earthquake. For both a 10,000 and 1,000 year return period, the mean response spectra of dataset B yield a higher value than strong-motion dataset A, resulting in a slightly different seismic design spectrum parameter. UHS at bedrock was utilized to determine some of the seismic design spectrum parameters, especially spectral acceleration at a period of 0.2 sec (S s ) and 1 sec (S 1 ). Spectral response acceleration at the ground surface for a 1,000 and 10,000 year return period as displayed in Figure 9 shows that at short period ranges and a 1 sec period, dataset A gives a 6% lower value compare to dataset B.
A previous study by  Table 3, it gives a 22-32% higher value in terms of the S DS coefficient, and 33-49% for the S S coefficient. Whereas, the S 1 coefficient tends to be 1-3% lower than Suntoko, H. et al. (2019). The different seismic design response coefficients create a significant impact to the civil design of the reactor and other buildings in the site area.
Seismic hazard analysis of a building structure cannot be carried out simply based on the maximum value of ground acceleration. Ground motion frequency content and dynamic soil properties influence the kind of response generated by the structure. The mean spectra on the ground surface reveals that the probability of ground motion exceedance level produces a significant difference to the resulted seismic response spectra.

Conclusions
Regional-scale seismic hazard analysis and site-specific data are essential to acquire a reliable seismic response design spectrum at the ground surface. Mean UHS resulted from a regional seismic hazard study has been utilized in this study. The mean UHS was assigned as the target spectrum when performing spectral matching and scaling. Eleven recorded strong-motions were selected based on the disaggregation analysis at the site. The calculation shows that, at the ground surface, peak ground acceleration (PGA) for a 10,000 and 1,000 year return period are in the range of 0.59-0.61g and 0.38-0.4g, respectively.
A significant finding to emerge from this study is that the generated design response spectrum differs from previous study results (Suntoko H., et. al 2019), which did not include the soil profile and detailed site seismic hazard analysis. The utilization of seismic hazard and site-specific data give a 22-32% higher value in terms of the S DS coefficient, 33-49% for the S S coefficient, and 1-3% lower of the S 1 coefficient compared with Suntoko (2019).
Besides, the selection criteria of input motion produce a different result of the ground surface's spectral response. Dataset A is comprised of eight selected strong-motions without considering the source mechanism or fault type. Meanwhile, dataset B consists of five strongmotions selected based on fault type, magnitude, and source distance. The calculation shows that dataset B gives a 6% higher peak ground acceleration compared to dataset A.