RAST-K v2—Three-Dimensional Nodal Diffusion Code for Pressurized Water Reactor Core Analysis

The RAST-K v2, a novel nodal diffusion code, was developed at the Ulsan National Institute of Science and Technology (UNIST) for designing the cores of pressurized water reactors (PWR) and performing analyses with high accuracy and computational performance by adopting state-of-the-art calculation models and various engineering features. It is a three-dimensional multi-group nodal diffusion code developed for the steady and transient states using microscopic cross-sections generated by the STREAM code for 37 isotopes. A depletion chain containing 22 actinides and 15 fission products and burnable absorbers was solved using the Chebyshev rational approximation method. A simplified one-dimensional single-channel thermal-hydraulic calculation was performed with various values for the thermal conductivity. Advanced features such as burnup adaptation and CRUD modeling capabilities are implemented for the multi-cycle analysis of commercial reactor power plants. The performance of RAST-K v2 has been validated with the measured data of PWRs operating in Korea. Furthermore, RAST-K v2 has been coupled with a sub-channel code (CTF), fuel performance code (FRAPCON), and water chemistry code for multiphysics analyses. In this paper, the calculation models and engineering features implemented in RAST-K v2 are described, and then the application status of RAST-K v2 is presented.


Introduction
Nodal diffusion codes are generally employed in the second part of the conventional two-step approach in the traditional procedure for nuclear reactor core analysis and design. Several nodal diffusion analysis codes such as SIMULATE [1], PARCS [2], MASTER [3,4], ASTRA [5], SCOPE [6], ROCS [7], and ANC [8] have been developed and used for light water reactor core design in recent years. Despite the existence of many nodal diffusion codes that exhibit high performance and accuracy, a novel nodal code was developed in the present study to incorporate state-of-the-art calculation methods. In particular, a new nodal diffusion code is needed to efficiently implement engineering features and allow coupling with other nuclear engineering codes for high-fidelity multiphysics simulations. The novel code developed here is expected to efficiently handle the transition of the regulatory environment, including the related reactor operation issues represented as axial offset anomalies (AOAs), and strengthen safety criteria.
The RAST-K v1 [9,10] was developed for evaluating the static and dynamic control abilities of control rods during dynamic control rod measurements (DCRMs). In particular, it was employed in a real-time core simulator for 21 pressurized water reactors (PWRs) operating in South Korea. In

Calculation Models
In this section, we explain in detail the calculation models and methods implemented in RAST-K v2, such as the nodal diffusion method, cross-section model, internal TH solver, and fuel depletion method. Since 2019, the capability to model hexagonal geometry has been newly implemented in RASTK-HEX [12][13][14], thereby enabling the RAST-K v2 to model both square and hexagonal lattices. This paper focuses only on the introduction of RAST-K v2 to model square lattices, and the RASTK-HEX will be introduced in the future.

Nodal Diffusion Analysis
To implement the non-linear nodal method in RAST-K v2, the multi-group diffusion equation is generally solved using the two-node nodal method with the coarse-mesh finite difference (CMFD) method [15]. Since the assumption of linear flux variation is not accurate for a large node size, the non-linear nodal method involves iterations entailing the solution of the CMFD and two-node nodal problems. The CMFD problem incorporates the global coupling of the nodes, whereas the two-node nodal problems incorporate local higher-order coupling.
In the two-node nodal problem, the interface current is obtained by employing the analytic nodal method (ANM) or the nodal expansion method (NEM) [16,17]. In the RAST-K v2, the formulation of the two-node nodal method is based on the unified nodal method (UNM) [18]. The flux solution of the previous CMFD problem is employed as the boundary condition during the two-node nodal calculation. The two-node problem is solved for each interface of every node and in all directions to provide an improved estimate of the net current. The interface current is then used to determine the corrective nodal coupling coefficient. The next CMFD matrix is formulated using the aforementioned relationship for the current.
The advantage of the UNM formulation is that the nodal basis functions can be optionally chosen for a specific condition that improves the nodal calculation stability and performance. As described 1.
Moderator temperature: 12 points, ranging from 293.6 K to the input moderator temperature + 25 K 3.
Geometry: 4 points for geometry changes (generally for control rod) STORA is a linking code that processes STN files into a cross-section file for use in RAST-K v2. The code reads the following data from STN files [11]: • State points: burnup, fuel and moderator temperature, boron concentration, control rod • Group-wise assembly/corner discontinuity factors • 1-group power distribution (gamma smeared power distribution) • 1-group form function for pin power reconstruction The group constants can be expressed in terms of the fuel temperature, moderator temperature, boron concentration, position of the control rod, and burnup, as shown below: where S is the state index, BU is the burnup point, Σ is any cross-section, Σ b is the base cross-section, S b is the base state point, dΣ is the deviation of the cross-section from the base state, and g is the given state. The cross-section is interpolated separately from S and BU. Figure 1 shows the burnup-dependent cross-section model.
where is the state index, is the burnup point, is any cross-section, is the base crosssection, is the base state point, is the deviation of the cross-section from the base state, and is the given state. The cross-section is interpolated separately from and . Figure 1 shows the burnup-dependent cross-section model.  The burnup calculations corresponding to the base state were performed using a fine burnup grid, whereas the restart calculations were performed only at selected burnup steps; i.e., the branch calculations were based on a coarse burnup grid. Therefore, the base state cross sections are given on a fine burnup grid, whereas the cross-sections for the branch states are provided on a coarse grid. The expressions for Σ b and dΣ are given as follows: Here, the weighting factors for each branch can be obtained using the 3D/2D Lagrange interpolation, leading to the following expression for dΣ:

