Abstract
The response of nuclear power plants (NPPs) to seismic events is affected by soil-structure interactions (SSI). In the present paper, a finite element (FE) model with transmitting boundaries is used to analyse the SSI effect on the response of NPP buildings subjected to vertically incident seismic excitation. Analysis parameters that affect the accuracy of the calculations, including the dimension of the domain and artificial boundary types, are investigated through a set of models. A numerical SSI analysis for the 10 MW High Temperature Gas Cooled Test Reactor (HTR-10) under seismic excitation was carried out using the developed model. The floor response spectra (FRS) produced by the SSI analysis are compared with a fixed-base model to investigate the SSI effect on the dynamic response of the reactor building. The results show that the FRS at foundation level are reduced and those at higher floor levels are altered significantly when taking SSI into account. The peak frequencies of the FRS are reduced due to the SSI, whereas the acceleration at high floor levels is increased at a certain frequency range. The seismic response of the primary system components, however, is reduced by the analysed SSI for the HTR-10 on the current soil site.
1. Introduction
The structural integrity of nuclear power plants (NPPs) built on medium or soft soil sites during an earthquake has been a focus of research in recent years [1–7]. Soil-Structure Interactions (SSI) may greatly amplify the seismic response of the NPP’s reactor building and increase the safety requirements [1, 2, 8]. Given their inherent safety features such as meltdown-proof safety, negative temperature reactivity coefficients, and passive decay heat removal [9, 10], High Temperature Gas Cooled Reactors (HTGRs) present advantages for construction on various soil sites when compared to Pressurized Water Reactors (PWRs). The seismic design of safety related nuclear structures and facilities on soil sites is built upon SSI analysis. Hence, SSI is of significant importance for nuclear safety and therefore the future development of commercial HTGR. The 10 MW High Temperature Gas Cooled Test Reactor (HTR-10) was constructed on a soil site in Beijing, China. Although SSI have been considered during the reactor design, a simplified model was adopted due to the lack of required computing capabilities. A more comprehensive SSI analysis will enhance the understanding of the dynamic response of the reactor building in seismic events.
A variety of SSI effects on the seismic response of NPPs have been investigated by recent studies. The floor response spectra (FRS) of the AP1000 nuclear island on five generic soil type sites as well as a hard rock site were compared by Tuñón-Sanjur et al. [1]. Chen and Maslenikov [3] calculated the SSI response of a nuclear reactor building and showed that FRS of the structure were highly sensitive to soil stiffness. Farahani et al. [4] carried out an SSI analysis of NPP buildings and concluded that the seismic responses were dependent on the subsurface profile. Saxena et al. [5] observed that slip and separation at the soil-foundation interface due to SSI have significant effect on the response of the reactor building. Politopoulos et al. [2] focused on a partially embedded NPP on a layered soil site and demonstrated that the base rocking excitation induced by SSI may amplify the seismic response.
Two groups of methods, namely, the substructure method and the direct method, are commonly used to study SSI. In the substructure method [11, 12], the SSI problem is divided into several subproblems which are solved separately based on the assumption of linear model and the solutions are superposed to give the complete results of the SSI problem. In the direct method, the combined soil-structure system is analysed in a single step. The interaction problem can be calculated using numerical methods such as the finite element method (FEM) [13], the boundary element method (BEM) [14], and their coupling procedures (BEM/FEM) [15, 16]. Comprehensive literature reviews can be found in papers and books [17–19]. Several computer codes designed for SSI analysis, such as SASSI [20] and CLASSI [21], have been developed. Compared to these codes, SSI analysis by general purpose commercial FE software provides advantages including modelling flexibility via the choice of element, improved nonlinear solvers, and postprocessing capabilities. When SSI is analysed using FEM, special handling, for instance, through transmitting boundaries, is required to simulate the wave propagation from the finite element mesh to the far field of infinite soil [22–25]. In this way, waves are artificially transmitted to the far field rather than reflected by the boundary. Among the various transmitting boundaries developed and implemented in SSI analysis, the most widely used one is the viscous boundary [22] which replaces the far field with viscous damping. Also widely used is the viscoelastic boundary [24] which employs springs to the viscous boundary to improve the accuracy. Since these artificial boundaries are approximate representation of infinite soil, modelling parameters such as size and shape of calculation domain and the artificial boundary type are expected to affect the reliability of the SSI analysis.
In this paper, a numerical model to simulate SSI was developed using commercial FE code. The model has been verified and the requirements for the modelling parameters are investigated. The seismic SSI of the HTR-10 NPP building is then analysed using the validated model. The FRS, which affect directly the design of nuclear components to withstand seismic events, are calculated to investigate the SSI effect.
2. Numerical Modelling
2.1. Dynamic Response Function in Frequency Domain
The SSI analysis was performed in frequency domain using a complex frequency response method. The calculation procedure is as follows.(1)The input acceleration of a seismic event, , in time domain is converted to frequency domain by Fourier transformation(2)The response function of the structure is calculated in frequency domain.(3)The response function in frequency domain is converted to time domain via inverse Fourier transformation in step (2) is derived from a three-dimensional FE model of the soil and the structure. The accuracy of the response function depends on the validity of the soil-structure model.
For a three-directional earthquake input, the acceleration input in formula (1) can be decomposed for each direction; that is, , so the frequency response function is written as
The component represents the response along direction when subjected to simple harmonic acceleration with unit amplitude in direction. The response corresponding to excitation in each of the three directions is decoupled and calculated separately.
2.2. FE Model for Calculation of Frequency Response Function
2.2.1. Modelling of Structure and Soil
A three-dimensional FE model representing the structure and soil with truncated boundaries as illustrated in Figure 1 is constructed to obtain the frequency response function in (3) for a soil-structure system subjected to vertically incident seismic excitation. As the nonvertically incident excitations have influence mainly on structures with large spans, such as dams, bridges, and pipes, the current method is considered adequate for the SSI analysis for NPPs. The subsurface soil is modelled as viscoelastic layered deposits with different shear modulus. The bottom of the model is set at the surface of the bedrock which has a shear wave velocity of 2400 m/s. The soil domain is meshed using structured hexahedral elements. The element size in direction is smaller than 1/5 of the local shear wave length for a cut-off frequency of 25 Hz according to the requirements of ASCE Standard 4-98 [26]. Input acceleration, that is, simple harmonic acceleration with unit amplitude, is applied to the nodes at the bottom surface of the model.
2.2.2. Artificial Boundary Conditions
In order to account for the infinite soil, viscous and viscoelastic boundaries have been implemented in the FE model in Figure 1.
For models with viscous boundaries, the wave propagation in the far field is modelled as viscous damping following Lysmer and Kuhlemeyer [22]. The equivalent damping constant is calculated aswhere is a viscous damping constant, is the soil density, and is the shear wave velocity. The subscriptions and correspond to shear and compression wave, respectively.
For models with viscoelastic boundaries, springs are added to the viscous model following Deeks and Randolph [24]. The spring constant for axisymmetric shear wave propagation problems is determined aswhere is the soil shear modulus and is the distance between the boundary and the symmetry axis. Liu et al. [27] reported that satisfactory results can be obtained using this damping-spring approach for both compression and shear waves’ propagation modelling. Therefore, this approach is adopted in this study, withThe coefficients and are set to be 2/3 and 4/3, and is approximately the distance between the truncated boundary and the structure, as recommended by Liu et al. [27].
The viscous damping and spring constant in (4) and (6) are discretized in dashpot and spring elements connected to the nodes of the truncated boundary following Liu et al. [27]. To simulate viscous boundary, 3 dashpots with orthogonal degrees of freedom are connected to each boundary node as illustrated in Figure 1. For viscoelastic boundary, 3 dashpots and 3 springs are connected to each boundary node. Fixed boundary conditions are applied to the dashpot and spring elements.
The earthquake motions at the truncated boundaries are applied as equivalent forces to the boundary nodes. For vertically incident seismic excitation, the equivalent force is calculated as [28]where is the damping constant of the dashpot elements and is the spring constant. are force, displacement, and velocity, respectively, corresponding to free-field motion. The damping constant (4) and spring constant (6) are used to calculate the force, , caused by the shear wave. Similarly, damping constant (4) and spring constant (6) are used to calculate the force, , caused by the compression wave.
In (7), are commonly derived based on one-dimensional theory of wave propagation. In this study, a method based on FE modelling is proposed to calculate these parameters directly. When free-field soil is subjected to one-dimensional input, the responses are the same at any position in the direction perpendicular to the excitation. Therefore, a long slice of soil represented by one layer of elements is used to model the free-field response as shown in Figure 2. Symmetry boundary conditions are applied to the side surfaces of the model by constraining the motion perpendicular to the excitation. Soil properties, mesh sizes in direction, and the input acceleration at the bottom of the model are the same as the soil-structure model in Figure 1.
In the first step of the analysis, free boundary conditions are applied to both ends of the model. The slice model is sufficiently long (>1000 m) in the direction of excitation so that the reflected waves from both ends are attenuated significantly to affect the middle of the soil slice model. The motion in this part of the model can be considered as free-field motion of infinite soil. The required , and can be calculated for the nodes in the middle cross section.
In the second step of the analysis, is applied to both ends of the slice to replace the free boundary in the first step. The boundary conditions on the sides and the input motion at the bottom of the model remain the same as these in the first step above. The reaction forces, , at the boundary nodes of the model are calculated and extracted.
From the two-step analysis, the free-field variables required to calculate the equivalent forces in (7) are obtained. This method facilitates the modelling of the artificial boundary.
3. Validation of the Numerical Model
3.1. Method and Model
The model validation is carried out by a simplified model of an NPP building with dimensions of 56 m × 50 m and 27 m height, as shown in Figure 3(a). The building is modelled using shell and beam elements with a size below 2 m. Distributed mass is added to the floor of the structure to represent the reactor and other major facilities in the building. The step profile soil type proposed by Tuñón-Sanjur et al. [1], which is based on a survey of 22 commercial NPPs in the United States, is used. The soil parameters for each layer are listed in Table 1. The damping ratios of all layers are assumed to be 5%, which is approximately an average value for the adopted soil profile. The thickness of the calculation domain in direction is 40 m, including layer 1 through 3 in Table 1. Vertical ( direction) size of each layer is also shown in Table 1. The horizontal ( and direction) element size is restricted to be less than 10 m. Sensitivity analysis was conducted and the results indicated that further refinement of the mesh in either horizontal or vertical direction would have negligible effects on the results.
(a)
(b)
The NPP model is subjected to one-directional excitation along axis. The FE model is shown in Figure 3(b). With this model, the components in (3) are calculated. Due to the symmetry of the NPP building, only half of the model is analysed.
In order to test the SSI model, 8 cases with varied sizes and artificial boundary types (see Table 2) were studied. Figure 4 is a diagram of the soil-structure calculation domain. To validate the models, an extended model with a size of 1600 × 800 m was also studied for comparison. Free boundary conditions were applied to the sides of the extended model. Most of the wave energy reflected from the truncated boundaries would be absorbed before reaching the building again. The results of the extended model can be regarded as more representative solution.
3.2. Results
3.2.1. Effect of Calculation Domain Size and Shape
Figures 5(a)–5(c) show the frequency response at the floor level (FL) +0.0 m in three orthogonal directions, direction (along the excitation), direction (perpendicular to the excitation), and in direction for cases 1, 2, 3, and 5 in Table 2. The accuracy of all 6 cases in Table 2 with viscous boundaries is summarised in Table 3.
(a) direction ()
(b) direction ()
(c) direction ()
As shown in Figure 5(a), significant differences are observed in the frequency response peaks corresponding to the natural frequencies of the soil-structure system. The smallest model in size (case 1, 96 × 70 m) is characterized by a peak value 52% less than the extended model. As the size of the model increases, the response curve changes to match closer that of the extended model. A model with a size of 280 × 160 m shows almost no difference in the peak values. Comparison between cases 1, 2, 3, 5, and 6 in Table 3 indicates that models with large size along the direction of excitation produce better results than small models. The maximum relative error is less than 10% over the whole frequency range when the size of the model along direction is 400 m; that is, the distance between the truncated boundary and the building is about 3 times the dimension of the building.
In and directions, as shown in Figures 5(b) and 5(c), the acceleration is smaller than that in direction. For models with dimensions of 96 × 70 m and 160 × 160 m, the responses in and directions differ from those of the extended model. A good agreement is observed between the response curve of the extended model and models with dimension greater than 280 m. Since motions perpendicular to the excitation are much smaller than that along excitation, their contribution to the overall seismic response is negligible.
Based on the results in Figure 5, when the distance between the truncated boundary and the structure is ~3 times the length of the structure, the response function is calculated with acceptable accuracy. This model size is adopted in the HTR-10 study described in Section 4.
Figures 6(a)–6(c) illustrate the frequency responses at FL +0.0 m in three orthogonal directions for cases 3 and 4 when excited in direction. Comparison between case 4 (280 × 260 m) and case 3 (280 × 160 m) shows that excessively enlarging the calculation domain in the direction perpendicular to excitation has limited effect on the accuracy of the response curve.
(a) direction ()
(b) direction ()
(c) direction ()
The dynamic response of the soil-structure model subjected to a simple harmonic acceleration with unit amplitude at a typical frequency of 4.7 Hz, corresponding to the fundamental frequency in direction, is shown in Figure 7 in order to illustrate the appropriate shape of soil model in the - plane. The dimension of the soil influenced by the structure along the excitation to that perpendicular to the excitation is about 2 : 1. In this sense, for a square-based structure, a soil model with an aspect ratio (i.e., model size along the excitation to that perpendicular to the excitation) of 2 is reasonable to represent the response in the soil. This aspect ratio is adopted to determine the model size perpendicular to the excitation in the following study.
3.2.2. Effect of the Artificial Boundary Type
To investigate the effect of viscoelastic boundary on the soil-structure model, the results of cases 7 and 8 are compared to those of cases 1 and 2. The amplification factors for viscoelastic and viscous boundaries of the 90 × 70 m models are presented in Figure 8. As the model dimension is not long enough in the direction of excitation, the result from case 1 is not acceptable (see Section 3.2.1). The same conclusion is valid for the model with viscoelastic boundaries. As the model dimensions increase, the difference between the two boundary types was observed to decrease (2.5% for 160 × 160 m model). As illustrated by Deeks and Randolph [24], the viscoelastic boundary improves the results related to rigid body displacement and response in the low frequency range. In the current study, because the natural frequencies of the NPP structures are usually relatively high (between 4 Hz and 6 Hz), the advantages of the viscoelastic boundary are not significant.
(a) direction ()
(b) direction ()
(c) direction ()
In summary, the artificial boundary type has little effect on the accuracy of SSI analysis for NPP buildings when using complex frequency response method. Both of viscous and viscoelastic artificial boundaries provide satisfactory results when compared with an extended model following the recommended dimension in Section 3.2.1. Since the viscous boundary provides implementation benefits, it is adopted in the HTR-10 study.
4. SSI Study for HTR-10
The numerical model as described in Sections 2 and 3 was used for seismic SSI analysis of HTR-10 NPP.
4.1. SSI Model of HTR-10
The HTR-10 NPP building consists of a reactor hall and an annex building. The reactor hall (30 m long, 24 m wide) rests on soil at elevation −17.0 m and rises to an elevation of +28.0 m. The 5 m high annex building (46 m long, 34 m wide) surrounds the reactor hall with a foundation at elevation of +0.0 m. The primary system cabin inside the reactor hall is a shaft structure stretching from the −17.0 m level to the +11.0 m level. The NPP is constructed on a soil site in Beijing, China. The parameters of the soil layers based on geological surveys are listed in Table 4. The soil layers above the bedrock are 52.1 m in thickness.
The HTR-10 building is subjected to a two-directional seismic input in the - plane. The standard spectrum recommended by NRC Regulatory Guides 1.60 [29] was normalised to 0.32 g zero period acceleration (ZPA) and adopted as design spectrum. Synthetic earthquake time histories developed from the design spectrum in strict conformance to the requirements in ASCE 4-98 [26] and Standard Review Plan [30] as shown in Figure 11 were applied to the ground surface as the design free-field motions. This earthquake motion was then deconvoluted to the bottom layer of soil (−52.1 m) using computer program SHAKE [31] to calculate the acceleration input used in (1).
The HTR-10 building is modelled using shell and beam elements as shown in Figure 9. The Reactor Pressure Vessel (RPV), Steam Generator (SG), and Hot Gas Conduct Vessel (HGCV) are simplified to a lumped-mass model which captures the essential features of the structures. Other facilities and components are modelled by distributed mass on the floor.
The responses corresponding to excitation in each direction are calculated separately according to Section 2.1. The soil-structure models used to calculate the frequency response function when subjected to excitation in and directions are shown in Figures 10(a) and 10(b), respectively. FE models with dimensions of 330 × 170 m and viscous boundaries are adopted for the SSI analysis (see Section 3). The equivalent forces at the truncated boundaries were obtained by the method in Section 2.2 using the soil properties in Table 4. The results in the two directions are combined using (1) through (3) to produce the response of the structure. The FRS with damping ratio of 5% were calculated at all floor levels and the locations of the supports of the RPV and SG.
(a) Excitation in direction
(b) Excitation in direction
(a) direction
(b) direction
In order to investigate the effect of SSI, a fixed-base model assuming that the NPP building is sitting on the rigid rock was also analysed and the FRS were generated for comparison. In the fixed-base model, the design free-field motions are directly applied to the nodes that connect to the rock, including the floor at −17.0 m, the embedded walls of the reactor hall, and the floor of the annex building at +0.0 m.
4.2. SSI Effect of HTR-10
In order to investigate the effects of SSI on the reactor building and the main components, the FRS were investigated at three different floor levels, including the bottom of the building (FL −17.0 m), top floor (FL +11.0 m), and RPV and SG supports (FL −0.29 m). The FRS at these floor levels are shown in Figure 12. The results of the SSI model (solid lines) are compared with those of the fixed-base model (dashed lines). At all floor levels, the spectrum curves in both and directions are clearly different for the SSI and fixed-base results.
(a) FL −17 m direction
(b) FL −17 m direction
(c) FL +11 m direction
(d) FL +11 m direction
(e) RPV and SG supports at FL −0.29 m direction
(f) RPV and SG supports at FL −0.29 m direction
Figures 12(a) and 12(b) are the FRS at FL −17.0 m in and directions, respectively. The dashed lines represent the fixed-base results at FL −17.0 m, which are the same as the design spectra (i.e., free-field response spectra at FL +0.0 m). The free-field response spectra at FL −17.0 m obtained from the soil slice model are also presented in Figures 12(a) and 12(b) as long dashed lines. The results show that the free-field response is amplified from −17 m (long dashed lines) to 0 m (dashed lines) by the soil layers. The embedded structure resting on the −17.0 m soil layer is subjected to lower level of excitation than that at 0 m. It can be concluded that using fixed-base model is conservative at lower floors for HTR-10.
The FRS for FL +11.0 m is illustrated in Figures 12(c) and 12(d), demonstrating that the FRS is substantially altered because of the SSI. The peak frequency is reduced from over 10 Hz to less than 5 Hz due to the SSI. The relatively high frequencies in the range from 6 Hz to 30 Hz are suppressed by the SSI. The responses in the relatively low frequencies ranging from 3 Hz to 5 Hz are amplified by a maximum of 2.5 times. Thus, the SSI effect should be considered in calculating FRS at floors above the ground surface, especially at high floor levels.
Figures 12(e) and 12(f) present the FRS at the supports of RPV and SG in the primary system cabin at FL −0.29 m. The acceleration in the high frequency range between 6 Hz and 30 Hz is reduced by the SSI for both and directions. At lower frequency ranges, the SSI model produces larger responses in both and directions. However, the increase is not as obvious as that in FL +11.0 m. Since the fundamental frequencies of the integrated structure of primary system components (i.e., RPV, SG, and HGCV) in both and directions are around 10 Hz, seismic responses of these components are reduced when taking SSI into consideration for the current soil sites.
5. Conclusions
In this study, a numerical method for seismic SSI analysis with commercial code is developed. Three-dimensional FE models composed of upper structure and subsurface soil with transmitting boundaries were solved in the frequency domain to determine the frequency response function of an NPP under vertically incident seismic excitation. A comprehensive investigation of this numerical method is carried out and the modelling parameters which produce reliable results were determined. Several recommendations are provided for the SSI analysis of NPP subjected to horizontal earthquake ground motion.(1)The distance between the truncated boundary of the soil model and the building should be at least 3 times the length of the structure in the direction of excitation.(2)Soil model with an aspect ratio (size along the excitation to that perpendicular to the excitation) of 2 : 1 produces satisfactory results for a square-based structure.(3)SSI analysis in the frequency domain with viscous boundary is adequate to produce satisfactory results for NPP.
The seismic response of the HTR-10 reactor building subjected to horizontal ground motion was analysed. The results show that the SSI do not amplify significantly the dynamic response of HTR-10 for the particular soil profile. The following conclusions are made by comparing the FRS results from SSI and fixed-base models.(1)The response at lower floor levels is reduced due to SSI.(2)The FRS curve is substantially changed at higher floor levels. The peak frequency of the FRS is reduced due to the soil flexibility and the spectral acceleration increases at certain frequency range. SSI analysis is essential for the seismic designs of the safety related components installed at these floors.
The above results are valid for a broad range of earthquakes with ground motion that can be represented by the standard spectrum of RG 1.60 [29]. The analyses performed in this study may be used for seismic design of NPPs based on codes and regulations. However, SSI effect under a wider range of seismic excitations, such as near-source ground motions characterized by large ground velocity and displacement pulses with long-periods [32] should be investigated separately in order to obtain a more comprehensive conclusion about the SSI effect on NPPs.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
This work is supported by the National Natural Science Foundation of China under Grant no. 11302116. The authors wish to express their thanks for the support.