Analysis of LBLOCA of APR1400 with 3D RPV model using TRACE

is very dif ﬁ cult to capture the multi-dimensional phenomena such as asymmetric ﬂ ow and temperature distributions with the one-dimensional (1D) model, obviously, due to its inherent limitation. In order to overcome such a limitation of the 1D representation


Introduction
Thermal-hydraulic (TH) system codes have aimed primarily at evaluating the TH behaviors and the performance of the engineered safety features of nuclear power plant systems during postulated accidents for both design and licensing purposes.Thus, the conservatism was of importance in the safety analysis in order to guarantee the safety of nuclear power plants [1].However, the perspectives have been changed with the advancement of the probabilistic safety assessment (PSA) method as it revealed that blindly posed conservatism could mislead the reactor safety [2].The objective of the system code analysis has been extended to produce best-estimate results for the reactor transients and accidents as realistic as possible.At present, the best-estimate analysis is conducted for selected design basis accidents, such as large-break loss-of-coolant-accidents (LBLOCAs), with the quantification of uncertainties in the analysis [3].
The best-estimate analysis has been conducted by using bestestimate TH system codes which typically describes the nuclear reactor system with one-dimensional (1D) representation.However, multi-dimensional phenomena such as asymmetric flow and temperature distributions could occur during the transients and accidents in complex components such as the reactor pressure vessel (RPV).It is very difficult to capture such multi-dimensional phenomena with the 1D model, obviously, due to its inherent limitation.In order to overcome such a limitation of the 1D representation, three-dimensional (3D) components for the multidimensional analysis capability have been developed in many state-of-the-art system codes such as MARS-KS [4], CATHARE [5], TRACE [6] and so on.TRACE is one of those best-estimate system analysis codes which has developed by United States Nuclear Regulatory Commission (USNRC) by consolidating the legacy codes of USNRC and implementing the state-of-the-art physical and numerical models [7].The code is capable of multi-dimensional analysis with its own 3D component, namely VESSEL.
Since the 3D component can simulate the multi-dimensional phenomena in a complex geometry, the typical application of the 3D component is the RPV.Especially, the reactor core in the RPV has been modeled with the 3D component in many accident analyses by sectioning and grouping the fuel assemblies based on the power distribution [8e11].Concerning the multi-dimensional effect, Clifford and Ferroukhi [12] conducted a series of calculations with 1D and 3D models for the reactor core using TRACE to figure out the effect of the 3D nodalization on reflood and CCFL (Counter-Current Flow Limitation).They found that, for large PWR (Pressurized Water Reactor) core, CCFL at the core upper plate could be simulated more realistically with the 3D model.In addition, their analysis of LBLOCA at a cold leg showed that the uniform quench front development was made in the 3D model during the reflood phase.They found that the uniform distribution of the quench front in the 3D model was attributed to the inter-channel crossflows which prevented the unbalanced level development featured in the 1D model.
The utilization of the 3D core model is also encouraged especially when the detailed TH analysis is required for multi-physics problems such as SLB (Steam Line Break) and RIA (Reactivity-Initiated Accident).The power behaviors of both accidents should be modeled by means of more detailed reactor kinetics model in order to obtain realistic analysis results.TRACE provides the spatial reactor kinetics simulation capability by coupling with a 3D kinetics code, PARCS (Purdue Advanced Reactor Core Simulator) [13] for aforementioned purpose.In the analysis of those accidents, detailed TH calculation is necessary to provide more realistic boundary conditions for the power calculation.Miglierini and Kozlowski [14] showed a good example of the RIA analysis with TRACE/PARCS for the VVER-1000 type PWR by postulating the control rod ejection at a specified assembly.They used the 3D model to consider the spatial configuration of the reactor core.The detailed core TH analysis is also of importance for the steam line break (SLB) accident, which accompanies asymmetric cooling and corresponding reactivity changes within the core.An international benchmark organized by USNRC and OECD/NEA (Organisation for Economic Co-operation and Development Nuclear Energy Agency) [15] indicated that the use of coupled 3D core model led to different results from the conventional 1D model.This benchmark revealed that the conservatism of the point kinetics model brought recriticality while the 3D model did not.
With those in mind, a standard 3D analysis model of APR1400 (Advanced Power Reactor 1400) has been developed in this study by using TRACE.The VESSEL component has been employed for modeling the entire RPV.Initially, the reactor core has been modeled by implementing all the assemblies in the active core considering the coupled analysis with PARCS.The zone-averaged core model with multiple hot pins has been also developed since such a detailed core representation is not necessary for LBLOCA analyses without the PARCS coupling.The consistency of both models has been confirmed by a comparative analysis for a doubleended guillotine break at a cold leg.