Internal TH Solver
The internal TH model used in RAST-K v2 entails a simplified TH solution method. In particular, it solves the axial heat convection and radial heat conduction equations only for PWRs [22]. The model assumes the following:

1.
The coolant flow is parallel to the channel; that is, the cross-flow is ignored.

3.
The power produced by the fuel rods within a node is deposited into the coolant in that node. 4.
The pressure drop across the core is assumed to be negligible.

5.
All water properties are evaluated at a single pressure.
Under the single-phase formulation, no boiling is considered in RAST-K v2. In addition, since the momentum equation is not solved under the assumption of constant pressure, only the mass continuity and energy conservation equations are solved with regard to the flow problem. The constitutive relations, which are required to close the field equations, are provided in the form of a polynomial at a given pressure. Based on the aforementioned assumptions, two processes are considered in RAST-K v2, namely, heat convection in the coolant and radial heat conduction in the fuel.
For each TH node, the radial temperature distribution within a fuel rod should be calculated to determine the Doppler temperature. The finite difference scheme is employed to obtain the radial temperature distribution and its temporal evolution. The theta method is used for the time differencing of both the fuel heat conduction and coolant heat convection. Axial heat conduction is neglected in the solution process because it is considerably smaller than radial heat conduction. At the fuel pin wall, the radial heat conduction and axial heat convection are coupled through the heat flux terms, and this poses a nonlinear problem. To resolve the nonlinearity, the heat conduction equation is decoupled Energies 2020, 13, 6324 5 of 21 using a past iterative value of the bulk temperature of the coolant. Consequently, heat conduction and heat convection problems can be solved separately.
In the radial direction, a node is divided into the pellet, gap, and cladding regions. A second-order accurate finite difference approximation was used for spatial differencing with the central differencing scheme. By collecting the difference equations for each region in the fuel pin, a linear system for the radial heat conduction solution can be obtained. The coefficient matrix is determined using the thermal properties of the fuel and heat transfer coefficients. The unknowns are the temperatures at each of the mesh points within the fuel rod in the radial direction, while the coolant temperature is the boundary condition. The properties of the coolant are known only at the bottom of the channel. The heat fluxes at the fuel wall were determined from the previous step of the fuel heat conduction calculations. Therefore, the TH properties of the channel can be obtained by sweeping from the bottom to the top of the channel.
Various thermal conductivity models can be selected for the internal TH model of the RAST-K v2. Briefly, the fuel and cladding thermal conductivities and gap conductance were obtained from the fitting function of conventional PWRs. The burnup-dependent thermal conductivity model from FRAPCON-4.0 [23] was implemented in RAST-K v2 to model the thermal conductivity degradation effect. Thus, fuel conductivity can be expressed as a function of the burnup, temperature, and gadolinia content. Similarly, the internal TH solver can consider the thickness of the ZrO 2 layer with regard to the temperature drop across the coolant and the outer cladding surface. In addition, the thermal conductivity fitting functions based on the experimental data for accident tolerance fuel (ATF) were also implemented.
The internal TH solver has been verified with other physics codes that can calculate the TH profile. For verifying the coolant temperature, a one-pin model is used for calculations with a fixed power shape during the hot state, and the result is compared with the sub-channel code CTF [24]. The maximum difference between the coolant temperature from the internal TH solver and CTF is smaller than 1 K at the same outlet temperature. For verifying the fuel temperature, a one-pin model is used for calculations with a fixed power shape and coolant temperature as the boundary condition. The maximum difference in the fuel centerline temperature from the internal TH solver and FRAPCON is smaller than 85 K owing to the different thermal conductivity information. At the same time, the maximum difference in the volume-averaged fuel temperature is smaller than 16 K. These calculations verify that the internal TH solver can provide a temperature profile even with a simple assumption for the hot state.

Fuel Cycle Analysis
For the fuel cycle analysis, RAST-K v2 adopts the micro-depletion method for important nuclides. The nuclide number densities of 22 actinides, 12 fission products, and some gadolinium isotopes are considered within the depletion chain of RAST-K v2. In PWR, the heavy nuclide chain is largely composed of U-235 chains and U-238 through Np-237 and Pu-238. To model the U-235 chain and the plutonium build-up from the alpha decay of curium, RAST-K v2 heavy nuclide chain includes 22 actinides from U-235 to Cm-244. To track the fission product nuclides that have large absorption cross-sections, the iodine/xenon chain and samarium chain are independently included in the fission product chain. The depletion chain in RAST-K v2 is internally verified by comparing it with depletion calculation using STREAM for various types of fuel assemblies with a 1% difference in the number density during burnup. The neutron transmutation equation can be expressed as the matrix form of the Bateman equation [25], which, in turn, is solved using the Chebyshev rational approximation method (CRAM) [26,27]. There are several techniques for reducing the time discretization error during burnup calculations [28][29][30]. RAST-K v2 adopts three schemes: the explicit Euler method (only predictor), implicit Euler method (Full predictor/corrector), and semi-predictor/corrector methods, and the user can choose any of the three schemes. Figures 2 and 3 show the depletion chain of the 22 heavy nuclides and 12 fission products used in the RAST-K v2 depletion calculation. approximation method (CRAM) [26,27]. There are several techniques for reducing the time discretization error during burnup calculations [28][29][30]. RAST-K v2 adopts three schemes: the explicit Euler method (only predictor), implicit Euler method (Full predictor/corrector), and semipredictor/corrector methods, and the user can choose any of the three schemes. Figures 2 and 3 show the depletion chain of the 22 heavy nuclides and 12 fission products used in the RAST-K v2 depletion calculation.    The information regarding the decay constants and fission yield, which forms the depletion matrix, is transferred via a cross-section file obtained from STREAM. The neodymium chain is included in the fission product chain to improve the accuracy of the Sm-149 number density. The fission yield of neodymium is calculated in STREAM by summing the cumulative fission yields and the portions generated from other nuclides. It should be noted that the I-135/Xe-135 chain is not included with other fission product chains. Xe-135 is the most crucial fission product required for tracking owing to its large absorption cross-section and relatively short half-life. Hence, the I-135 and Xe-135 isotopes are treated separately from other depletion chains and are explicitly tracked.
In the calculation, the chain of gadolinium isotopes (Gd-154, Gd-155, Gd-156, Gd-157, and Gd-158) is replaced by an effective gadolinium isotope [31][32][33]. Owing to the effect of the absorption reaction rate along the burnup, which arises from the spatial self-shielding effect, it is better to solve for the effective gadolinium isotope, instead of solving for a direct chain of all gadolinium isotopes. The number densities and microscopic absorption cross-sections for the effective gadolinium isotopes are defined as follows: The information regarding the decay constants and fission yield, which forms the depletion matrix, is transferred via a cross-section file obtained from STREAM. The neodymium chain is included in the fission product chain to improve the accuracy of the Sm-149 number density. The fission yield of neodymium is calculated in STREAM by summing the cumulative fission yields and the portions generated from other nuclides. It should be noted that the I-135/Xe-135 chain is not included with other fission product chains. Xe-135 is the most crucial fission product required for tracking owing to its large absorption cross-section and relatively short half-life. Hence, the I-135 and Xe-135 isotopes are treated separately from other depletion chains and are explicitly tracked.

