Methodology for Pressurized Thermal Shock Analysis in Nuclear Power Plant

1. Introduction

The relevance of the fracture mechanics in the technology of the nuclear power plant is mainly connected to the risk of a catastrophic brittle rupture of the reactor pressure vessel. There are no feasible countermeasures that can mitigate the effects of such an event that impair the capability to maintain the core covered even in the case of properly functioning of the emergency systems.

The origin of the problem is related to the aggressive environment in which the vessel operates for long term (e.g. more than 40 years), characterized by high neutron flux during normal operation. Over time, the vessel steel becomes progressively more brittle in the region adjacent to the core. If a vessel had a preexisting flaw of critical size and certain severe system transients occurred, this flaw could propagate rapidly through the vessel, resulting in a through-wall crack. The severe transients that can lead the nuclear power plant in such conditions, known as Pressurized Thermal Shock (PTS), are characterized by rapid cooling (i.e., thermal shock) of the a part of the internal reactor pressure vessel surface that may be combined with repressurization can create locally a sudden increase of the stresses inside the vessel wall and lead to the suddenly growth of the flaw inside the vessel thickness.

Based on the long operational experience from nuclear power plants equipped with reactor pressure vessel all over the world, it is possible to conclude that the simultaneous occurrence of critical-size flaws, embrittled vessel, and a severe PTS transient is a very low probability event. Moreover, additional studies performed at utilities and regulatory authorities levels have shown that the RPV can operate well beyond the original design life (40 years) because of the large safety margin adopted in the design phase.

A better understanding and knowledge of the materials behavior, improvement in simulating in a more realistic way the plant systems and operational characteristics and a better evaluation of the loads on the RPV wall during the PTS scenarios, have shown that the analysis performed during the 80’s were overly conservative, based on the tools and knowledge available at that time.

Nowadays the use of best estimate approach in the analyses, combined with tools for the uncertainty evaluation is taking more consideration to reduce the safety margins, even from the regulatory point of view. The US NRC has started the process to revise the technical base of the PTS analysis for a more risk-informed oriented approach. This change has the aim to remove the un-quantified conservatisms in all the steps of the PTS analysis, from the selection of the transients, the adopted codes and the criteria for conducting the analysis itself thus allow a more realistic prediction.

This change will not affect the safety, because beside the operational experience, several analysis performed by thermal hydraulic, fracture mechanics and Probabilistic Safety Assessment (PSA) point of view, have shown that the reactor fleet has little probability of exceeding the limits on the frequency of reactor vessel failure established from NRC guidelines on core damage frequency and large early release frequency through the period of license extension. These calculations demonstrate that, even through the period of license extension, the likelihood of vessel failure attributable to PTS is extremely low (≈10-8/year) for all domestic pressurized water reactors.

Different analytical approaches have been developed for the evaluation of the safety margin for the brittle crack propagation in the rector pressure vessel under PTS conditions. Due to the different disciplines involved in the analysis: thermal-hydraulics, structural mechanics and fracture mechanics, different specialized computer codes are adopted for solving single part of the problem.

The aims of this chapter is to present all the steps of a typical PTS analysis base on the methodology developed at University of Pisa with discussion and example calculation results for each tool adopted and their use, based on a more realistic best estimate approach.

This methodology starts with the analysis of the selected scenario by mean a System Thermal-Hydraulic (SYS-TH) code such as RELAP5[2][3], RELAP5-3D[1], CATHARE2 [4][6], etc. for the analysis of the global behavior of the plant and for the evaluation of the primary side pressure and fluid temperature at the down-comer inlet.

For a more deep investigation of the cooling load on the rector pressure vessel internal surface at small scale, a Computational Fluid Dynamics (CFD) code is used. The calculated temperature profile in the down-comer region is transferred to a Finite Element (FE) structural mechanics code for the evaluation of the stresses inside the RPV wall. The stresses induced by the pressure in the primary side are also evaluated.

The stress intensity factor at crack tip is evaluated by mean the weight function method based on a simple integration of the stresses along the crack border multiplied by the weight function. The values obtained are compared with the critical stress intensity factor typical of the reactor pressure vessel base material for the evaluation of the safety margin.

2. Origin of the problem

The internal components of PWR vessels that are closest to the core (baffles, formers and core envelope in Solution Annealed 304 stainless steel[7], bolts in Cold Worked 316 stainless steel[8], etc.)[9] are highly irradiated; the most irradiated areas of some of these components may be exposed to doses reaching around 80 dpa after 40 years of operation. This neutron irradiation changes their microstructure and their mechanical properties, so they harden, lose ductility and toughness, suffer irradiation creep [9][10][11][12]. In addition, these changes seem to be the basis of increased sensitivity to stress corrosion [13].