TRACE model for APR1400 for LBLOCA analysis
The standard 3D analysis model of TRACE has been developed on the basis of a thermal-hydraulic input for MARS-KS [16].As depicted in Fig. 1, the primary side consists of a RPV, a pressurizer (PZR), four reactor coolant pumps (RCPs), two steam generators (SG), and piping systems.An additional volume using a BREAK component is connected to the top of the PZR to give the pressure boundary condition during the steady-state calculation.The connection between the PZR and the pressure boundary is closed in the middle of the steady-state calculation to figure out whether the system can maintain the system pressure without imposing the boundary pressure.The PZR level is controlled during the steadystate calculation by a FILL component connected at the bottom of the PZR.It is expected to have zero flow at the level controller when the steady-state conditions are achieved.The model postulates the double-ended guillotine break at a cold leg.The broken cold leg is modeled by using three VALVE components, one for the cold leg and remaining two for the breaks.During the steady-state calculation, the VALVE components for the cold leg opens and the flow path to the break is closed.At the beginning of the accident, the VALVE component for the cold leg is closed to model the doubleend guillotine break.At the same time, two VALVE components for the break open simultaneously to establish the break flow for both ends of the broken cold leg.
The secondary side consists of two SGs, feedwater systems, and pressure boundary for steam.The operation of the main steam-line isolation valves (MSIVs) is implemented by closing both pressure boundary connections with the turbine trip signal.There are two FILL connections for each SG to model the feedwater to the economizer and downcomer, respectively.The feedwater injections are stopped once the turbine trip signal is initiated.

Description of 3D RPV model
TRACE uses three orthogonal-coordinate velocity-component momentum equations for the VESSEL component to describe the 3D flow behavior in either Cartesian or Cylindrical geometry [7].The mass and energy conservation equations for the VESSEL component were derived to cope with the change in the convective terms by three separated velocity components.All the field equations for the VESSEL component can have additional source terms which allows one or more 1D connections to any interface of a cell in the VESSEL component.This feature allows the connections to the hot/cold legs, ECCS (Emergency Core Cooling System), guide tubes, and so on.
The 3D RPV model used in this study is depicted in Fig. 2, which includes the ECCS consisting of the safety injection tanks (SITs) and safety injection pumps (SIPs).The 3D RPV model has been developed by modeling the RPV with single VESSEL component.As listed in Tables 1 and 2, the 3D RPV model consists of five radial rings, six azimuthal sections, and 35 axial levels.The first to third rings are reference to the reactor core and the fourth ring is designed to describe the core shroud bypass channel.The fifth ring models the downcomer.The azimuthal sections are used to describe the separated connections of two hot legs and four cold legs to the RPV.The axial nodalization was decided with reference to the internal structures listed in Table 2.The hot legs are connected to the upper plenum volumes, located at the azimuthal section 1 and 4 of the 4 th ring in the 30 th axial level.The cold leg connections are made to the rest of the azimuthal sections in the 5 th ring at the same axial level.The guide tubes in the core are modeled with the PIPE components separately.The hydraulic connection between each guide tube and the RPV is modeled by using crossflow connection indicated in blue dots in the figure in Table 2.In total, 1,050 hydraulic cells are used in the 3D RPV model.
The direct vessel injection (DVI) lines are connected to the 5 th ring of the 32 th axial level, where the connections are 2.1 m higher than the cold leg connections.The SIT of each DVI line is modeled by employing two individual flow paths: one for the standpipe and the other for the fluidic device.This model had been validated against the SIT blowdown test [17].The SIT operates once the upper downcomer pressure decreases below 4.03 MPa.Four FILL connections are employed to model the SIP injections.The SIPs are modeled to inject the ECC water based on the pressure dependent mass flow rate table with 40 s delay after the safety injection actuation signal (SIAS).