Engineering Features
In the calculation, the chain of gadolinium isotopes (Gd-154, Gd-155, Gd-156, Gd-157, and Gd-158) is replaced by an effective gadolinium isotope [31][32][33]. Owing to the effect of the absorption reaction rate along the burnup, which arises from the spatial self-shielding effect, it is better to solve for the effective gadolinium isotope, instead of solving for a direct chain of all gadolinium isotopes. The number densities and microscopic absorption cross-sections for the effective gadolinium isotopes are defined as follows: Energies 2020, 13, 6324

Engineering Features
In this section, we describe the various miscellaneous features implemented in the RAST-K v2 for multi-cycle simulation.

Multi-Cycle Simulation
To perform a multi-cycle core follow calculation, the restart, fuel shuffling, and rotation functions should be capable of reloading the fuel material and burnup information from the previous cycle. In the restart file, the material number densities, burnup information, and calculation options are stored, and this file is read before the core calculation. The restart information can be folded and unfolded according to the core symmetry with periodic or reflective boundary conditions. Owing to the lack of core depletion history, the jump-in function is implemented to approximate the material information of the burned fuel by interpolating the number density based on the assembly burnup.
During the core follow calculation of commercial PWRs, the critical state is maintained using the boron concentration. Thus, the boron concentration indicates the target eigenvalue. The critical boron concentration (CBC) was estimated as follows: where i is the CBC feedback iteration index, k i is the eigenvalue at i, and C is a conversion factor from the eigenvalue to the boron concentration. Normally, the conversion factor is set to 10 4 assuming that the reactivity of the 1 ppm boron concentration can be controlled to 10 pcm reactivity. When the CBC oscillation occurs, the conversion factor is obtained by linear extrapolation using the CBC and eigenvalues at the (i − 1) th and i th iterations. The control rod and power level can also be employed for estimating the critical state. The equilibrium xenon can be evaluated based on the depletion chains for I-135 and Xe-135 at the equilibrium state with a given power level. The core mass flow rate can be determined for the target average moderator temperature. RAST-K v2 is capable of executing branch calculations, whereby nodal calculations can be performed by changing the options or core conditions at a specific time step during steady-state calculations. It is possible to change the calculation mode, calculation options such as those regarding xenon, samarium, and gadolinium, feedback options, and core conditions such as the boron concentration, temperature, and control rod position. The branch calculation capability facilitates the computation of the reactivity worth of temperature change, boron, xenon, samarium, gadolinium, and control rods, allowing the safety margin to be easily assessed during the core design.

Pin Power Reconstruction
Since RAST-K v2 employs nodal methods based on the homogenized nodal cross-section, the resulting intranodal flux distribution cannot provide any local heterogeneity within an assembly. Thus, the pin power information is obtained through a pin power reconstruction process. Currently, only two-dimensional reconstruction is available to obtain the pin power distribution on each axial plane. The pin power reconstruction model employed in RAST-K v2 assumes the separability of the global flux (homogeneous intranodal flux) and local flux shapes (heterogeneous form functions), similar to most other nodal codes [34,35].
The use of a spatially homogenized cross-section implies that the nodal model is capable of providing only smoothly varying flux distributions within each assembly. Using the separability assumption, the pin-by-pin power distribution within an assembly can be approximated via the product of the global homogenized power distribution and a local heterogeneous form function. The form function should account for assembly heterogeneities caused by water holes, burnable absorber pins, Energies 2020, 13, 6324 8 of 21 enrichment variations, etc. In the RAST-K v2, the power form functions of all types of fuel assemblies are obtained from the prior lattice calculations. Further, they are stored in the interpolation table as a function of the state as a cross-section.
The homogeneous intranodal flux is obtained using analytic functions, each of which is a solution of the two-dimensional neutron diffusion equation, as the basic function of the expansion. The determination of the intranodal flux becomes a Dirichlet problem in which a two-dimensional partial differential equation is solved exactly with the function values specified at the boundary. To solve this problem, four surface average currents corresponding to two-nodal sweeping and four-corner flux are employed with regard to the boundary conditions. The corner flux can be calculated by solving the corner point balance equation and multiplying it with the corner discontinuity factor.
The pin power reconstruction module of the RAST-K v2 was verified using STREAM transport and SIMULATE-3 nodal calculations with pin power reconstruction. For the hot full-power OPR-1000 core calculation, the root mean square (RMS) difference of the pin power from the RAST-K v2 is 7.2% and 5.7%, compared with the STREAM core and SIMULATE-3 calculations, respectively.

