Use of the simplified PN equations and transport corrections in the WIMS/PANTHER Embedded Supercell Method

The WIMS/PANTHER Embedded Supercell Method (ESM) provides a significant improvement in prediction accuracy in radial power distributions for pressurised water reactors compared to the standard “two-step” approach, without the need for a significant increase in computational resource. A companion paper presents validation of the ESM as previously presented, using PANTHER pin-by-pin diffusion to correct interface errors arising from the standard two-step approach. However, in principle any reference method can be used to solve the embedded supercells and correct the basic solution – this is a significant advantage of the ESM. A paper presented at PHYSOR 2016 demonstrated that use of diffusion theory introduces significant error relative to transport theory only in the high energy range, due to the discontinuous fission source when modelling the interface between two types of pincell. This paper investigates further improvement of the PANTHER solution through use of simplified PN in the fast energy groups. This solution method can be implemented without the need for any significant change to the calculation route and further improves agreement with a transport reference for a small computational cost. The remaining error in the solution is examined by derivation of an effective diffusion coefficient from a heterogeneous transport reference. This suggests a simple characterisation of the discrepancy, which can be easily corrected. Results are presented in comparison to WIMS for supercells and the KAIST benchmark. For the KAIST small core benchmark assembly power errors relative to a WIMS fine group transport reference are shown to be less than 0.5 %. Results are also presented in comparison to Monte Carlo for the Watts Bar benchmark.


INTRODUCTION
The ESM is a recent development in the WIMS [1] and PANTHER [2] codes which uses detailed local calculations to correct the basic "two-step" solution. The corrections are calculated by solving 2D supercells containing four quarter assemblies. The detailed local method used to date in the ESM is PANTHER pin-by-pin diffusion with homogenised pincells. The method is aimed particularly at capturing the severe flux gradients across UOX/MOX interfaces in PWR cores, but it is used at all interfaces including at the core edge. A companion paper at this conference ( [3]) presents validation of the ESM against plant measurements, for a range of operating PWRs, and a previous paper ( [4]) gives more detail of the methodology.
A paper presented at PHYSOR 2016 [5] demonstrated that use of diffusion theory introduces significant error relative to transport only in the high energy range, due to the discontinuous fission source when modelling the interface between two types of pincell. This paper presents a method for correction of this high energy transport effect by solution of the simplified spherical harmonics (SPN) equations, following previous work which developed a nodal SPN multi-group transport model for use in PANTHER [6].
Solutions are first investigated for a simple problem in order to demonstrate the fundamental effects. The method is then applied to solve more realistic problems as part of the ESM.

METHOD OF CORRECTION USING SPN
The SPN equations are an approximation to the full neutron transport equation, and have been shown to be successful in realistic problems, capturing effects beyond diffusion theory while being significantly less costly to solve than full transport. SPN solvers have been implemented in a number of nodal codes [7][8][9][10].
The approach taken here is not to perform a full multi-group SPN solution, but to solve the SPN equations for a single energy group in isolation. This is done for each of the fast groups only (where diffusion has been observed to be in error relative to transport [5]). These solutions are then used to derive transport corrections for use in the standard multi-group nodal diffusion solution. This method is appealing because it simplifies significantly the SPN solution and allows the continued use of existing nodal diffusion solvers.
The 1D even-parity transport equation can be solved with an angular expansion using normalised Legendre polynomials, resulting in a system of equations coupled by a tri-diagonal matrix: where A is the tri-diagonal geometric matrix linking Legendre leakage terms, ψ + is the column matrix of even-order Legendre polynomials, and q is the column matrix of source terms, with only q 0 non-zero.
Finite difference is used to solve these equations for each of the fast energy groups that contain the fission source -groups 1 and 2 in the six group structure used here. Using a pin-wise mesh, finite difference gives a sufficiently accurate solution for these fast groups.
The zeroth order solution of equation 1, for the scalar flux, is then influenced directly only by the second order leakage term (because A is tri-diagonal). This zeroth order solution is simply the solution to the standard diffusion equation, modified by the contribution from higher order terms which capture transport effects.
Rearrangement of the zeroth term of equation 1 allows the diffusion equation for the scalar flux to be restated with a source term now modified by this second order leakage: The contribution from the second order leakage can therefore be used as a correction to the source and the diffusion scalar flux solved for as usual.
The fission source and cross sections are taken from an existing diffusion solution, and mesh-and group-dependent corrections are calculated for inclusion in the multi-group nodal diffusion solver. A high order SPN expansion is used here, truncating after N = 21, since the computational penalty associated with increasing the order is minor.