Hardening (or embrittlement) starts at the nanometer level as the high energy neutrons are absorbed by the material causing lattice defects which cluster. The mechanisms proposed for the radiation damage are many, but on a fundamental level a single neutron scattering event can be considered.If a neutron of sufficient energy scatters off a nucleus, the nucleus itself is displaced. The atom associated with the nucleus finds itself embedded into the structure elsewhere in a high-energy, interstitial site. It is termed a self-interstitial as the matrix and interstitial atoms are in principle the same. The site the atom previously occupied is now empty: it is a vacancy. In this way, self-interstitial-vacancy pairs are formed.

Neutron scattering events are not isolated. On average, each displaced atom might then go on to displace further atoms, and likewise the neutron that caused the first displacement might go on to displace further atoms. This means that there is a local cascade of displacements, known as a displacement spike, within which there is a large amount of disorder in the structure.

Both the interstitial atoms and vacancies can diffuse through the lattice, but the interstitial atoms are more mobile. Both interstitials and vacancies are eventually removed from the lattice (when they reach sinks such as dislocations or grain boundaries). However, they are also always being generated by the neutron radiation. Thus steady-state populations of interstitials and vacancies are formed.

The majority of the self-interstitial quickly cluster to form small, disc-shaped features that are identical to small dislocation loops. Along with self-interstitial, these loops are very mobile. Diffusion of self-interstitial and loops within the cascade region causes additional recombination prior to their rapid long-range migration (unless they are strongly trapped by other defects or solutes). Although they are less mobile than the self-interstitial, vacancies also eventually diffuse.

In summary, displacement cascades produce a range of sub-nm clusters (defects, solutes, and defect-solute complexes) that directly contribute to irradiation hardening. Expressing damage exposure, or neutron dose, in terms of displacements-per-atom (dpa) partially accounts for the effect of the neutron energy spectrum on the generation of cascade defects and the net residual defect production scales with dpa.

The important thing to know, however, is that the neutron fluence causes the material to loose fracture toughness and in addition causes a shift upwards in the nil ductility transition temperature (RTNDT). It is this shift in the RTNDT value that is the heart of the problem. Consider Figure 1 which illustrates Charpy-V-Notch (CVN) results for ferritic steel before and after irradiation. The CVN test is an old approach to measure the fracture toughness of a material that is still in use today. Improved techniques are available but the CVN results can be used to clearly show the impact of neutron fluence.

The upper curve on Figure 1 represents the CVN energy curve for the un-irradiated specimens. The test measures the energy required to break a specific specimen at a given temperature. If the material is ductile (tough), the energy required to break the specimen is high. As the material is cooled, it loses fracture toughness (becomes more brittle) as illustrated by the curves in Figure 1. The upper curve shows relatively high fracture toughness for temperatures greater than 20° C. Since RPV temperatures are not expected to drop to this level, this material will remain tough during an overcooling event [14].

The risk is mainly associated to the presence of defect in the welding lines that face the core region, because at the edge of this defect an intensification of the stresses can occurs in case of a fast overcooling phase can be generated during the plant life by injection of cold water by the emergency systems.

From a safety point of view, the material properties of the RPV have to be regularly and carefully checked in order to evaluate the embrittlement level of the material by mean the analysis of specimens placed close to the inner RPV surface that can be analyzed by mean the classical CVN test.

Figure 1.

Effect of radiation damage on the CVN transition characterization of ferritic steels

Beside this process a complex analysis is required aimed at identifying those scenarios that can lead an overcooling on the RPV internal surface where a crack flaw is supposed to exist.

The overcooling transients are usually very complex. It is often not possible to define in advance conservative or limiting conditions for all system parameters. Engineering judgment might not be sufficient to decide whether an accident under consideration is, by itself, a PTS event or along with other consequences can lead to a PTS event that may potentially threaten RPV integrity. Therefore thermal hydraulic analyses are often necessary for choosing, from a number of accidents, those initiating events and scenarios that can be identified as limiting cases within the considered group of events.The calculation period of a transient should be long enough to reach stabilized conditions or at least to overreach the critical time from the point of view of RPV integrity [16].

A wide variety of transients can contribute to the risk of vessel failure. These transients include reactor system overcooling attributable to a LOCA or a stuck-open primary side relief valve, a component failure that results in an uncontrolled release of steam from the secondary side (e.g., Main Steam Line Break (MSLB) or stuck-open secondary side relief valve), or a control system failure that results in overfilling the steam generators. Combinations of failures are also of concern and have to be taken into account during the analysis.