Parameter Edits
As the reactivity of commercial PWRs is controlled by boron, the B-10 in the coolant undergoes a depletion during the cycle. A B-10 depletion model was implemented in the RAST-K v2 to compare the measured CBC values with the calculated values directly. The micro depletion solver does not directly track the B-10 nuclide and thus, the RAST-K v2 only corrects the calculated CBC value through the B-10 depletion model.
An ex-core detector signal can be obtained using the results of core calculations, such as the results of the flux and cross-section of the detector material. The pre-generated detector response functions of OPR-1000 and Westinghouse 3-loop-type reactors are stored in RAST-K v2. The ex-core detector signal is calculated by multiplying the detector response function with the flux from the core calculation.
During cross-section generation via lattice calculation, the in-core detector cross-section and box-to-pin flux ratio can be generated for all branch cases. Thus, the information is processed using the same cross-section feedback routine. By multiplying the cross-section of the detector material, the box-to-pin flux ratio, and the assembly flux from the core calculation, the in-core detector signal can be obtained. Thus, depletion calculations for the detector material can be performed. RAST-K v2 has depletion chains for rhodium, cobalt, vanadium, and silver as the detector materials, which are usually employed as in-core self-powered neutron detectors (SPNDs).
To calculate the departure from nucleate boiling ratio (DNBR) during the steady and transient states, the critical heat flux (CHF) can be calculated based on the W-3 CHF correlation [36], which is widely used in the sub-channel code for PWRs. The CHF is calculated using information such as the enthalpy, equivalent diameter, and local steam quality, which is obtained during the internal TH calculation.

Burnup Adaptation
It is crucial to predict the reactor parameters accurately during reactor operation, even though the operating conditions may change owing to unpredictable practical conditions. It is difficult to quantify the effects of unpredictable conditions that should be considered in the calculations, and the discrepancy between the prediction data and actual data will cause errors and biases in modeling the next reactor cycle. The burnup adaptation model is implemented in the RAST-K v2 to compensate for the effect of past unpredictable core conditions [37].
The burnup of the nodes in each fuel assembly was modified using the shape function and burnup multiplier. The shape function in the burnup adaptation model is sinusoidal, and the number densities of major isotopes used in RAST-K v2 can also be modified with the corrected burnup. The multipliers of the burnup and shape function are searched to minimize the error of the CBC and axial shape index (ASI) and the error of the two-dimensional radial power distribution with respect to the reference data. Prior research has revealed that the assembly power and CBC are considerably more sensitive to the burnup multiplier while the ASI is more sensitive to the shape function multiplier. The reference data from the measurement were utilized for the burnup adaptation, and the numerical results of a typical OPR-1000 core analysis were successfully demonstrated.

CRUD Modeling
The deposits of corrosion products, referred to as chalk river unidentified deposits (CRUD), have a significant impact on reactor operation. This can lead to a CRUD-induced power shift (CIPS), which is an operational issue, and CRUD-induced local corrosion (CILC), which is a safety issue. It is difficult to predict CRUD deposition by considering the local thermohydraulic and chemical processes. Hence, RAST-K v2 includes a simple CRUD modeling capability that considers the power shift and fuel temperature increment owing to CRUD deposition.
The simple CRUD modeling function of RAST-K v2 expresses the thickness of the CRUD as a function of the axial shape function. It is assumed that the CRUD thickness is 80 µm at a fuel assembly burnup of 60 MWd/kg. Figure 4 shows the axial shape function used to calculate the thickness of the CRUD. These fitting functions were obtained from the measurement data.  The change in the thickness of the CRUD can be expressed as follows: where is the axial shape function, as shown in Figure 4, and is the average fuel assembly burnup. The threshold of the CRUD deposition can be set to prevent the deposition of CRUD from reaching the range corresponding to the low burnup condition. The increment in the fuel temperature can be expressed as follows: where ′′ is the average heat flux in each TH node and is the thermal conductivity of the CRUD, which is 1.5 W/m-K. The CRUD deposition affects the neutron cross-section in two different ways: an increase in the fuel rod temperature and an increase in the concentration of boron deposited on the outer surface of the cladding. The boron deposition rate in the CRUD can be expressed as follows: where is the burnup step index, is the number density of boron in the CRUD at the burnup step, ∆ is the volume of newly deposited CRUD, is the volume of the node, is the porosity of the CRUD, and ∆ is the time interval. To consider the reduction of the CBC by feeding the CRUD layer, the boron number density in the CRUD is corrected as follows:  The change in the thickness of the CRUD can be expressed as follows: where S CRUD is the axial shape function, as shown in Figure 4, and BU is the average fuel assembly burnup. The threshold of the CRUD deposition can be set to prevent the deposition of CRUD from reaching the range corresponding to the low burnup condition. The increment in the fuel temperature can be expressed as follows: where q is the average heat flux in each TH node and k CRUD is the thermal conductivity of the CRUD, which is 1.5 W/m-K. The CRUD deposition affects the neutron cross-section in two different ways: an increase in the fuel rod temperature and an increase in the concentration of boron deposited on the outer surface of the cladding. The boron deposition rate in the CRUD can be expressed as follows: where n is the burnup step index, N n B is the number density of boron in the CRUD at the burnup step, ∆V n CRUD is the volume of newly deposited CRUD, V node is the volume of the node, P CRUD is the porosity of the CRUD, and ∆t n is the time interval. To consider the reduction of the CBC by feeding the CRUD layer, the boron number density in the CRUD is corrected as follows: where N 0 B is the unadjusted boron density in the CRUD and CBC n is the CBC at the burnup step n. In addition, the boron number density in the CRUD can undergo depletion during core depletion, and the change in the boron number density in the CRUD is expressed as follows: where σ n a φ n is the absorption reaction rate of boron at the burnup step n. If the measured power is provided, the amount of CRUD deposition can be precisely predicted by adjusting the measured ASI. Figure 5 shows the CBC and ASI during core depletion obtained using the CRUD modeling capability of RAST-K v2. The axial offset anomaly (AOA) has occurred during the (N − 1) th cycle in practice. As demonstrated by the behavior of the ASI of the (N − 1) th cycle, the axial power shape cannot be predicted by RAST-K v2 without CRUD modeling. Thus, it will be difficult to model the N th cycle owing to the incorrect depletion history of the previous cycle contained in the restart file. Using the CRUD modeling function, the decrease in the CBC owing to the CRUD is shown for the (N − 1) th cycle, and the core depletion calculation for the N th cycle is performed more accurately. Hence, the ASI of the N th cycle is close to the measured values.

