A Comparison of Advanced Boiling Water Reactor Simulations between Serpent/CTF and Polaris/DYN3D: Steady State Operational Characteristics and Burnup Evolution

Abstract

High fidelity modelling for nuclear power plant analysis is becoming more common due to advances in modelling software and the availability of high-performance computers. However, to design, develop and regulate new light water nuclear reactors there are, up until now, limited requirements for high fidelity methods due to the already well-established computational methods already being widely accepted. This article explores the additional detail which can be obtained when using high fidelity methods through Monte Carlo/Sub-channel analysis compared to industrial methods of cross-section/nodal analysis using the Advanced Boiling Water Reactor as a case study. This case study was chosen due to the challenges in modelling two phase flow and the high levels of heterogeneity within the fuel assembly design. The article investigates how to implement such an approach, from a bottom-up procedure, by analysing each stage of the modelling process.

1. Introduction

Modelling and simulation in the topic of nuclear reactor cores has made very strong progress over the last few decades, mainly sponsored through large international programmes like NURESIM [1] or through national efforts like Consortium for Advanced Simulation of Light Water Reactors (CASL) [2]. However, not all of the new efforts seem to be suitable for near term industrial application; thus, the UK has recently delivered their own industrial nuclear reactor core simulation [3] approach. In addition, all the new approaches have been targeted towards pressurized water reactor (PWR) development as these are the current leading reactor technology worldwide. We intend in this publication to have a deeper look into the application of high-fidelity methods for boiling water reactor (BWR) core modelling to determine the requirements which will be used to help improve future industrial BWR simulation technology.

The relevant key safety parameters and their evaluation in industrial application in nuclear reactors have not significantly changed over the past thirty years, despite significant advancements in computation power within this time frame. The main reason for this is that new high-fidelity methods, such as Monte Carlo analysis are yet to efficiently show greater depths of understanding when compared to the traditional methods which are significantly less computationally expensive [4]. Presently, there is no motivation to change the current methodology from an industrial perspective as regulatory bodies are satisfied with the detail being provided because the designs being commissioned are often well understood.

Today’s industrial methods determine few-group neutronics cross section sets for the core simulator via 2D, multi-group transport solvers in the lattice codes [5,6,7] for each fuel assembly, or several levels for axially varying fuel assemblies. These cross sections are then provided to a core simulator software package such as a 2D-1D nodal diffusion solver [8,9,10] which is coupled with an internal heat transfer module consisting of thermal fluid dynamics and a fuel rod model. These well-established solvers have the advantage of being extremely time efficient and have low computational requirements, which allows for full core burnup analysis, fuel re-shuffling and transient analyses to be performed with reasonable accuracy in acceptable time. The current industrial standards of different countries are listed in

Table 1. A description of potential cross-section generation and system tool codes used for several countries.

Country	Cross-Section Generation	Core Simulator
UK	WIMS [11]	Panther [12]
France	APOLLO-2 [13]	MANTA/SMART [14] CRONOS-2 [15]
Sweden/USA	CASMO	SIMULATE [16]
USA	Scale-Polaris [17]	PARCS/RELAP [18]

for light water reactors (LWRs); however, the methodology within these simulation packages very similar.

By obtaining the power and temperature distributions, these results can be later coupled to external system codes and fuel performance codes to provide sufficient regulatory information for the licensing process. The benefit of these methods is that they are widely accepted as they have been developed alongside years of experimental data. One of the main drawbacks of the current safety evaluation methodology is that nuclear regulators place a high emphasis on the most challenged pin under steady state and transient conditions. Individual pins are not specifically modelled within the current industrial standard methodology as these methods use a representative pin per node model. Pi power reconstruction is implemented offline which requires large safety margins due to the estimations of pin powers without appropriate thermal hydraulics feedback [3]. This conservative offline approach is widely used to ensure the safety of the nuclear installations.

In recent years, the USA has created the largest nuclear modelling and simulation project, CASL, which was provided with $242M by the US government over ten years to advance the country’s nuclear reactors modelling capability [19,20]. In addition to this, the USA has access to the most powerful super computers in the world [21] which have been utilised to produce high-fidelity models of their current reactor designs [22]. One of the drawbacks of these methods is the massive computational power and time required, although the consortium has made a significant effort to reduce these requirements by simplifying the process through the use of MPACT/CTF coupling, where MPACT is a 2D-1D method of characteristics neutron transport solver and CTF is the Nuclear Regulatory Commission’s subchannel code [23]. The high computational requirements and lack of access to the CASL software rules out most of these procedures for current industrial application, as well as for smaller nations. The complexity of the modelling procedure within CASL places only very few people in the world to be able to use these toolkits.

High fidelity modelling and its application is becoming increasingly popular within the research community [1], with Monte Carlo and Computational Fluid Dynamics (CFD) or subchannel models recently gaining momentum. It is often argued that these methods are required for non-conventional nuclear systems, due to the complexity of some designs and the absence of specially designed tools, which has also contributed to this area of research [24]. As computational power becomes more easily accessible, the modelling industry is likely to incorporate these methods more. Recent developments have seen a short transient being performed in 3D by coupling Serpent [25], a Monte Carlo neutron transport solver and SUBCHANFLOW [26], a subchannel code this project produced a pin by pin power profile of a small core during a reactivity control insertion [27]. The computational power required was 1000 cores, which is currently seen as inaccessible to most industrial institutions; however, two decades ago, this achievement would have been thought to be impossible.

The UK governments Department for Business, Energy and Industrial Strategy (BEIS) has supported the development of the UK nuclear industry under the £180M multidisciplinary Nuclear Innovation Program. This program aims to support the UK’s civil nuclear industry and reduce the overall cost of nuclear power plants by 30% which will support the decarbonisation of the country [28]. A part of this work aims to increase the UK’s skillset and ability to provide high fidelity modelling techniques and it is, therefore, important to understand the additional benefits that using high fidelity modelling provides over the conventional industrial standard methods.

The main focal point for high fidelity modelling has been on PWRs, which are the most common deployed reactor technology globally. There are several benefits with these designs as their geometry is relatively simple, with symmetrical fuel assemblies, an almost uniform water density and the fuel assemblies have uniform material properties. This has meant that the validation of the modelling procedure has been relatively straightforward. The second most deployed reactor type is the BWR, which are significantly more challenging to model, due to their two-phase flow and the complex fuel assembly design to promote a varying power profile. As most modelling and simulation toolkits have focused on the development of PWRs (for commercial application), the inclusion of BWRs has not received as much attention and simplified thermal-hydraulic models have been adapted to support two phase flow regimes. This poses the question: Would there be any benefit for the BWR community to embark on a high-fidelity modelling and simulation programme as we have seen for PWRs? Such a program could further the understanding of the transition of the power profile which is vital to understand when determining the most challenging pins for safety cases.