For a wide discussion of the PTS issue in PWR Nuclear Power Plant (NPP) see[15], hereafter only a short overview is reported.

Beside the main categories of transients highlighted, some phenomena were deemed to be most important to down-comer conditions during PTS events:

• natural circulation;

• Emergency Core Cooling Systems (ECCS) injection (mixing and condensation following ECCS injection);

• flow stagnation in case of primary system pressurization.

Natural circulation and flow stagnation are important because if loop mass flow continues (or restarts during a transient), warm water at the average coolant system temperature will be flushed through the reactor vessel down-comer, increasing the down-comer fluid temperature. In contrast, if the loop flow is stagnant, the cold ECCS water will not be mixed with water flowing from other parts of the reactor system and the down-comer temperature will be colder in comparison with the natural circulation case. Integral system response is important because the ECCS injection behavior (flow rates, timing, and to some extent temperatures) are functions of the overall system behavior. System pressurization is itself the primary phenomenon in the PTS analysis. The phenomena listed above were considered because of their potential impact on the down-comer conditions, in particular, the mixing phenomena occurring in the down-comer have the capability to mitigate the cooling effect of the cold water injected by the emergency systems [16].

3. UNIPI methodology for PTS analysis

The objective of the PTS analysis is to determine the safety margin for the RPV operability. The safety margin is obtained comparing the stress intensity factors at crack edge calculated for the identified spectrum of overcooling scenarios with the critical stress intensity factor of the RPV material obtained from the CVN test. In each one of the selected scenarios, there must be enough margin to be sure that the vessel can withstand the selected loads conditions.

University of Pisa developed a methodology concerning the use of a chain of codes to quantify this margin by mean a deterministic approach to the PTS issue,see[17][18].

A preliminary analysis of the plant configuration and logics is required in order to identify the spectrum of the scenarios that lead to an overcooling of the down-comer region.

The methodology starts with the thermal hydraulic analysis of the Nuclear Power Plant (NPP) using a SYS-TH code such as RELAP5, CATHARE2, or equivalent, during a selected transient scenario. The goal of this step is to evaluate the plant response and to calculate the cooling load induced on the internal RPV wall surface by the Emergency Core Coolant (ECC) injection or by the cooling plug following a MSLB initiating event, to calculate the primary circuit pressure and to provide boundary conditions for the next step.

If the transient evolves in single phase, a more detailed analysis of the mixing phenomena occurring in the down-comer region can be performed by mean a CFD code. The result of this step needed for the PTS analysis is the temperature distribution inside the down-comer.

Figure 2.

UNIPI Methodology for PTS Analysis

The thermal load to be applied to the FE model for the stress analysis can be extrapolate from the CFD result considering the temperature profile at the interface between down-comer fluid model and RPV wall. This step is accomplished using suitable subroutines developed for this purpose. The stresses due to mechanical load such self-weight, pretension in bolts and internal pressure are also accounted in the Ansys FE calculation. In the last step of the analysis, the Stress Intensity Factor (SIF) KI is calculated by means the Weight Function method, once the stresses generated by the loads identified before are known. The KI has to be compared with the critical SIF (KIc) of the material for the evaluation of the safety margin for the RPV operability. In the following paragraphs an example application to the methodology is provided.

In the next paragraphs a more detailed analysis will be provided for each steps of the analysis previously identified.

4. Analysis at system level

The thermal-hydraulic analysis has the task to determine accident sequences where the temperature differences in adjacent parts of the RPV-wall are large and last for a longer period of time. A catalogue on all relevant load conditions having the potential to cause such temperature differences has to be drawn up for each reactor plant individually, since each plant is characterized by system-specific equipment features.

Under the conditions described above, a high degree of detailing is required from a thermal-hydraulic analysis. Usually, the down-comer is being subdivided in thermal-hydraulic analyses into one or several vertically arranged flow areas, the so-called parallel channels. With regard to the accuracy required for a detailed brittle fracture analysis, this subdivision is not sufficient. This becomes obvious when considering the cooling mechanisms which have been observed in the corresponding test facilities.

Two different cooling mechanisms, the plume and stripe cooling, are considered in connection with the determination of thermal loads. Stripe cooling occurs when the cold-leg emergency core cooling takes place at a time when the water level in the down-comer has fallen below the opening of the cold-leg nozzle. Here, those loss of coolant accidents are taken into consideration which either show a sufficiently large leak cross section or where the emergency core cooling has limited availability.