Reactor core model
Each axial plane of the reactor core is modeled by three radial rings and six azimuthal sections.Thus, the plane consists of 18 zones and the size of each zone is decided by grouping of the fuel assemblies considering the power distribution, as listed and depicted in Table 3 and Fig. 3, respectively.Every fuel assembly is modeled separately in the reactor core of the standard APR1400 model to cope with the detailed multi-physics analysis using TRACE-PARCS coupling.APR1400 has 241 fuel assemblies in the core.However, the assemblies numbered from 113 to 129 in Fig. 3 are divided in half in the model for grouping and thus, 17 additional heat structures are used to describe the halved fuel assemblies.The fuel assemblies denoted by red number in Fig. 3 feature the highest assembly power and they contain the hot pins with the highest pin power in the reactor core.As the behavior of the hot pin is of interest for LBLOCA analysis, additional eight heat structures are employed to model the hot pins.Thus, in total, 266 heat structures are modeled for the reactor core and this model is named as the full-assembly model.
The full-assembly model described above describes the reactor core in detail to consider the coupled analysis with the 3D neutronics code, PARCS.However, such a detailed model is not necessary when the analysis is conducted with the point kinetics model.Especially, the model will require lots of computational costs for reflood calculation during the LBLOCA analysis where the fine meshing is made to track the quench front location precisely.However, the result of the figure of merit is expected to show insignificant difference because the peak cladding temperature (PCT) is the figure of merit in the LBLOCA analysis and it is likely to occurs at the hot pin in general.Thus, an alternative core model has been developed by averaging the assemblies in each zone into one heat structure component, as depicted in Fig. 4, for better efficiency during the LBLOCA analysis.Since the number of cells in each axial plane of the core region is 18 and each cell should include an averaged fuel assembly, 18 heat structures are used to model the averaged fuel assemblies in the reactor core.Considering additional 8 heat structures for the hot pins, in total, 26 heat structures are employed to describe the fuel assemblies in the reactor core of the simplified 3D model and this model is named as the zone-averaged model.As listed in Table 3, the averaged peaking of each cell in the zone-averaged model is made by dividing the sum of peaking by the number of assemblies in each section.

Physical models for LBLOCA analysis
For reflood calculation, the fine mesh and axial conduction models are activated for the entire fuel assemblies in both 3D models.The fine-mesh rezoning is restricted by three permanent nodes per each axial node.Considering the oxidation of the   cladding during the accident, the Cathcart-Pawel model [18] is used as a metal-water reaction model.The dynamic gap conductance model derived on the basis of the model in FRAPCON 3.4 [19] is activated to take into account the variation of the gap geometry during the accident.The Ransom-Trapp model [20], which is the default critical flow model of TRACE, is applied to the break to model the break flow and the discharge coefficient of the model is set to 1.0.

Steady-state
A null-transient calculation with each model has been conducted for 5,000 s in problem time using TRACE V5.0 Patch 6 in order to achieve steady-state conditions.Table 4 summarizes the plant parameters from the steady-state calculations.Both models predict the reference plant conditions appropriately, except for the steam flow rate which has a relative error of 2.7 %.When comparing the heat transferred from the fuel rods to the coolant based on the enthalpy difference and core mass flow, it was found that the fullassembly model with 266 heat structures resulted in a smaller heat transfer of 3,892 MW, whereas the zone-averaged model with 26 heat structures resulted in 3,948 MW.This indicates that there is a loss of the thermal power in both models and, especially, the fullassembly model loses the power more.In order to figure out the root cause of such losses in the thermal power, additional investigation has been carried out by conducting the steady-state calculations with six different cases which have different number of the heat structures for the fuel rods maintaining the same total power and peaking conditions to the reference model.As depicted in Fig. 5, the results show that the error of the thermal power increases in proportion to the number of heat structures, and the trend of this errors could be correlated with a specific function.By performing additional calculations with 68 and 110 heat structures and comparing both results with the derived function, it has been confirmed that the function could predict the losses of additional cases.This indicates the dependency of the error with respect to the number of heat structures.The same discrepancy in the thermal power is indicated when a constant power is applied.Thus, it is concluded that the error in the thermal power is dependent on the number of associated heat structures.However, the impact of the loss of the thermal power should be very limited in the LBLOCA analysis because the thermal power decreases drastically as soon as the LBLOCA initiates.