The overarching aim is to be able to understand where the current methods fall short of the high-fidelity methods, which will allow for an insight into where the industrial methods can be improved, without becoming too computationally expensive. This leads to the research question: Can we judge and observe the effect of increased complexity of BWR cores by using cutting edge coupled methods? This would allow us to identify the requirement for model improvements for standard coupled core simulators for BWR application. This work focuses on the Advanced Boiling Water Reactor (ABWR) as this is one of the most complex advanced reactor design types due to the presence of two-phase flow and axial heterogeneity. The ABWR has passed the UK’s Generic Design Assessment (GDA) [29]; however, GE-Hitachi have recently pulled out of committing to building the ABWRs due to concerns over the costs and financing of the project [30].

To model the ABWR, a standard industrial method is used, which takes a four-level approach to cross-section generation using Scale-Polaris (SP) [31] which are then applied to the core simulator package DYN3D [10]. The high-fidelity method couples the Monte Carlo transport solver Serpent [25] to the subchannel code Cobra-TF (CTF) [32]. The benefit of this method is that Serpent can provide an accurate pin power representation and CTF is then able to perform a thermal-hydraulic analysis on the fuel assemblies.

2. Model Description

As the ABWR is a commercial product the exact arrangements of their fuel assemblies as well as the loading pattern is not available within open literature so this work uses a patent filed by Global Nuclear Fuel-Americas in 2014 [33] which contains several fuel assembly types. Global Nuclear Fuel America is a GE-led joint venture with Hitachi limited [34] and, therefore, the patent does not provide proprietary information; however, the 10 × 10 assembly design provides a reasonable starting point. The patent provides seven assembly types; however, this study considers only two, thus reducing the modelling complexity and effort. Minor adjustments to the design filed within the patent have been made to provide a more realistic fuel assembly and the pin layout and pin compositions are shown in Figure 1 and Figure 2, respectively.

The initial water density within the assemblies is provided using benchmark data for boiling water reactor [35] and further data used can be found within

Table A1. Thermal hydraulic data used within the model.

Thermal Hydraulic Data	Value	Unit
Pressure at the outlet	7.171	MPa
Power per assembly	4502	KW
Mass flow rate per assembly	16	kg·s−1
Inlet temperature	278.3	°C
Gas gap conductance	2942	W·m−2·K−1

,

Table A2. Serpent’s initial coolant channel data.

Hight (cm)	Assembly Model Density (g·cm−3)	Initial Temperature (K)
0	0.7552	600
15.24	0.7552	600
30.48	0.7500	600
45.72	0.7309	600
60.96	0.6956	600
76.2	0.6433	600
91.44	0.5903	600
106.68	0.5395	600
121.92	0.4993	600
137.16	0.4598	600
152.4	0.4194	600
167.64	0.3958	600
182.88	0.3670	600
198.12	0.3357	600
213.36	0.3225	600
228.6	0.3002	600
243.84	0.2764	600
259.08	0.2717	600
274.32	0.2506	600
289.56	0.2398	600
304.8	0.2380	600
320.04	0.2159	600
335.28	0.2164	600
350.52	0.2095	600
365.76	0.1961	600
381	0.2020	600

and

Table A3. Material information used.

Part		Description		Value (cm)
Fuel pin		Fuel pin diameter		0.876
		Gas gap outside diameter		0.894
		Fuel Cladding outside diameter		1.026
		Fuel pin pitch		1.295
Water rod		Water rod water diameter		2.322
		Water rod cladding outside diameter		2.522
Assembly cladding		Inside cuboidal dimensions		13.42
		Outside cuboidal dimensions		13.82
Water surrounding the assembly		Cuboidal dimensions		15.539
Steel/water mixture below the assembly		Height		50
Outside of the reflector		Diameter		710
Part	Density (g·cm−3)	Temperature (K) If Not Fixed, Starting Temperature Provided	Fixed Density	Fixed Temperature
Water rod	0.7552	600	Yes	Yes
Water rod cladding	6.57	600	Yes	Yes
Outside cladding	6.57	600	Yes	Yes
Outside water	6.57	600	Yes	Yes
Natural uranium ends	10.424	900	Yes	Yes
Fuel Pellets	10.424	900	Yes	No
Fuel cladding	6.57	600	Yes	No
Subchannel coolant	Varies axially	Varies axially	No	No

within the Appendix A.

The full core model is achieved by arranging the assemblies with four different burnup stages (0, 5, 10 and 30 GWd/TU) which are slightly adopted from the ABWRs equilibrium core arrangement to provide a critical core configuration. The full core arrangement and the cycle number is depicted in Figure 3, which is taken from the ABWRs GDA [36]. Similarly, to the ABWR design, this analysis uses an assembly supercell arrangement, where the bottom right-hand corner of Figure 1 is always at the center of each four-assembly group.

3. Codes and Toolkits Used

3.1. Scale Polaris

SCALE 6.2 is a modelling and simulation tool developed by Oak Ridge Nuclear Laboratory (ORNL, USA) for criticality analysis, reactor and lattice physics calculations, sensitivity, and uncertainty analysis. Polaris has recently introduced the SCALE 6.2 module into the code for lattice physics calculations. Polaris performs transport calculations using Method of Characteristics (MOC). Self-shielding cross-sections are obtained with Embedded Self-Shielding Method (ESSM) which uses Bondarenko interpolation approach [31]. Polaris can perform lattice calculations using ENDF/B-VII.1 [37] nuclear data library with 252- or 56-energy groups.

3.2. DYN3D

DYN3D is three-dimensional code designed for performing coupled neutronic/thermal hydraulic steady state and transient analysis of LWRs. Helmholtz-Zentrum Dresden-Rossendorf (HZDR) have developed the code for more than 20 years [10]. DYN3D was verified and validated against various benchmarks, mostly for Russian VVER-type reactors and is a NURESIM reference code. BWR analysis was performed for the BWR Turbine Trip Benchmark developed by OECD/NEA [38]. DYN3D is widely used in academia for research purposes and also provides 16 licenses to industrial stakeholders in 7 countries [10].