The stripe cooling causes the biggest thermal load by far, because an only moderately heated cold water stripe of relatively small width cools down the RPV-wall with a large temperature difference to its surrounding area.

However, the significance of cold water stripes for the determination of thermal loads is restricted, since corresponding experimental analyses show that cold water stripes already become detached from the RPV-wall with a relatively low mass flow rate.

Depending on the constructive layout of the cold-leg nozzle, 10 kg/s are for example sufficient to detach a water stripe from the vessel surface (see Figure 3). The injection rate for a high-pressure emergency core cooling system in the conventional PWR NPP is higher than this value; therefore, the stripe cooling is most significant for the area of the cold-leg nozzle.

Figure 3.

Water stripes in the down-comer with large and small mass flow rate mechanisms

The width of the cold water stripe depends on the flow velocity of the draining cold water stripe. Experimental investigations show that this velocity is determined by the flow phenomenon "critical flow”. Dependent on the mass flow of the emergency cooling water, a water level in the cold-leg appears at which the flow velocity exactly corresponds to the critical flow velocity. It turns out that for the relevant mass flows stripe widths of about 10 to 30 cm are to be expected in the cold-leg nozzle. Within this width, the nozzle is cooled down locally.

Cold water plumes are formed if the cold-leg emergency coolant is injected into a down-comer filled with water. Such a situation arises, e.g., from loss-of-coolant accidents with a smaller leak cross section in the hot leg.

Since these plumes stay for a longer period until final mixing, i.e. about one to two hours, the temperature differences can act on the RPV-wall correspondingly long. In contrast to the stripe cooling, the concentration of cold emergency coolant here is lower. By admixture with the water surrounding the plume, there is a permanent exchange. The plume width is dependent on the injected mass flow and the exchange with its environment.

Such phenomena require a more detailed analysis at small scale level in order to better identify the position, the size of the overcooled region, the duration of the cooling phase and the cooling rate.

In order to reach this goal the use of the CFD code is envisaged taking the boundary conditions from the system thermal hydraulic code results. The only limitation in the use of such codes is the large computational resources required.Because of this, the analysis is restricted to the time period during the maximum cooling phase.

The boundary conditions to the CFD code from the system thermal hydraulic code are transferred following a scheme reported in Figure 4. The mass flow rates and the fluid temperature of the injected fluid in the cold legs and down-comer calculated by the SYS-TH code are imposed to the CFD model.

Figure 4.

Coupling technique between SYS-TH and CFD codes

5. Analysis at small scale by means a CFD code

The second step of the methodology foresees the analysis at small scale by mean a CFD code for a more detailed calculation of the profile temperature inside the down-comer.

The CFD computer codes solve the Navier-Stokes equations optionally two or three-dimensionally, having the potential to reach the necessary degree of detailing.

With the increasing speed of modern computers, CFD techniques are becoming more widely used and may provide the best tool for computing the thermal fluid mixing. An example application is the emergency injection in the cold-leg by the emergency water reported in Figure 5.

Figure 5.

Mixing phenomena occurring during emergency injection

These approaches, however, still suffer from considerable user effect and the need for best practice guidelines continues. Challenges for the single phase CFD user include the turbulence modeling approach which must be able to adequately simulate the various mixing regimes which each have their own unique geometry and driving forces. In addition, the wall modeling approach is important. Mixed convection in the down-comer region is expected and typical CFD wall treatments do not account for this phenomenon.

All of the issues associated with CFD for the single phase PTS issue are compounded in the multi-phase problem by the relative immaturity of the multi-phase CFD techniques. For the near future, system analysis codes will still provide the overall system behavior and experimental results will provide the best source of information on the details of the multi-phase behavior related to PTS [14].

For reliable simulation of PTS related mixing processes the CFD methods must be validated to determine how well the CFD model, defined by the detail level of model geometry, the mesh and the used numerical and physical models can simulate the relevant physical processes and produce the needed data. The final target data of thermal hydraulic analysis are the pressure and temperature fields on structures needed as an input for structural analysis. The CFD model should be able to model the complex mixing and stratification processes in the cold legs and the down-comer of the pressure vessel as well as the heat transfer between fluid and structures accurately enough to reproducethis data. In the Figure 6 an example result of the CFD calculation of the temperature profile inside the down-comer is provided.

Figure 6.

Example result of the down-comer CFD result

The temperature profile and pressure profiles calculated with the CFD code are the boundary condition for setting up the FE simulation for the evaluation of the stresses inside the RPV wall.

6. Structural mechanics analysis