Spent Nuclear Fuel Analysis
The prediction of isotope inventory is important for the following reasons: (1) to ensure the safety of the management and transportation of spent nuclear fuel (SNF) and (2) to predict the radiation source terms such as activity, decay heat, and neutron and gamma intensities. In addition, the need for the analysis of SNF is increasing because the saturation issue of spent fuel pools is gaining importance in South Korea. Therefore, an SNF analysis module was developed based on the Lagrange non-linear interpolation method, which is employed for cross-section feedback by using the history indices and power correction factor. Lagrange non-linear interpolation is performed to predict the isotope inventory in the SNF based on the number density file pre-generated by the history depletion branch calculation in the STREAM. History indices (moderator temperature, fuel temperature, and boron concentration) were used to consider the whole core simulation conditions. The power correction factor was employed to adjust the power levels during fuel depletion. Moreover, the SNF module supports the calculation of source terms and decay heat for the assembly and pin levels.
The inventory calculation module has been validated with 58 pin samples discharged from five different reactors (Takahama-3, Calvert Cliffs-1, Turkey Point-3, TMI-1, and Obrigheim-1) [38]. A total of 1234 isotopes were compared with the corresponding measurements. Of these, 1034 isotopes (i.e., 84% of the total comparison data) matched the measurement data within a 20% relative error. The decay heat calculation module was validated using 58 fuel assembly samples discharged from five different reactors (Ringhals-2, Ringhals-3, Turkey Point-3, Point Beach-2, and San Onofre-1). The samples span a discharge burnup range from 19 to 51 MWd/kg, such that high burnup fuel assembly samples are also contained in the comparison samples. The relative errors of decay between the calculated and measured values lie within ±4% in the case of assembly-wise calculation and ±4.3% in the case of pin-wise calculation.
The application of SNF analysis in a commercial reactor was elucidated in [39]. A Westinghouse 2-loop reactor, which is one of the reactor types employed in South Korea, was used for comparison. Three operating scenarios were validated using STREAM-SNF: (1) assembly depletion during the N th , (N + 1) th , and (N + 2) th cycles, (2) assembly depletion during the N th and (N + 2) th cycles, and (3) assembly depletion only during the N th cycle. STREAM-SNF was employed for a code-to-code comparison, which showed that the results of the isotope inventory obtained using RAST-K v2 agree with those obtained using the STREAM-SNF within a relative error of 5%.

Applications
RAST-K v2 has been utilized for nuclear design, core analysis, coupling with other physics codes, and machine learning. In this section, the application status of RAST-K v2 is described.

Verifications and Validation
RAST-K v2 has been verified and validated using the operation history of Korean PWRs. A total of 114 cycles of 13 reactor cores, including the OPR-1000, APR-1400, Westinghouse 2-loop reactor (WH-2L), and Westinghouse 3-loop reactor (WH-3L) were employed for the core follow calculation. The core depletion calculation results are compared with the results obtained from commercial nuclear design code and measurement data for the verification and validation (V&V) of the steady-state depletion capability of the RAST-K v2. Several nuclear design codes are used in Korea's nuclear power plants. The DIT/ROCS [7,40] is used for the design of the initial core of the OPR-1000 and ARP-1400 reactors, while PARAGON/ANC [41] is used for the design of the reloaded core of the OPR-1000. The KARMA/ASTRA [5] is used for designing the reloaded core of the APR-1400, and PHOENIX-P/ANC [42] is used for designing the core of the Westinghouse-type reactor.
The CBC, ASI, and fuel assembly power distribution are compared. For the group constant generation, STREAM uses the combined ENDF/B-VII.1 and JENDL4.0 libraries for the neutron cross-section and depletion libraries. The two-dimensional fuel assembly depletion calculation is performed with the on-the-fly kappa library, semi-predictor/corrector method, and quadratic depletion for gadolinium isotopes. The resonance elastic up-scattering approximation of U-238 is used, and the critical spectrum is searched. The assembly discontinuity factor (ADF) is calculated using the homogeneous ADF method, and the radial reflector discontinuity factor preserves the fuel and reflector homogenized flux ratios. The fuel assembly pattern, which has an asymmetric pattern of burnable absorber rods, in a WH-2L reactor was developed from the cross-section generation by STREAM [43].
The RAST-K v2 calculation is performed with a three-dimensional quarter-core model, including four nodes for one assembly and 24 to 46 axial nodes. Equilibrium xenon calculations and CBC searches were conducted for every burnup point. The jump-in option is used in the absence of nuclear design information. The depletion calculation is performed based on the predictor/corrector method. The internal TH calculation is carried out by considering the effects of fuel thermal conductivity degradation and the generation of cladding ZrO 2 as a function of burnup. Table 1 presents a summary of the V&V error statistics of RAST-K v2. The CBC and ASI were compared for each burnup point, and the radial fuel assembly power distribution was compared at the beginning, middle, and end of cycles. The calculation results of the jump-in cycle and subsequent cycle are excluded from the comparison. The difference in the CBC obtained via RAST-K v2 and the measurements is smaller than 40 ppm, and this difference is smaller than 70 ppm for all reactor types when RAST-K v2 and the nuclear design code are compared. The absolute root mean square (RMS) error of the ASI between RAST-K v2 and measurements is smaller than 0.013 while that between RAST-K v2 and the design code is smaller than 0.017. The RMS error of the fuel assembly power between RAST-K v2 and the measurement is less than 3%, while that between RAST-K v2 and the design code is less than 5%. In summary, the CBC, ASI, and assembly power during steady-state depletion are validated within a maximum difference of 70 ppm, 0.017, and 5%, respectively.