The nodal expansion method, implemented within DYN3D for neutron kinetics, solves the neutron diffusion equations for Cartesian and hexagonal geometry in two- or multi-group approximation. Balanced equations for mass, energy and momentum are being solved using code’s internal thermal-hydraulics model for one- and two-phase flow in the reactor core linked with a fuel rod model for heat transfer.

3.3. Serpent

Serpent is a continuous energy Monte Carlo neutron transport code developed at the VTT Technical Research Centre of Finland. Serpent provides a pin-by-pin power profile across the core. To provide an accurate comparison with SP, a beta version of the ENDF/VII.1 data library was provided by the Serpent development team as Polaris is unable to use the ENDF/VII.0 data library which has been validated with Serpent. The Serpent version used within this article was version 2.1.30 and all burnup calculations were implemented using the Transmutation Trajectory Analysis (TTA) burnup method.

3.4. CTF

CTF is a sub-channel code produced by North Carolina State University (and now ORNL) within the framework of the CASL project [39]. The CTF project aims to advance the 1980’s Cobra-TF to a higher standard by providing significantly more thermal-hydraulic data required for regulation than its predecessor. The CTF version used within this work is version 3.7.0.

3.5. Serpent and CTF Coupling Procedure

Serpent provides an axial pin power profile across the geometry based upon twenty-five axial dependent positions along the length of the assembly and one for each pin in radial direction. This allows for a full 3D power profile to be provided along the fuel pins to CTF, which uses a total amount of nodal heights of 250, to allow for an accurate representation of the flow regime which must capture the development of the two-phase flow across the geometry. CTF does not directly output the coolant temperature yet provides the nodal enthalpy and pressure which are then input into a Python module [40,41] produced by The International Association for the Properties of Water and Steam (IAPWS) to provide the coolant temperature. CTF’s 250 nodes are then averaged across 30 axial heights to return the temperatures and densities back to Serpent. Within the full core analysis, the crossflow between assemblies is ignored, which is acceptable within ABWRs due to their wrapped assembly structure. A brief description of the coupling procedure is highlighted within Figure 4 to demonstrate in interplay between both codes and the developed script.

4. Methodology

The first investigation explores the preparation of cross sections and, thus, the neutron transport solver and burnup modules of Serpent and Scale/POLARIS (SP). To achieve this, four 2D cut outs of the assemblies are taken at intervals which would be used to define the axial variations within the nodal method as highlighted in pink in Figure 5.

Figure 5 highlights that the fuel assembly levels are placed at standard positions where the materials compositions change within the assemblies. One small discrepancy between SP and Serpent is the natural uranium cap on F10, which is modelled as purely fuel within SP whereas Serpent explicitly models the natural uranium.

The cross-sections from SP are then exported to DYN3D and the cross-section preparation as well as the correct transfer are verified via a comparison of DYN3D and SP. The next stage aims to understand the differences between Serpent and SP by investigating the 2D cross-section preparation under reflective boundary conditions for each of the fuel levels which are depleted in stages up until 60 GWd/TU to represent the maximum burnup experienced during operational conditions. This will highlight any differences in the burnup procedures within Serpent and SP and determine if any errors are being introduced by using Serpent’s beta ENDF/VII.1 data library and any discrepancies when the cross-section set is used within DYN3D.

The second investigation evaluates a fuel assembly, and this takes part in two stages, firstly purely neutronic and secondly a coupled neutronic and thermal hydraulic procedure. Within this study reflective boundary conditions are used on the sides of the assembly with vacuum boundaries at the top and bottom to provide the correct power profile. The purely neutronic study explores how high-fidelity methods vary from the industrial standard in the areas shown in

Table 2. Single assembly progressive test comparisons performed.

Model Description	Water Density Distribution	Aims
Uniform fuel	Uniform	To understand any differences in a standard model, like Pressurised Water Reactors (PWR)s, first validation of the modelling in both systems
Uniform fuel	Non-uniform	To understand how the water density averaging in levels in the industrial method affects the modelling approach
Four-level fuel	Uniform	To understand how the high flux gradients of the fuel behaves between the models
Four-level fuel	Non-uniform	To understand how the high flux gradients of the fuel is represented between the models and the water density combined

.

The coupling procedure between neutronics and thermal hydraulics is automatic in the nodal diffusion solver; however, this will take several iterations between Serpent and CTF to determine the final power profile and the related feedbacks which will be measured by the effective multiplication factors (k_eff) variation between iterations. By using either method, the correct operational water density distribution in axial direction will be obtained, which is important to provide vital safety parameters within the design.

The third investigation aims to model the full ABWR core under steady state operating conditions under normal operation. Within this case all boundary conditions are vacuum at the outer surfaces and Polaris’s assembly discontinuity factors (ADFs) [42] feature was used to provide ADFs at the boundary between fuel assemblies and between the exterior fuel assemblies and the reflectors. ADFs aim to provide a more realistic representation of the heterogeneous flux distribution at the boundary between nodes, as this is unaccounted for within the homogenized nodal diffusion approximation and can lead to errors at high flux gradient boundaries, in the thermal neutron flux, especially. One of the challenges within this study arises from obtaining the materials at each cycle for Serpent; this requires a highly detailed fuel assembly model which accounts for spatial self-shielding and each material to be modelled separately. An evaluation of the burnup curves over the fuel cycle is performed to provide an understanding of the effect of the different levels of detail the water density in combination with the material levels has on the burnup curve. Following this, the full core is modelled using both solutions and the cores power distribution is compared on a 25-node level to determine the power variations.

5. Results and Discussion

5.1. Single Assembly 2D Levels Results

This investigation aimed to understand the differences within the softwares solvers, and the master cross-section used by taking 2D cross-sections for each level of each fuel assembly, these were then burnt up to 60 GWd/TU using a constant fuel temperature of 900 K. This will provide an understanding of the way that each code handles burnup sequences and comparing them to one another. The results are presented for infinite multiplication factor (k_inf) in Figure 6 and Figure 7 with results for DYN3D added to test the few group cross section preparation and transfer procedure for SP to DYN3D as well as the burnup module.