LBLOCA calculation
The event chronologies listed in Table 5 represent the major system responses from both models in the early stage of the postulated accident.As soon as the accident initiated, the void fraction was increased in the reactor core instantaneously because of drastic flashing in the core and loss of the coolant through the break.Thus, the reactor power dropped down to the decay heat level due to large negative reactivity by the coolant density change, as depicted in Fig. 6.The RCP trip occurred with the accident initiation as the loss of offsite power (LOOP) was assumed as a coincident event.All the RCPs underwent the coastdown after the RCP trip, except that the one in the broken side accelerated in the early stage of the accident, as shown in Fig. 7.This discrepancy was made by the established break flow during the blowdown phase depicted in Fig. 8.During the blowdown phase, most of the coolant was discharged through the broken cold leg.As indicated in Fig. 8, larger amount of broken flow was discharged through the break at vessel side.However, about 33 % of total break flow was established through the break at pump side during the blowdown phase and this resulted in the further acceleration of the RCP at the broken cold leg.The primary pressure decreased rapidly by the loss of inventory in the reactor system and reached the setpoints of low pressurizer pressure trip (12.48 MPa) and SIT injection (4.03 MPa)  at 9.0 s and 14.7 s, respectively, as depicted in Fig. 9.The SIT injection was initiated automatically when the pressure went below the setpoint.However, the SIP injection did not start until the SIT injection because of the LOOP.Since the LOOP was assumed as a coincident event, the emergency diesel generators (EDGs) started its operation with the accident initiation and a delay of 40 s was assumed to consider the time required for the full functionality of the EDGs.In addition, the failure of an EDG was assumed as the single failure and an EDG is connected to two SIPs in the reference APR1400 design.Thus, SIP-1 and SIP-3 were assumed to be available during the calculation and actual injections using both SIPs started at 40 s from the accident initiation.General plant behaviors obtained by both models indicated no significant difference.When comparing the predicted maximum cladding temperatures depicted in Fig. 10, both models revealed similar behaviors as well.The predicted PCTs were 997.68 K and 998.96K for full-assembly model and zone-averaged model, respectively.Thus, it was concluded that both models are consistent with each other and the zone-averaged 3D model with 26 heat structures could be employed for LBLOCA analysis for better efficiency as listed in Table 6.Following the sequence of the accident, further evaluation has been made focusing on the TH behaviors within the RPV.The details are described for blowdown and reflood phases, respectively, based on the results of the full-assembly model with 266 heat structures.

Blowdown period
In Fig. 11, two-dimensional void and fuel temperature distributions of RPV consisting of intact and broken axial planes is depicted.The figure clearly shows that it took no more than 3 s for the active core to be dried out.Most of the assemblies were heated up because of the loss of coolant in the reactor core.Especially, the excessive heat up was observed in the second ring which has the highest average peaking.During this period, the steam flow was established upward to the upper plenum as depicted in Fig. 12.This upward steam flow interfered the water in the upper head to drop into the active core.However, as the water inventory was continuously exhausted, the steam flow became subsequently diminished and depleted within 5 s as in Fig. 13.Since then, the water in the upper head started to drop into the active core.As depicted in     Fig. 14, the downward flow of water has been made twice during the blowdown phase.As show in Fig. 15, the flow comes from the upper head through control element assemblies (CEA), at first.However, the flow was diminished in 5 s as the water inventory of the upper head was depleted.Since then, flows back from the hot legs were established to the upper plenum as depicted in Fig. 16, and this yielded the downward flow to the active core.Meanwhile, the figure clearly shows that more flow was established from the intact side.This asymmetry was attributed to the difference in the flows from both hot legs, as shown in Fig. 17.The figure revealed that the hot leg of the intact side had more flow than the broken side.Since the hot leg of the broken side was directly connected to the break, the coolant in the broken side could flow over the utubes more than one in the intact side.Thus, the larger flow back to the upper plenum could be established from the intact hot leg, as shown in Fig. 18, and this resulted in an asymmetric downward flow distribution to the core.This asymmetric downward flow into the core brought an asymmetric cooling, and, as a result, fuel assemblies in the intact side core has relatively lower temperature distribution compared to the broken side, as depicted in Fig. 19.