Internal Loose Coupling
High-fidelity core analysis can be achieved through multiphysics coupling calculations. By coupling the CTF subchannel analysis code [24], it is possible to perform multi-channel modeling Energies 2020, 13, 6324 13 of 21 with cross-flow and obtain detailed coolant information to solve the equations for the conservation of mass, momentum, and energy for the two-fluid model. The FRAPCON fuel rod analysis code [23] calculates the fuel rod performance during core depletion based on the pellet-to-cladding heat transfer, mechanical deformation, pellet-to-cladding mechanical interaction, fission gas release, cladding oxidation, hydrogen pickup, pellet burnup, and power radial distribution models. The code uses the axial distribution of the linear power, coolant pressure, and outer surface temperature of the cladding as the feedback parameters to be iterated with the neutronics and thermal-hydraulics codes.
RAST-K v2 is coupled with CTF and FRAPCON to perform high-fidelity core depletion calculations. As CTF is a TH subchannel analysis code, its channel-centered model is changed to the pin-centered model to ensure consistency during coupling. The pin power information from the RAST-K v2 pin power reconstruction function is transferred to CTF and FRAPCON as input variables. The CTF solutions, such as the coolant temperature and pressure, are transferred to FRAPCON as a boundary condition for solving the heat conduction equation. Next, the coolant temperature from CTF and fuel temperature from FRAPCON are merged into node-wise information from pin-wise information and transferred to RAST-K v2 for cross-section feedback. For the multi-cycle calculation, the restart capability of FRAPCON is implemented by developing an interface program called FRAPI [44]. Figure 6 shows the flowchart and data exchange of RAST-K v2, CTF, and FRAPCON. The depletion calculation of RAST-K v2 is performed using a predictor/corrector scheme for higher-order prediction with linear extrapolation and quadratic interpolation [28]. Otherwise, the fuel depletion of FRAPCON is performed using only the predictor scheme with a given power. In the coupling, the fuel depletion process between RAST-K v2 and FRAPCON is matched by using the restart capability of FRAPI, so that the fuel temperature calculation of FRAPCON can consider the pin burnup. The multi-physics simulation for the core follow calculation was performed using four cycles of the OPR-1000 type reactor. In the equilibrium cycle, the core design parameters such as CBC, ASI, peaking factor, and reactivity worth obtained from multiphysics calculations are compared with those from the standalone RAST-K v2. Table 2 shows the difference in design parameters between the multiphysics version using RAST-K v2, CTF, and FRAPCON and the standalone RAST-K v2. The differences in CBC, ASI, and FA power were smaller than 10 ppm, 0.01, and 1.5%, respectively. The local pin-wise fuel and coolant temperatures can be obtained from the multiphysics simulation, but the temperature distribution is merged with the RAST-K v2 nodal mesh size. Hence, the fuel temperature for the fuel rods within a node is the same, and the pin-by-pin power reconstruction for The multi-physics simulation for the core follow calculation was performed using four cycles of the OPR-1000 type reactor. In the equilibrium cycle, the core design parameters such as CBC, ASI, peaking factor, and reactivity worth obtained from multiphysics calculations are compared with those from the standalone RAST-K v2. Table 2 shows the difference in design parameters between the multiphysics version using RAST-K v2, CTF, and FRAPCON and the standalone RAST-K v2.
The differences in CBC, ASI, and FA power were smaller than 10 ppm, 0.01, and 1.5%, respectively. The local pin-wise fuel and coolant temperatures can be obtained from the multiphysics simulation, but the temperature distribution is merged with the RAST-K v2 nodal mesh size. Hence, the fuel temperature for the fuel rods within a node is the same, and the pin-by-pin power reconstruction for pin-wise and node-wise heat transfer models is similar, resulting in a pin peaking factor difference lower than 1%. We note that since the fuel temperature averaging effect evidently decreases in the pin peaking factor owing to the Doppler effect, the calculated pin peaking factor value is conservative. The differences in each reactivity worth, such as the fuel temperature coefficient (FTC), moderator temperature coefficient (MTC), and control rod worth are lower than 4%, 8%, and 0.5%, respectively. In other words, the internal TH solver based on the equivalent pin model can accurately predict the core depletion results. A detailed analysis of the numerical results obtained from the multiphysics simulation will be performed in the future. To track the effect of CRUD deposition during core depletion, RAST-K v2 is coupled with the subchannel TH code VIPRE [45] and water chemistry code BOA [46]. The input file of VIPRE is automatically generated based on the node power distribution, which is calculated from RAST-K v2. After simulating VIPRE and BOA, the output file of BOA, which includes the CRUD thickness and boron mass in the CRUD, is read by RAST-K v2. The boron mass in the CRUD is used for the cross-section feedback, and the CRUD thickness is used to correct the increment in the fuel temperature. Owing to the absence of the source codes of VIPRE and BOA, the coupling is performed by handling the input and output files.
The coupled version of RAST-K v2, VIPRE, and BOA was utilized to evaluate the cycle depletion of OPR-1000, which has an AOA. Figure 7 shows the ASI values obtained from the calculation performed using the combined RAST-K v2, VIPRE, and BOA, nuclear design code (PARAGON/ANC), and measurements. As shown in the figure, it is possible to track the ASI by employing BOA during the AOA. Hence, the N th cycle depletion can be initiated by predicting the depletion history of the fuel in the (N − 1) th cycle.