The two assemblies were very similar, with only two fuel pin types (X6 and X9) materials being changed; however, it impacts on k_inf on all levels leading to different results. Each pin within the assembly plays a role in controlling the power distribution over time and the bottom level is heavily poisoned using Gd₂O₃ to reduce the initial reactivity which is caused by the higher water density present. In general, the results coincide very well for all 8 shown levels keeping in mind the different master libraries used for the calculations. The initial test is to determine if the cross-sections within DYN3D matched those which were produced in SP and this was confirmed within the figures as these values were identical between the codes. The maximum discrepancies between SP and Serpent occurs at the peak reactivity positions when the Gd₂O₃ is depleted and the speed at which this happens is dependent on the method used. There is a trend that Serpent reaches a peak reactivity at an earlier time to SP, with the maximum discrepancy being 671 pcm, in level three of assembly two which is caused by the large amount of Gd₂O₃ within this level. The bottom level of the model plays the highest role in the behaviour of the fuel assembly and will determine how quickly the power profile moves axially. Within this case the maximum variation is 315 pcm in assembly one, which occurs after at 12 GWD/TU, which is prior to the peak reactivity. The peak reactivity at 18 GWD/TU is in a very good agreement between the codes with a variation of 38 pcm. From Figure 6 and Figure 7 it is clear that Serpent and SP are in good agreement with one another and that the beta data library in Serpent seems to be compatible with SP’s validated library.

5.2. Single Assembly 3D Neutronic Evaluation

The 3D single assembly evaluation aims to progressively increase the complexity of the design to determine at which points the assumptions within the software methods start to affect the results. The initial study assumes a uniform fuel composition based on level one of fuel assembly two for the whole height of the fuel, with a uniform water density but real boundary conditions; this represents the simplest case study with the axial power profile depicted in Figure 8. The power profiles are in very good agreement with one another with this simple test.

The next stage introduces the water density variation across the axial height while the material is kept constant through the whole axial assembly; thus, the first case is created with real influence of a thermal hydraulics parameter but still without the coupling of neutronics and thermal hydraulics. Within the industrial method, the water density must be averaged across the heights represented by that cross-section, which deviates from the high-fidelity method which can incorporate each stage of the water density. It must be mentioned here that the water density in the core simulator is based on the water density of the four levels of the cross-section preparation while in later simulations the real water density is finally determined through the thermal hydraulics module of the core simulator.

Within Figure 9, it becomes obvious that the averaging process (based on four levels) plays a significant role within the peak power location as DYN3D’s power peak is at a different location to Serpent’s. DYN3D also shows that following a new water density input, there is a change in gradient of the power profile; this is due to the averaging process creating a virtual boundary between these levels, and the diffusion equation must determine the flux distribution between these levels and, hence, the power change. The results in the analysis are concerning as the assumptions for the purely neutronic model does seem to break down once the water distribution is added. It should be noted that as DYN3D is meant to be used as a coupled solver and a further evaluation of these results is required once DYN3D produces its own water density profile across the geometry as the boundaries seen within Figure 9 will be reduced in this case as the water density variation between nodes becomes smaller.

The criticality of the system (k_eff) provides an insight into the behaviour of the neutronic solver and is often used as a method of validation. Within this case two tests were performed using Serpent, one with and one without the thermal scattering libraries linked to the coolant material, these libraries increase the accuracy when bound atoms such as water [43] are used. In the case of the industrial method, the thermal scattering libraries are already applied within SP when the cross-sections are produced.

Table 3. Comparison of the reactivity coefficients of the single assembly with a uniform fuel material.

Model	DYN3D keff	Serpent keff	Serpent Error (±) 95% Confidence [pcm]	Variation [pcm]
Uniform water density no thermal scattering library	1.04872	1.04227	12	−590
Uniform water density thermal scattering library	1.04872	1.04953	11	74
Non-uniform water density	1.03798	1.02540	11	−1182

highlights the very prominent role of the thermal scattering libraries within the Monte Carlo simulations accuracy, as the discrepancy reduces from 590 pcm to 74 pcm. Due to the differences in the power profile, it is difficult to determine if this is an error reduction as the systems behave in very different ways; however, the inclusion of the scattering libraries does play a significant role within the systems’ behaviour. This highlights one drawback when using the high-fidelity method, as the current version of Serpent does not support the inclusion of thermal scattering libraries when using the interface cards feature which is required for the coupling procedure to assign the water density and temperature axially. Thus, this study has highlighted that the high-fidelity method will include an error when modelling BWR’s due to the modelling methodology not currently being capable of accurately predicting the criticality effect of the thermal scattering in water.

The next stage represents the real fuel assembly, with four axial material compositions (see Figure 5) with each using a new cross section set to account for the change in materials while the water density is kept uniform. Nodal diffusion methods have been shown to work accurately when passing through a uniform material although this method breaks down at high flux gradient boundaries which could occur between the fuel levels due to the inclusion of fixed burnable poison or between high/low moderated regions. Figure 10 investigates the severity of the flux gradients between the fuel levels. In this case, the codes are within very good agreement of one another, which implies that these flux gradients caused by the material compositions are on a magnitude that can be treated to an acceptable accuracy through diffusion approximation; thus, no major differences can be noticed between the solvers.

The final study aims to determine the combined effect of all the previous cases, with four levels and a non-uniform water distribution based on the levels, without direct neutronic/thermal hydraulic coupling.

Applying the full set (four levels with different materials and different water densities) leads to a similar power peak shift within Figure 11, as with the earlier study and this was determined to be caused by the non-uniform water density (Figure 9). An emphasis on the importance of the thermal scattering libraries is evident from

Table 4. Four leveled fuel assembly reactivity coefficients comparison.

Model	DYN3D keff	Serpent keff	Serpent Pcm Error (±) 95% Confidence	Pcm Variation
Four-level assembly uniform water density no thermal scattering library	1.03108	1.04924	11	1679
Four-level assembly uniform water density thermal scattering library	1.03108	1.03172	11	60
Four-level assembly non-uniform water density	1.02430	1.01703	17	−698

. However, the total discrepancies for the four leveled approach are slightly lower than with the uniform fuel case with different water densities. Overall, up until now, this study has already highlighted several key aspects. Initially the high-fidelity method is unable to use the thermal scattering libraries when a non-uniform water density is used which has been shown to create quite large deviations under a uniform water density. Secondly, due to the averaging process used for the cross-section preparation for DYN3D, the power profile is unable to match that of Serpent’s for a purely neutronic analysis with the peak power shifting using the current industrial standard method.

To reduce this shift in the peak, a serial analysis was performed to see if a refined averaging process will provide a closer match in the power profile. To achieve this, level one’s cross sections were split into section by using the two different water density averaging schemes shown in Figure 12 to reduce the strong step change in water density between level one and two.