APPLICATION TO A LINE OF PINCELLS
A simple line of pincells is solved to test SPN in correcting a smeared pin-by-pin diffusion solution, as used in the ESM. Previous work studied a similar problem in order to isolate transport effects at an interface between UOX and MOX pins [5].
In this work the two pincell types are both UOX and differ only in their enrichment. This simplification allows us to focus on correction of the transport effect due to the fission source discontinuity at the interface. The line is 40 pins long, with 20 of each pin type and reflective boundary conditions.
The line is solved using WIMS-CACTUS (method of characteristics) with heterogeneous pincells and this reference is compared to PANTHER smeared pin-by-pin diffusion with and without SPN corrections. The use of six energy groups was previously shown to be sufficient in the ESM [5], and the same six group structure is used here in both codes for consistency, in order to concentrate on transport effects. Indeed the same six group pin data is prepared for both WIMS and PANTHER solutions, before homogenisation is performed to generate data for PANTHER. Pincell discontinuity factors are derived for each pin type using a flux based on the in/out partial currents at the pincell edge, as described in [5]. Figure 1 shows the pin power and reaction rate differences by pincell, for the basic PANTHER diffusion solution and that with SPN corrections, both in comparison to the heterogeneous CACTUS reference solution. Both thermal (U235 absorption) and fast (U238 fission) reactions are shown to illustrate the effect.
The uncorrected PANTHER diffusion solution is unable to capture the flux shape in the fast groups due to the step change in the fission source. Application of SPN is successful in correcting this transport effect. The corrections are applied in the fast groups only, but affect the full multi-group solution, as can be seen by the changes in U235 absorption and in pin powers (which are largely determined by the thermal flux).
The result is pin power agreement with CACTUS to within a fraction of a percent. Thus for this problem, we can conclude that the only significant inter-pin transport effect arises from the fission source in the fast groups and is correctable using SPN.

COMPARISON OF EFFECTIVE DIFFUSION COEFFICIENTS
Various recent work has emphasised the importance of having an accurate method for the calculation of the diffusion coefficient used in a whole core nodal solution (see e.g. [11,12]). This section addresses the same issue, albeit empirically.
The difference between diffusion and transport theory can be examined by comparison of effective diffusion coefficients derived from the pincell average fluxes and interface current: where δx is the pincell width. This can be done in PANTHER and CACTUS, and with a smeared diffusion solution simply returns the diffusion coefficient. Figure 2 shows that there is a dramatic reduction in fast D eff at the interface between pincell types in the line of pincells described above. This feature is not present in the diffusion solution but is well captured by the SPN solution.
It can also be seen that there are small differences between the three solutions away from the interface. Firstly, the difference between diffusion theory and SPN is not only local to the interface where there is a discontinuity in the fission source. Secondly, the SPN solution does not perfectly reproduce D eff from transport.  boundaries away from the interface. A highly converged Monte Carlo solution using Serpent ( [13]) is also used to provide an independent reference (but using the same underlying nuclear data, from JEF-2.2, as is used in all calculations for this paper). Comparison of the SPN diffusion coefficient to that from a transport solution shows a small underprediction. This discrepancy is believed to be due to pincell homogenisation, but a simple group, material and condition-independent scaling is successful as a correction. E.g. for comparison to a 22g WIMS reference a scaling of 1.5% would be derived from the fifth row of table I.