A common approach in the calculation of the stress at the crack edge is the simulation of a part of the vessel obtained by mean symmetry consideration on the geometry and load condition where the crack is supposed to be. This approach in followed basically for saving computational resources, an example of this approach can be found in the ref. [19].

In the methodology developed at University of Pisa, the approach followed is to model the full geometry of the RPV without modeling the crack inside the wall. The reason of this approach is due to a more precise analysis of the local stress inside the RPV wall avoiding any simplification due to the fact that the phenomena occurring in a PTS scenario are intrinsically not symmetric. This choice is supported even from the adoption for the calculation of the stress intensity factor of the weight function method (described in the next paragraph) that needs the stress calculated in the undamaged structure (see Figure 7).

Figure 7.

FE model of the RPV (example)

The thermal load calculated in the previous step by mean the CFD code has to be applied FE model. At this point, a problem arises because the CFD model is based on a finite volume mesh, while the structural mechanics model is based on FE mesh.Two techniques (shortly described hereafter) can be adopted, to transfer the temperature time history from one model to another:

• The first foresees the evaluation of the temperaturetime history of the FE RPV nodes internal surface interpolating the temperature of the CFD down-comer model closest nodes (external down-comer surface, Figure 8) and to calculate the profile inside the wall thicknesses solving the heat conduction equation by mean the FE code itself. The FE model implements all the material properties of the RPV for solving the thermal and stress calculation inside the thickness.

• A second approach foresees the modeling of the RPV wall and down-comer with the CFD model in order to solve the conductivity problem directly. This is called in literature as “conjugate heat transfer calculation”.In this approach the temperature time trend in each of the nodes of the FE model is obtained interpolating the temperature of the closest nodes in the RPV CFD model. This second technique is more time consuming and requires a more sophisticated subroutines for performing the transfer (from RPV-CFD model to RPV-FE model) compared to the previous one because the number of nodes to manage is in the order of magnitude of millions.

An example of this transfer technique is reported in Figure 9.

Figure 8.

Section of the RPV

Figure 9.

Temperature profile in the FE model of RPV

Once the pressure and the temperature profiles are implemented in the FE model, the stress analysis can be executed. An example result is shown in Figure 10.

Figure 10.

Von Mises stresses

The knowledge of the stress profile inside the RPV wall is the basis for the application of the Weight Function (WF) method that is described in the next paragraph for the calculation of the stress intensity factor.

7. Fracture mechanics

Most numerical methods require a separate calculation of the stress intensity factor for each given stress distribution and each crack length. The weight function procedure developed by Bückner[20] simplifies the determination of stress intensity factors. If the weight function is known for a crack in a component, the stress intensity factor can be obtained by multiplying this function by the stress distribution and integrating it along the crack length. The weight function does not depend on the special stress distribution, but only on the geometry of the component.

If σ(x) is the normal stress distribution in the uncracked component along the prospective crack line of an edge crack (see Figure 11), the stress intensity factor is given by the expression (1):

Figure 11.

One dimensional crack model scheme

$K_I=\underset{0}{\overset{a}{{\displaystyle \int }}}\sigma \left(x\right)h_I\left(x,a\right)dx$E1

The integration has to be performed over the crack length. The WF h(x,a) does not depend on the special stress distribution, but only on the geometry of the component.

The general procedures for the determination of weight functions are described below for the weight function component h_I. The relation of Rice (see[20]) allows to determine the weight function from the crack opening displacement V_r(x,a) under any arbitrarily chosen loading and the corresponding stress intensity factor K_Ir(a) according to:

$h_I\left(x,a\right)=\frac{E^{\text{'}}}{K_{Ir}\left(a\right)}\frac{\partial V_r(x,a)}{\partial a}$E2

(E' = E for plane stress and E' = E/(1-ν²) for plane strain conditions), where the subscript r stands for the reference loading case. It is convenient to use σ_r(x) = σ₀ = constant for the reference stress distribution.

One possibility to derive the weight function with eq.(2) is the evaluation of numerically determined crack opening profiles which may be obtained by Boundary Collocation Method (BCM) computations. For more details on this method see [21].

Once the weight function is defined it can be implemented in a simple spreadsheet program that can easily compute the integral for the stress intensity calculation taking into account the stresses evaluated in the undamaged structure by mean the FE code described in the previous paragraph.

Figure 12.

Safety margin evaluation for RPV brittle rupture

The results of the fracture mechanics analysis is a curve representing the value of the stress intensity factor in the selected scenario that is compared with the critical stress intensity factor curve of the material obtained from the analysis of the specimens by mean the CVN test.

In Figure 12 an example of this comparison is shown. The distance between the two curves gives the safety margin for the RPV operability.