External Loose Coupling
The intimate interaction of neutronics, thermal-hydraulics, mechanics, materials, and other sciences is of great significance in reactor core safety analysis. Numerical simulation of nuclear reactor cores is accomplished by linking existing engineering codes using various strategies generally divided into loose and tight coupling with respect to data exchange [47,48]. The loose coupling methods simulate each multiphysics phenomenon separately and exchange information between relevant equations through the boundary or source terms associated with the entire reactor core domain. The tight coupling methods implement data exchange on a subdomain-specific level, allowing more efficient acceleration techniques to create a cohesive application with much better performance. While the tight coupling strategy is superior in computational efficiency, the loose coupling is much more extensive as it is more manageable for code development.
In view of the full-scale multi-physics reactor core simulation, we have adapted the RAST-K v2 code to the external loose coupling interface to be linked with a multiple-purpose multiphysics system such as MPCORE [49,50]. In contrast to the internal coupling considered in Section 4.2, where the leading neutronics code incorporates TH and fuel performance as sub-modules, the external coupling aims to link engineering codes on an equal basis by a common software framework to facilitate the autonomy of the multiphysics modules. Consequently, for external coupling, any multiphysics module can be replaced by a more efficient or accurate one depending on the purpose, with some programming efforts applied to the new module interface, while other modules remain intact. On the other hand, the internal coupling of the master code cannot be dismissed, principally admitting improvement only through updates.
The external loose coupling interface consists of the following functions: allocating and deallocating the internal memory space of the module, saving and loading the internal memory for restarting the module, calculating the initial steady state, performing the next time step, and updating the feedback parameters. The last function transfers data from a module's internal memory space to the common memory space allocated by a multiphysics system framework shared by all linking modules. A detailed description of the interface functions and loose coupling algorithms of the RAST-K v2 code implemented for MPCORE has been presented in a previous study [49,50]. The RAST-K v2 configuration of the MPCORE multi-physics system was designed for safety analysis and uncertainty quantification of PWR depletion and transient cores.

Machine Learning
RAST-K v2 has also been employed to check the efficacy of the application of machine learning (ML) algorithms in the field of nuclear reactors [51][52][53][54]. The application of RAST-K v2 in ML includes the development of a loading pattern (LP) optimization algorithm based on the simulated annealing (SA) method, development of artificial intelligence (AI) models that predict nuclear reactor parameters, and development of a nuclear reactor core diagnostics system based on AI models. The utilization of RAST-K v2 in ML research is advantageous for several reasons: 1.
RAST-K v2 can provide an accurate solution for nuclear reactor cores with only a small difference compared to the measured data. Hence, RAST-K v2 can be used to develop ML models, which can be further used in a real nuclear power plant (NPP).

2.
Various reactor core models, even the ones whose operating conditions are unacceptable for safety reasons, can be modeled using RAST-K v2. Data on abnormal core conditions such as CRUD deposition can be employed to train ML models in a nuclear diagnostics system. 3.
RAST-K v2 can generate a substantial amount of nuclear data owing to the implemented acceleration methods, and such data can be used for efficiently training ML models within the stipulated time. 4.
RAST-K v2 can provide various nuclear reactor parameters for prioritizing LPs in SA algorithms or mimicking the core monitoring system of NPPs to develop a diagnostic system.
The safety and operating margins of core cycles are highly dependent on the method based on which fuel assemblies are loaded in the core. Therefore, the optimum LP, which maximizes the operating margin and also satisfies the safety criteria, is manually searched for in every core cycle. To find the optimum LP for which the core cycle is on beyond human limitations, an optimization algorithm based on the SA method is currently under development [51]. The algorithm randomly selects, estimates, and ranks an LP according to the objective function, which is defined based on several nuclear reactor parameters. The objective of each LP is determined by embedding RAST-K v2 in the LP optimization algorithm.
ML and artificial intelligence (AI) technologies are among the most extensively employed technologies owing to their effectiveness in various fields. Recent research has led to the development of ML models that predict nuclear reactor parameters such as cycle length, peaking factor, and maximum CBC [52,53]. The development of ML models requires a substantial amount of data regarding various LPs and nuclear reactor parameters corresponding to each LP for training. RAST-K v2 has been utilized in data generation with approximately one million core calculations, and the corresponding dataset was used to train various ML models. In a recent study, a convolutional neural network (CNN) model was built and trained. It was observed that the ML model exhibited high prediction accuracy for cycle length, peaking factor, and maximum CBC with average errors of 0.04%, 0.5%, and 0.1%, respectively. Figure 8 shows the accuracy of the cycle length predicted using the ML model. network (CNN) model was built and trained. It was observed that the ML model exhibited high prediction accuracy for cycle length, peaking factor, and maximum CBC with average errors of 0.04%, 0.5%, and 0.1%, respectively. Figure 8 shows the accuracy of the cycle length predicted using the ML model. The recent trends observed in light water reactors, such as high thermal power, obsolescence, long cycles, and water chemistry management, entail abnormal core operation conditions. Early diagnosis of abnormal core states such as CIPS and follow-up measures can reduce the costs incurred by such states. Thus, the ML algorithm has been applied in reactor core diagnostics systems to improve the current reactor core diagnostics, which is based on the operator's proficiency. The purpose of the ML model is to classify the core operating state using core operation data as reactor operators do in the main control room (MCR) in NPPs. In a recent study, a framework for developing ML models was established, and the feasibility of using ML in the core diagnostics system was confirmed [54]. In this study, RAST-K v2 was employed to generate the operation data, that is, incore instrument (ICI) signals, as shown in Figure 9. The data generation system AUTOGEN operates in three steps: 1. Generation of normal and abnormal core models by randomly perturbing the input parameters such as the positions of the individual control rod and multiplication factor of the CRUD accumulation model, 2. Running RAST-K v2 for all generated input files, The recent trends observed in light water reactors, such as high thermal power, obsolescence, long cycles, and water chemistry management, entail abnormal core operation conditions. Early diagnosis of abnormal core states such as CIPS and follow-up measures can reduce the costs incurred by such states. Thus, the ML algorithm has been applied in reactor core diagnostics systems to improve the current reactor core diagnostics, which is based on the operator's proficiency. The purpose of the ML model is to classify the core operating state using core operation data as reactor operators do in the main control room (MCR) in NPPs. In a recent study, a framework for developing ML models was established, and the feasibility of using ML in the core diagnostics system was confirmed [54]. In this study, RAST-K v2 was employed to generate the operation data, that is, in-core instrument (ICI) signals, as shown in Figure 9. The data generation system AUTOGEN operates in three steps:

1.
Generation of normal and abnormal core models by randomly perturbing the input parameters such as the positions of the individual control rod and multiplication factor of the CRUD accumulation model, 2.
Running RAST-K v2 for all generated input files, 3.
Extraction of target output parameters from output files and writing them in a single train dataset file with "csv" format.
Various reactor core abnormal conditions such as CRUD deposition and control rod mislocation were modeled by randomly perturbing the input parameters of RAST-K v2. The data were used to train various ML models: Naïve Bayes (NB), support vector machine (SVM), deep neural network (DNN), and random forest (RF). ML models, especially the DNN and RF models, classify the nuclear core abnormal conditions with a high prediction accuracy of over 88% depending on the core power and refer to the detailed analysis of ML methods. A detailed analysis of ML will be presented in the future. 3. Extraction of target output parameters from output files and writing them in a single train dataset file with "csv" format. Various reactor core abnormal conditions such as CRUD deposition and control rod mislocation were modeled by randomly perturbing the input parameters of RAST-K v2. The data were used to train various ML models: Naïve Bayes (NB), support vector machine (SVM), deep neural network (DNN), and random forest (RF). ML models, especially the DNN and RF models, classify the nuclear core abnormal conditions with a high prediction accuracy of over 88% depending on the core power and refer to the detailed analysis of ML methods. A detailed analysis of ML will be presented in the future.

Conclusions
This paper presents the calculation models, engineering features, and application status of the newly developed nodal diffusion code RAST-K v2. State-of-the-art neutronics, TH, and depletion solvers are adopted in this code. A non-linear nodal method based on UNM is used with multi-group cross-sections generated by STREAM for a specific case matrix. The single-phase energy-mass conservation equation is solved for a closed channel for TH feedback, and the micro-depletion calculation is performed using the CRAM solver for fuel cycle analysis. To support users in conveniently designing and analyzing reactor cores, several miscellaneous engineering features are implemented in RAST-K v2. In particular, burnup adaptation and simple CRUD modeling techniques have been developed to consider unpredictable core conditions such as AOA due to CRUD. In addition, during the core follow calculation, backend SNF analysis can be performed based on the same method with cross-section interpolation using the history information. The practical operational history of Korean PWRs is utilized for the V&V of RAST-K v2. According to the V&V, RAST-K v2 results are in good agreement with measured data and commercial nuclear design codes with respect to CBC, ASI, and fuel assembly power, with a maximum difference of 70 ppm, 0.017, and 5%, respectively. Recently, high-fidelity core depletion simulations have been achieved via multiphysics coupling with CTF and FRAPCON. Further, RAST-K v2 can be coupled with VIPRE

Conclusions
This paper presents the calculation models, engineering features, and application status of the newly developed nodal diffusion code RAST-K v2. State-of-the-art neutronics, TH, and depletion solvers are adopted in this code. A non-linear nodal method based on UNM is used with multi-group cross-sections generated by STREAM for a specific case matrix. The single-phase energy-mass conservation equation is solved for a closed channel for TH feedback, and the micro-depletion calculation is performed using the CRAM solver for fuel cycle analysis. To support users in conveniently designing and analyzing reactor cores, several miscellaneous engineering features are implemented in RAST-K v2. In particular, burnup adaptation and simple CRUD modeling techniques have been developed to consider unpredictable core conditions such as AOA due to CRUD. In addition, during the core follow calculation, backend SNF analysis can be performed based on the same method with cross-section interpolation using the history information. The practical operational history of Korean PWRs is utilized for the V&V of RAST-K v2. According to the V&V, RAST-K v2 results are in good agreement with measured data and commercial nuclear design codes with respect to CBC, ASI, and fuel assembly power, with a maximum difference of 70 ppm, 0.017, and 5%, respectively. Recently, high-fidelity core depletion simulations have been achieved via multiphysics coupling with CTF and FRAPCON. Further, RAST-K v2 can be coupled with VIPRE and BOA to model the CRUD deposition. Notably, RAST-K v2 has been validated as it has unique functions such as CRUD modeling. Thus, it is utilized in the area of machine learning for determining optimum loading patterns and developing core diagnostics systems.