The new power distributions in Figure 13 highlights that case 1 and 2 both provide a power profile which more accurately represents that of Serpent with a broadening of the power peak in DYN3D. It is, therefore, confirmed that the power profiles do experience a shift due to the averaging process used within the nodal analysis, which especially within the lower region of the core which could produce inaccuracy in identifying the peak power location. This has highlighted that a refinement study could be considered within nodal methods to achieve the correct location of the power peaks in the case the problem still appears in the fully coupled solution. This would be less important in PWRs due to the less severe axial variations in the water density.

The effective multiplication factor is also improved by splitting the bottom as depicted in

Table 5. K_eff of the three different assembly models.

	Serpent	DYN3D Original	DYN3D Case 1	DYN3D Case 2
keff	1.0170	1.0243	1.0205	1.0189
pcm variation from Serpent	-	698	332	184

. Further evaluations within this article use the original DYN3D model, to best represent the industrial method applied. This has also implied that the use of the thermal scattering libraries is likely to have minimal affect and that the main driver for the errors comes from the nodalisation process of the water density during the cross-section preparation process.

5.3. Single Assembly Coupled Neutronic and Thermal-Hydraulic Evaluation

The core simulator DYN3D comes with a complete coupled neutronics and thermal-hydraulics package which provides steady state and transient opportunities by using a coupling scheme which is determined to have converged once the users defined k_eff variation between iterations is achieved [44]. Similarly, the coupling process between Serpent and CTF focused on the variation in k_eff between iterations due to thermal feedback effects, where these results are presented in

Table 6. Serpent/CTF variations between iteration numbers.

Iteration Number	keff	Serpent Pcm Error (±) 95% Confidence	Pcm Variation from Previous Stage
1	1.01703	17	NA
2	0.989339	18.1	2752
3	0.987905	18.5	147
4	0.989163	17.7	−129
5	0.988656	17.4	52
6	0.988832	16.8	−18
7	0.988786	16.4	5

.

After seven iterations the variations of k_eff the model has converged within the statistical error of the Monte Carlo calculation. A further analysis of the pcm variation per step is provided within Figure 14, where the largest deviation occurs after the first stage as the water density and temperature deviated from the initial assumption. After stage four, the variations between stages becomes minimal, thus hinting that the system has converged. Earlier efforts highlighted that, despite a low variation in k_eff, the power profile can still vary significantly; to determine if this was the case, Figure 15 shows that, similarly to k_eff, the power variations are minimal after iteration 5, which represents step 4 in Figure 14. Following the convergence being met, the full radial power distribution is presented within Figure 16 where the pin powers remain within their original position and just small deviations in the power are visible as the model converges.

Figure 17 shows the axial power profile delivered by the coupled high-fidelity calculation is within good agreement with the DYN3D results once the simulations have converged. A similar pattern where Serpent’s power profile is slightly lower than DYN3D’s has been observed in the previous chapter; however, this deviation is reduced in the coupled results. This result provides confidence that the 2D-1D approach in DYN3D is able to accurately predict the axial power profile to a high degree of accuracy for this fuel assembly.

From a safety perspective, the fuel temperature is one of the key parameters for the design process together with the related cladding temperature, as these values provide the safety margin for the design of the fuel itself, but also for the core design. CTF provides a new in-built fuel performance module which provides the fuel surface temperature and centreline temperature which can be obtained for each pin via the pin powers provided by Serpent. DYN3D does not obtain any details within the node, since the approach is based on one representative fuel rod per node, so the fuel centreline temperature is approximated by assuming an even power distribution across all pins which provides the average fuel centreline temperature. To obtain the maximum fuel temperature of the assembly, DYN3D would use an offline pin power reconstruction procedure based on the pin-power map of the lattice calculation. The results can be refined by coupling to a fuel performance code, such as ENIGMA [45] or TRANSURANUS [46,47]; however, this is beyond the scope of this project, so only the average temperatures are considered here.

The slightly shifted power peak in Serpent leads to a similarly shifted fuel temperature peak as highlighted in Figure 18. The average fuel temperatures between the codes are in very good agreement with one another, showing that the industrial method can provide a good representation of these average results. DYN3D’s slightly elevated temperature would also provide a small safety margin within this case—thus, the result of the industrial code is conservative.

The density of the coolant is one of the key parameters, since within a BWR steam is formed inside the fuel assembly and later is used within the turbine. However, the heat transfer ability of steam is lower than water and the steam bubbles reduce the moderation of neutrons; thus, the coolant density curve is of high importance for any BWR simulation. The average density can be seen in Figure 19 where the horizontal lines in water density are the positions of the grid spacers. The grid spacers accumulate water due to the restricted flow and induce a slight pressure change; this effect is ignored in DYN3D, which cannot incorporate grid spacers. The grid spacers also create a pressure drop across the axial height, which would slightly reduce the mass flow rate which could account for the slightly higher coolant temperature. The average water density profiles are nearly identical between the two codes, which shows that the high-fidelity method does provide addition detail, while the average value does not differ significantly. This is an impressive outcome keeping in mind the fully modelling each channel and the inclusion of crossflow within the assembly and the inclusion of grid spacers in the high-fidelity solution.

The coolant temperature is also important to understand the phase change location of the water, which significantly reduces the heat transfer properties of the coolant. Within Figure 20, CTF shows an accelerated heating up of the water compared to DYN3D. This is caused by the slightly shifted power profile; however, this effect should not be as significant as the figure suggests which is partly based on the very much stretched x-axis spanning only over 25 K. It should be noted that DYN3D has a directional change in the water temperature at ~144 cm; this is the location of the transition between fuel level one and two. It seems as if DYN3D is unable to accurately determine the temperature across this nodal boundary; this is likely to be caused by the large reduction in power across these boundary levels which causes this instability. To try and reduce this effect, both DYN3D’s internal two-phase flow method and the ASME IFC-67 thermodynamic steam tables (IFC-67 were replaced by IAPWS-IF97 in 1997 [48]) were compared and showed very little difference with the results.

Overall, the codes are in good agreement with one another, with their only being a 4 K temperature difference in the coolant temperature, which has negligible effect on the coolant density compared to the void content.

The coupled neutronic and thermal-hydraulic performance of DYN3D on fuel assembly basis has shown to provide a high level of agreement within the critical safety parameters which have been evaluated here. Despite the purely neutronic version showing a deviation in the peak power, this has been resolved once both code suites have provided their converged fluid flow regimes. The largest deviation came from DYN3D’s coolant temperature, which showed quite a variation compared to CTF; however, the total deviation was only 4 K, which is unlikely to be significant for safety procedures. Thus, based on these results, we can claim a substantial basis of good agreement on fuel assembly level between both approaches.

5.4. Full Core Pure Neutronic Evaluation