Reflood period
The blowdown quench ended up within 20 s after the break.Since then, the reactor core became dried out again due to the decay heat.Meanwhile, the injection of cooling water from ECCSs subsequently cooled down the downcomer, and, consequently, it refilled the lower plenum as depicted in Fig. 20.The core was still empty during this period and thus, the entire core became heated up as depicted in Fig. 10.With the continuous operation of ECCSs, the core reflood starts at around 50 s after the break.During the reflood, the subsequent recovery of the core inventory was made as given in Fig. 21.The contact of the reflooding water accompanied boiling at the heated surface of the fuels and the consequent quenching was developed as shown in Fig. 22.One should be noticed that the development of the quench front was not made uniformly.Most of the regions in the third ring showed rapid development while the second ring resulted in delayed quenching.This is because the third ring has the lowest average peaking and most of the stored heat in the third ring was sufficiently removed during the blowdown, but the second ring was not, as depicted in Fig. 19.Consequently, this has made the difference that the third ring is easier to be quenched whereas the second ring is not.In addition, the results depicted in Figs.23e25 indicate that there was an azimuthal quenching difference as well.All the ring sections consistently featured asymmetric quenching behavior such that the intact side showed rapid development whereas the broken side showed most delayed quenching.Such differences in the reflood quenching as well as the nonuniform blowdown quenching are highly connected to the peaking of each fuel assembly and the average peaking of each zone.In addition, such nonuniform behaviors during the accident could be captured by only with the 3D component.Thus, it is recommended to employ a multidimensional RPV model with a detailed fuel description for a realistic safety analysis with the consideration of the RPV design characteristics such as radial power distribution and hot pin locations.

Conclusion
The 3D component allows to capture the multi-dimensional phenomena such as asymmetric flow and temperature distributions in accident analyses.In order to perform a realistic safety analysis with the consideration of the multi-dimensional phenomena, a standard 3D analysis model of APR1400 for TRACE has been developed in this study.The model employs a VESSEL    component in TRACE for modeling the entire RPV.The standard model, namely full-assembly model, has been developed with the assembly-wise core description considering the coupled analysis with PARCS for SLB and RIA.Additional model, namely zoneaveraged model, has been developed on the basis of the standard model in order to improve the calculation performance in LBLOCA analyses.A comparative LBLOCA analysis reveals that general plant behaviors obtained by both models showed no significant difference.Based on the comparative analysis result, it is concluded that both models are consistent with each other, and thereby, the zoneaveraged model could be used for LBLOCA analysis with better computational performance.The detailed results from the standard

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Fig. 3 .Fig. 4 .
Fig. 3. Full-assembly model with the location of hot pins (in red).(For interpretation of the references to colour/colour in this figure legend, the reader is referred to the Web version of this article.)

Fig. 5 .
Fig. 5.The tendency of power balance error according to the number of heat structures.

Fig. 11 .
Fig. 11.Distribution of void fraction and fuel temperature: 2.4 s after break.

Fig. 18 .
Fig. 18.Distribution of void fraction and mixture flow across the hot legs: 13.4 s after break.

Fig. 19 .
Fig. 19.Distribution of void fraction and fuel temperature: 21.4 s after break.

Fig. 20 .
Fig. 20.Distribution of void fraction and fuel temperature: 42.5 s after break.

Table 3
Zone-averaged peaking in each section.

Table 4
Steady-state calculation results of zone-averaged and full-assembly models.

Table 5
Event chronology of each model during LBLOCA analysis.

Table 6
Computational performance comparison for LBLOCA analysis.