IMPLEMENTATION IN THE EMBEDDED ROUTE
Inspection of the SPN corrections in a line of pincells shows a major effect at the interface, as expected, and also a much smaller, but steady, long-range effect, which shows that the diffusion solution is deficient even in a region where the flux behaves as a buckling-only solution. The embedded supercell methodology will capture the effect at the interface, but will miss the longrange effect because the supercell only extends over a half assembly on either side of the interface.
Wherever there is buckling present, the linkage to higher-order Legendre terms is active, via the tridiagonal coupling matrix, A. The second order contribution to the scalar equation can be written as a correction to the leakage term, or in other words a correction to the diffusion coefficient: The ratio ψ 2 /ψ 0 can be determined from the equation for the second order angular flux, assuming a buckling-only solution. This correction can therefore be applied in the whole core solution, with the correction based solely on the local buckling and no knowledge of the angular flux moments. As part of the SPN solution in ESM, the overall D eff correction in table I is also applied.
In 1D SPN is equivalent to PN, but this is no longer the case in 2D, therefore it is also natural to ask whether the method is as successful for a simple 2D problem. Figure 3 shows a comparison of pin powers for PANTHER diffusion ESM and SPN-corrected solutions with WIMS-CACTUS, for a simple 2x2 assembly chequerboard supercell with varying enrichment. The reference solution has in this case been generated using a 22-group CACTUS calculation, as will be used in the next section. The PANTHER pin solution is also now using pincell-average data condensed from 22 group lattice data, as is standard in the ESM [5]. I.e. no DFs are applied at the pin level to account for within-pincell heterogeneity. All these changes have little effect on the level of agreement observed, which is essentially preserved from the pincell line studied above.

APPLICATION TO THE KAIST BENCHMARK
The KAIST 1A benchmark ( [14]) is a 2D small core with UOX and MOX fuel, burnable poisons and an explicit baffle-water reflector, which has been used extensively to test the ESM. The core layout is shown in figure 4. Figure 4: The KAIST 1A benchmark core layout. Table II shows results for this benchmark, for the basic 'single assembly' PANTHER diffusion solution (SA), the ESM using pin-by-pin diffusion (DFN) and the ESM using SPN corrections, in comparison to a 22 group WIMS-CACTUS reference solution. In reference [5] the emphasis was on the improvement in pin powers that ESM offered over SA, this being a very challenging core with MOX. In this paper the emphasis is more on the global power shape, which is why table II shows errors for each assembly. SPN shows a significant improvement in global power shape over DFN, especially in assembly powers but also in detector reaction rates.

APPLICATION TO THE WATTS BAR BENCHMARK
Finally, the SPN method is applied to generate preliminary solutions to the Watts Bar 1 benchmark, a realistic 3D whole core [15]. The beginning of cycle 1 hot zero power state is modelled, with the control bank approximately 1/4 inserted.
A preliminary reference solution is generated using a full 3D Serpent calculation, converged such that reported pin power relative errors are less than 1 %. This is a detailed reference, modelled as defined by the benchmark specification, and including aspects not modelled by PANTHER, such as the barrel and neutron shield pads.   Table III: k eff , absolute pin power, relative assembly power and relative detector response percentage differences versus a Serpent 3D reference for Watts Bar 1 BOC1 HZP, and relative CPU times, for the three PANTHER models described in the text.

CONCLUSIONS
A method of deriving SPN transport corrections has been implemented in the WIMS-PANTHER ESM, and it has been demonstrated that this allows close reproduction of a transport solution for a line of pincells with very small remaining discrepancies.
The SPN equations are solved and the resulting corrections applied in the fast energy groups only, greatly simplifying the method while still successfully capturing the relevant transport effects.
It has also been shown that the effective diffusion coefficient is valuable for understanding of transport effects; to achieve high accuracy using diffusion theory it is not enough to replicate the local reactivity effects. An effective diffusion coefficient can be derived from any transport reference, allowing a direct assessment of the accuracy of the diffusion coefficient used in a whole core diffusion solution.
The SPN method has then been applied to more realistic problems and shown to further improve agreement with a transport reference, for only a small increase in computational cost relative to the ESM.