The first part of this investigation was to obtain a detailed set of assembly materials required for the full core analysis; this procedure took place by depleting both fuel assemblies to the required burnup cycle stages using a uniform fuel temperature of 900K and the benchmarked water density. This task was a non-trivial one in Serpent as each fuel pin containing gadolinium had to be subdivided into ten radial regions and all pins were axially divided into thirty regions to account for changes in radial and axial depletion. The fuel assembly multiplication factor variation with the burn-up are presented in Figure 21 and Figure 22. The fuel assembly burnup curve of the modelled BWR fuel assembly is significantly different to the curves regularly seen for PWR assemblies with burnable poisons. Both codes deliver comparable trends showing a first peak about 5 GWd/TU, a dip between 12 and 27 GWd/TU and a second peak around 32 GWd/TU, but with clear differences in the height of the peaks and especially in the form of the dip. To deepen the understanding of the effect of the level structure and the constant fuel temperature in neutronics a fully coupled DYN3D burnup calculation is added in orange. First, there is the initial bias which is just confirming the differences seen in Table 6. Second, there is the strong dip in the DYN3D result around 12 to 13 GWd/TU, where the results significantly differ from the Serpent results. This seems to be an artefact of the level structure. Thirdly, there is the stretching of the whole burnup curve over time which could be an effect of the real fuel temperature and water density distribution which has typically an influence on the burnup behaviour. Finally, the proposal for a next, future step to apply a history correction [49] in the DYN3D calculation to eliminate or reduce the effects created by the assumptions of water density and fuel temperature in the lattice calculation.

From Figure 21 and Figure 22 and

Table 7. Pcm variation at each burnup stage for both fuel assemblies.

Assembly Number	A1				A2
Burnup (GWd/TU)	0	5	10	30	0	5	10	30
Serpent	0.9089	1.0827	1.0204	1.0098	1.0154	1.0933	1.0228	1.0265
DYN3D	0.9208	1.0839	1.0050	1.0238	1.0243	1.0955	1.0079	1.0425
Pcm variation	−1414	−97	1505	−1353	−860	−183	1441	−1496

, the fuel assemblies k_eff across the fuel lifecycle show large discrepancies between Serpent and DYN3D; however, the trends in the curves are similar. The initial pcm deviation at day zero was earlier demonstrated to be strongly influenced by the thermal scattering library in Serpent, which is currently not available in the required setup, another cause can be in the water density averaging across the axial levels of DYN3D. This could be further reduced by adding an additional level of cross-sections to reduce the errors incorporated by the averaging process as shown before. From a modelling perspective, the additional natural uranium region on the F10 fuel pins in Serpent could account for the slight initial reduction in k_eff due to the parasitic absorption within these regions and, over time, this could see ²³⁹Pu being bred, which could account for the higher k_eff later in the cycle, but the total effect of this is unlikely to be significant. To account for this, ²³⁹Pu in the industrial methods would require more cross-section productions between levels 2 and 3, and this level would only have a minor effect on the systems’ behaviour. Similarity to the second investigation, the original water density model was used within DYN3D, which has caused an overestimation in the bottom of level two, thus depleting the uranium faster over the full burnup curve; this is likely to be the driving factor for the differences between the two softwares. Overall, the accumulated differences within the model all contribute towards the significantly different results shown by DYN3D and Serpent, which could lead to an overestimation or underestimation of the reactivity coefficient of the equilibrium core.

Due to the lack of symmetry within the assembly, which is created by the water holes and the assembly map, a nonsymmetrical power distribution is created which Serpent can account for, while this not possible in the nodal approach of DYN3D due to the homogenisation process. To determine the severity of this effect, the pin power profile for assembly one is determined across the burnup cycle as shown in Figure 23. It should be noted that the supercells within the ABWR often contain cross shaped control rods between the assemblies, leading to asymmetric water gaps around each fuel assembly. If the control rods are inserted, this would cause additional power distribution imbalances, which would also be missed by nodal methods; these are likely to be more severe, yet more localised than the pin power distributions.

From Figure 23, the fresh fuel power profile is non-symmetrical from the design of the fuel assemblies and due to the water holes, there is only corner to corner symmetry whereas with a PWR there is usually full quarter symmetry. Through the burnup procedure, the power profile becomes highly non-symmetrical; this effect becomes quite severe at 30 GWd/TU with the power in the top left and bottom right showing 40% discrepancies. To account for the non-symmetrical power distribution, SP provides the option to create ADFs for each burnup and branching steps as well as for each assembly side separately. From Figure 23, it becomes obvious that a side ADF will produce an average boundary flux, where the true unbalance within the assembly is more corner based. Some studies have been performed to provide corner ADF’s for hexagonal assembly LWR’s [50]; however, this study has identified that a similar treatment is likely to be beneficial for modelling BWR’s due to the highly non-symmetrical nature of the fuel assemblies. ADF’s within SP is a simplification for assemblies as SP assumes that the neighbouring assembly is of the same type as that assembly being modelled; however, this is currently the only available approach to account for the non-symmetry of the problem inside the homogenization area of a full fuel assembly. ADFs are also used for assemblies adjacent to the reflector within the full core model, for which SP has its own method of producing ADFs. It can be expected that calculation of the ADFs using, for example, two-dimensional models (i.e., super cells) can potentially lead to further improvements. However, this approach would require much higher efforts for generation of the cross-section libraries due to the larger sizes of the super cells in comparison with a single fuel assembly.

Following on from Figure 11, Figure 24 shows the axial power profiles for the assembly one at different burnup stages. Within the second study the axial power profile of Serpent was slightly lower than that of DYN3D and this was deemed to be due to the assumptions of the water density in the lattice calculation which is forwarded into the nodal calculation. This work has been further investigated in Figure 24, as the depletion of the ²³⁵U will vary across each material zone axially, with this model using 30 axial depletion regions in Serpent and 26 in DYN3D. As Serpent’s power profile starts at a lower point, the materials within the lower section are depleted faster than in DYN3D and this explains the transitional discrepancies of the power profile between the burnup stages two and three within Figure 21. The different speeds at which the power profile moves between the stages is critical for a safety analysis as these determine where the hottest pin or pin section is at any point in time. This study has identified that the pure neutronic assumptions within DYN3D and Serpent cause large variances of the speeds at which the power profiles move in BWR’s due to the high levels of axial heterogeneity within the design which are not considered within current high-fidelity methods for PWRs. To fully investigate this, a coupled high-fidelity model is required for each assembly at each burnup step to determine if these match the coupled DYN3D burnup curve which is already included within Figure 21. This fully coupled calculation result with DYN3D is achieved using the water density as calculated through the thermal hydraulic module and the fuel temperature as calculated through the fuel rod model. These results help to deepen the understanding of the effect of the approximations which are introduced through the level structure of the lattice calculations. Overall, while the peak powers between the pure neutronic models do not differ significantly, their exact locations vary between the codes, while in all three cases the coupled solution shows reduced peaks which indicates that all results seem to be conservative. However, the identified differences between the codes does pose the question as to whether the used methods can accurately determine the time and location of the hottest pin to the required degree of accuracy through the life of a fuel assembly. This would be very helpful to justify to the quality of the safety analysis and to reduce over conservative safety margins. Nevertheless, the history of reactor operation with only very limited number of fuel failures confirms the approach used up to now as sufficiently accurate for reactor operation so long as the boundaries are not pushed to a too large extent.

The full core simulations of the ABWR were performed using Serpent’s MPI feature on the University of Liverpool’s High-Performance Computer cluster, Barkla. Due to the memory requirements the simulation used 80 Intel Xeon Gold Skylake processors. The simulation used 8M neutron populations for 1000 active cycles and 300 inactive cycles, after the fission source convergence was shown to be achieved with this configuration. The limitation of this simulation is like in other research institutes, where the queue times become significant due to other users and runtime limits for simulations are set to a maximum of three days.

The k_eff of the full core neutronic model in Serpent was 1.02394 (±8 pcm) which provides only a 48 pcm variation from DYN3D, despite the differences in the material compositions which shows that, in large systems, the differences can be easily overlooked due to error cancellation if only comparing a single variable. The previous steps have shown significant differences between the two softwares outputs, which would have been missed if only the k_eff of the full core was evaluated.

From Figure 25 it is obvious that the power profile is not completely symmetrical with Serpent, with the top right of the figure having a slightly higher power than the opposing corner. This is due to the non-symmetrical power distribution in each fuel assembly as highlighted in Figure 23.

Figure 26 shows the difference in each assembly power as the % of the total full core power, where the negative values represent a higher value in DYN3D compared to Serpent. Within this analysis DYN3D overestimates the highest power assemblies, which is a more conservative approach as these are the assemblies with the hottest pins. This can be explained due to DYN3D showing a larger value of k_eff for the central assemblies within the burnup stage, which implies the loading pattern provided this conservativeness and not the modelling procedure. The largest discrepancies occur near the exterior boundaries; this is within the vicinity of the assemblies where assembly burnup of 30 GWd/TU are positioned, and these have the second largest deviation in error from the assembly burnup models, and, since the evaluation is based on percentages of power, the error in low power areas is higher. It should be noted that the use of the reflector ADFs in DYN3D significantly reduced the discrepancies within the assemblies adjacent to reflector. Errors also occur within the centre of the core; this can be explained due to the largest variation burnup assembly on 10 GWd/TU being mainly positioned within this vicinity.

The radial power differences in the models can all be directly related to the input assembly compositions, a cascading effect of the modelling variances. Figure 26 emphasises the maximum deviation in the top right to 0.05% whereas the opposing assembly is only 0.043. Overall, the purely neutronic assembly power distribution is within reasonable agreement between the two codes and both codes identified the same fuel assembly as the limiting one with the highest power.

The purely neutronic normalised axial power profile is displayed within Figure 27 and the axial shift in power is like the single assembly analysis. There are now more factors within the production of the full core axial power profile compared to Figure 24 which has highlighted that the overall shift in power during burnup with each assembly with its burnup stage contributing towards the axial power profile. These effects of each fuel assembly are weighted to with the total power distribution per assembly shown in Figure 25 which has determined that the highest power assemblies are in burnup stages one and two. In this stage the fuel assembly simulations show that the lowest variations in their power profiles between the software. This allows Figure 27 to obtain the familiar shape and a good agreement between the codes; however, this could change if the core behaviour would be followed over the full cycle. Beyond the axial height of 150 cm, the power peaking shown within DYN3D is provided by the assemblies within burnup stage four, which within DYN3D will provide a higher flux at high axial points within the design. Despite this area not being critical from a safety perspective, the flux distribution is, therefore, slightly distorted due to the burnup discrepancies identified in Figure 24.

Overall, the purely neutronic code comparison for the full core is within good agreement between the two codes within the critical areas of interest within a snapshot in time at the core start up. In this instance the hottest power pins are within assemblies in burnup stages one and two, which have almost similar power profiles. This is very much dependent on the exact burnup stage as this investigation has highlighted that DYN3D’s axial power profile lags slightly behind Serpent’s over burnup and this gives evidence to the thought that the hottest power pins further into the burnup cycle are likely to shift, which would occur at different rates between the codes. The analysis has highlighted that the detailed behaviour of both modelling methods is significantly different between the models following the initial burnup stages, which has caused an overall cascading error through the results.

5.5. Full Core Neutronic and Thermal-Hydraulic Evaluation

The coupling procedures for the iteration steps for Serpent/CTF are shown in

Table 8. Full core K_eff iterations between Serpent and CTF.

Stage	Keff	Pcm Error 95% Confidence	Pcm Variation from Previous Stage
0	1.02394	8	NA
1	0.99218	8	−3126
2	0.99412	8	196
3	0.99587	8	177
4	0.99939	8	354
5	1.00765	8	820
6	1.00938	9	170
7	1.01227	8	283
8	1.01110	8	−114
9	1.01241	8	128
10	1.01125	9	−113

which shows that the coupling procedure does not converge in the same manner as the single assembly with a visual representation in Figure 28.

Up until stage four, the effective multiplication factor seems to be converging and then the results suddenly start to become unstable. The first area to determine why the system has not converged by investigating the radial power profile in Serpent, as shown in Figure 29 which shows that following the first iteration the build-up of unsymmetrical power profiles which shifts the power to the top right of the core layout after the second iteration. The effect of this is that the highest power assemblies become less favourable to fission within stage three because of the increase Doppler broadening and the lower water density. This affect then has the exactly opposite effect on the proceeding iteration, which identifies that this is the cause of the instability of the coupling procedure. This affect is thought to be more severe than the single assembly because the water and fuel temperatures variation is much wider in the full core system compared to the single assembly. One of the drawbacks of this, is that the hottest power assembly is also obtaining higher power values per each iteration which implies that the results will be inaccurate.

Stage ten of the coupling process had the same power distribution arrangement but not as high peak values as stage four within Figure 29. Due to convergence not being achieved, it would not be promising to go into a deeper comparison of any radial results; thus, only a short investigation will be performed on the integral values to see what has been achieved through convergence of k_eff as an integral value and to demonstrate what kind of additional information can be delivered using the fully resolved fluid dynamics with CTF.

The fully coupled axial power profile comparison in Figure 30 shows that Serpent’s higher axial profile trend continues, with Serpent having a slightly broader peak power than DYN3D, which can be explained by the 10 GWd/TU assembly power profiles power contribution.

The coupled results of the full core could not be used as a direct comparison due to the not achieved convergence on the coupled high-fidelity solver. However, this section has shown that the input material compositions of the equilibrium core varied significantly between the two codes which implies that even at the start of life there is a difference between DYN3D and Serpent to be expected. This could be investigated further by reducing the number of variables between the models to determine the main cause; however, this study has highlighted that the most likely cause is due to the different water density representation used within the different codes. Further analysis is required with regard to how to best apply high-fidelity methods with BWRs, where the lack of symmetry and high levels of heterogeneity have shown on a full core level that these could have significant implications for the long-term operation of the design and its fuel assemblies and following the potential for further optimisation of BWR fuel assembly design and the related operation.

5.6. Discussion of Possible Ways to Improve the Results

As it was demonstrated in the previous sections, the results of the currently applied industrial techniques can have significant discrepancies in comparison with the high-fidelity coupled neutronics-thermal hydraulics simulations. Undoubtedly, the results of the industrial codes used in present study can be further improved by applying more advanced techniques for ADF generation or, for example, application of the SPH factors (which were not used in the current study) or the combination of both mentioned techniques. However, all these techniques would require additional significant investments into the developing of the appropriate computational models, verification and validation and will not overcome the limit of nodal codes delivering only coarse mesh results instead of fuel pin related data. Therefore, this way can become computationally expensive and not so attractive for application by industry since it will still require the offline pin power reconstruction.

To overcome mentioned above limitations, another approach can be proposed. As it has been demonstrated in the present study, high fidelity coupled pin-by-pin neutronics and thermal hydraulic simulations can be very interesting, but unfortunately costly and difficult to implement and use in day-to-day industrial work. Therefore, another method was proposed [51,52] to overcome these challenges while keeping the computational time on the level acceptable for industrial application. In this technique, high fidelity methods would be applied only to the assemblies where high errors are expected (for example, assemblies adjusting to the reflector or control rods) or the highest loads will be needed to determine safety parameters. These assemblies could be resolved down to the pin level (using the boundary conditions from the full core nodal solution), while the fuel assemblies in the less challenging surrounding (i.e., central assemblies) would be still simulated on the nodal level. The idea of this approach is schematically shown in Figure 31.

Considering the multi physics nature of the reactor simulations, the application of the Monte Carlo codes as a neutronics solver is expected to still be time consuming even when applied on the level of a single assembly. Therefore, less general but faster neutron transport solvers should be used for this task (see, for example, [53,54]) and coupled to thermal hydraulics solvers. This approach, in our view, has strong potential to resolve the problems outlined in the present study without involvement of the ADF and SPH factors and other approximations aiming to reduce the naturally inherent error in diffusion approximation.

6. Conclusions

This article aimed to determine the variations between the industrial toolkit SCALE-POLARIS/DYN3D and a high-fidelity method of a coupled Serpent/CTF approach for the application in BWRs to determine where the industrial code struggles to find the details which were highlighted in the high-fidelity methods. This approach used a bottom-up procedure, where code to code comparisons were made from the 2D fuel assembly to, finally, a full core comparison. The depletion methods within both codes for 2D arrangement were in near perfect agreement between SCALE-POLARIS and Serpent, despite the different methods used within the software and the use of Serpent’s beta version of the ENDF/VII.1 data library.

Stage two investigated the performance of a single assembly with reflective boundary conditions by using the cross-sections from four axial levels from SP in DYN3D. From a purely neutronic perspective, the SCALE-POLARIS/DYN3D was not able to capture the axial power distribution in as much detail as required for the highly heterogeneous problem and this was identified to be due to the assumptions when splitting the axial geometry into the traditionally used four levels. Further improvements with the axial power profile shape and a reduction in K_eff variation were achieved by dividing the bottom level into more sections. A hot full power analysis was performed, and SCALE-POLARIS/DYN3D’s average fuel temperature and water density were in a good agreement with the high-fidelity methods and the axial power profile shape became closer to the high-fidelity methods. It was noted that SCALE-POLARIS/DYN3D did not produce the same water temperatures, yet matched the water density, which implies that the code used different water property data. Using the SCALE-POLARIS/DYN3D, a new fully coupled approach for the burnup case has been delivered for the investigated fuel assembly which significantly differs from both pure neutronic results.

Stage three investigated a full core supercell arrangement based upon the ABWR design. This stage highlighted that the burnup procedures of the individual assemblies provided large discrepancies between the two software tools and the axial power profile transition within Serpent varied significantly compared to DYN3D. These discrepancies are attributed to the long-term implications of the different power shapes between the methods, which will have an accelerated affect due to depletion. SCALE-POLARIS/DYN3D overestimated the power within the hottest fuel assemblies, which provided a slightly more conservative approach towards the modelling high fidelity approach; however, this was determined to be due to the fortunate arrangement of the core layout and this trend would not hold. Unfortunately, the third study could not provide a converged fully coupled Serpent/CTF model. Due to the large variations within the burnup procedure, this work has identified that the assemblies used within the equilibrium core layout could be significantly different when using industrial methods compared to high fidelity methods. At the time of writing, there is no solution to resolve the coupling errors of the full core analysis; however, it is assumed that these discrepancies of the full core model will be magnified and will have an impact on the safety parameter assumptions we make within BWR safety analysis.

Finally, a possible way to overcome the mentioned problems with the help of the novel multiscale and multi physics approach is proposed to support the future development of a more accurate industrial approach with acceptable computational demand.

7. Further Work

This work has identified that the materials used within an equilibrium core can be significantly different to those predicted by industrial methods. Due to limitations in the convergence of the high-fidelity method, the severity of these impacts could not be investigated in all detail; however, the power level profile is likely to be broadened from early indications. This could have a positive effect, due to the reduced peaking of the fuel axial power profile, but this must be established with further work. Thus, it would be ideal to apply the newly developed tools of the CASL project (even if this will require massive computational resources) to this kind of analysis to create reference solutions for the upgrading of industrial tools for BWR analysis.