© T. Okamura et al., Published by EDP Sciences, 2021

1 Introduction

The importance of nuclear fuel cycle simulation has grown in recent years as we have faced a transitional period of rebuilding a future that uses nuclear energy in response to drastic changes in the environment related to energy. This dynamic transition is occurring not only in Japan but throughout the world following the Paris Agreement on climate change in 2015. Although nuclear energy is thought to be a promising zero-carbon energy source, nuclear energy has problems from both technical and social perspectives, including nuclear waste, technology inheritance, and consensus building.

The nuclear fuel cycle promoted by developed countries is an extremely complex system that combines multiple processes including nuclear reactors, reprocessing, vitrification, storage, transportation, and disposal. This makes it necessary to precisely analyze nuclear energy scenarios by taking into account uncertainties in social trends and cross-sectoral perspectives in rational nuclear R&D and policy making. In addition, in order to achieve sustainable nuclear energy utilization and advanced nuclear energy systems, radioactive waste disposal, which is one of the biggest problems in nuclear energy, needs to be solved. The author previously reported that the amount and properties of radioactive waste depend on upstream conditions in the nuclear fuel cycle such as the nuclear reactor operation, reprocessing, and storage [1,2]. Integrated nuclear fuel cycle simulation spanning from the front- to the back-end is therefore required in order to design advanced nuclear energy systems and perform quantitative analysis of future scenarios.

About 30 codes for simulating the nuclear fuel cycle have previously been developed by various research institutes, and some international benchmark studies have been conducted in the last decade [3–5]. Table 1 shows the nuclear fuel cycle simulation codes of each institute in the report published by the Nuclear Energy Agency (NEA, under the Organization for Economic Co-operation and Development) in 2012 [5], as well as the fuel cycle code ANICCA, developed by SCK CEN [6]. It can be said that each of the codes excels at tracking actinide nuclides and performing nuclear fuel cycle simulation, particularly in relation to the front-end and reactor. However, the conditions of back-end simulation are considered to be limited using these codes because they do not track the necessary nuclides to analyze the diverse and variable scenarios. Although FAMILY-21 [7] and NMB3.0 [8] were developed to simulate the Japanese nuclear fuel cycle, they are specialized for mass balance analysis of actinide nuclides for mainly using fast breeder reactors (FBR) and accelerator driven nuclear transmutation system (ADS). Therefore, nuclear fuel cycle codes that flexibly simulates not only the front-end but also the back-end has not been developed.

In order to realize an integrated nuclear fuel cycle simulator that covers all stages from mining to nuclide migration after disposal, Tokyo Institute of Technology (Tokyo Tech) and the Japan Atomic Energy Agency (JAEA) began jointly developing Nuclear Material Balance version 4.0 (NMB4.0) in 2019. NMB4.0 is designed to provide flexibility and an all-in-one packaged nuclear fuel cycle simulation for users as open code. Various ideas have been implemented to improve user-friendliness. The main technical features of NMB4.0 are as follows.

• 179 nuclides, which were selected to enable nuclear fuel cycle simulation including not only the front-end but also the back-end, are tracked.

• The Okamura explicit method (OEM) is implemented for depletion calculations including short half-life nuclides with very low calculation cost and sufficient accuracy.

• The amount of waste and footprint of the geological repository, which are related to back-end simulations, are simulated using multiple approaches from waste production treatment to thermal analysis of the geological repository.

• The accuracy of depletion calculations agrees with ORIGEN [9], and the back-end simulation can be calculated with an accuracy within 1% of ORIGEN and COMSOL Multiphysics® (COMSOL) [10].

The main objective of this paper is to introduce the updated implementation of the integrated nuclear fuel cycle simulation and verify NMB4.0 through benchmark studies. Section 2 gives an overview of NMB4.0, including general functions and the three core calculation modules (front-end, reactor, and back-end). Section 3 shows the results of verifying NMB4.0 against other simulation codes. Verification of the front-end and reactor simulations was conducted with reference to scenarios 1 to 3 in the NEA benchmark study. In addition, verification of back-end simulation was performed compared with ORIGEN and COMSOL, which are general methods for static nuclear fuel cycle simulation.

2 Description of NMB

2.1 Basic structure

Figure 1 shows an overview of NMB4.0, which consists of three modules: front-end, reactor, and back-end. NMB4.0 was developed for dynamical and integrated analyses of the nuclear fuel cycle, and performs the depletion, material balance, thermal analysis of the geological repository, and nuclide migration from the repository to the environment under the given reactor operating conditions and reprocessing parameters. The electric generation, material balance, waste properties, footprint of the geological repository, and nuclide migration are then obtained as output data.

NMB4.0 runs in Microsoft Excel®, the programming language is mainly VBA, and some parts of the calculations are written in C++ for faster calculation. One of the targets of developing NMB4.0 is to distribute the nuclear fuel cycle simulator to many stakeholders, including policymakers, industry leaders, researchers, and students, in order to provide an opportunity for thinking about future nuclear energy utilization. Hence, NMB4.0 uses Excel VBA, which is widely and easily used by many people.

NMB4.0 is a kind of “tag-based simulator” that provides users with flexible nuclear fuel cycle simulation. Various nuclear fuel cycle scenarios can be realized by linking each process with user-defined tags. The simulation period can basically be defined by the user. For example, the simulation period is 10, 100, or 500 years.

Fig. 1 Overview of NMB4.0.

2.2 Tracked nuclides

There were 26 actinide nuclides selected as major nuclides with a half-life or more than 2 days which are associated with the burnup/decay chain, and 153 fission product (FP) nuclides were selected from the 1200 nuclides handled by ORIGEN. The 153 FP nuclides were selected to satisfy two criteria. First, nuclides were selected to reproduce ORIGEN by more than 99.9% in 5 evaluation indexes that are important for back-end analysis (mass, decay heat, radioactivity, toxicity, prohibited elements into vitrification). The burnup and decay calculations were carried out under 4 conditions such as Light Water Reactor (LWR)-UO₂, LWR-MOX, fast breeder reactors (FBRs)-MOX and FBR-Blanket. Second, the FP nuclides also include those required to configure a simple burnup chain with the same calculation accuracy as ORIGEN in the NMB code.

The mass of nuclides, which are not tracked in this code, is conserved into daughter or granddaughter nuclides through searching the decay chain of the nuclide with reference to nuclear data library. For this reason, the mass of nuclides not tracked in the code are taken into account. Further details and list of selected FP nuclides were reported in [18].

2.3 Nuclear data library

As shown in Figure 2, the nuclear data libraries used in NMB4.0 were created using the “Nuclear Data Library Editing System” framework, which is provided as an external code. This structure enables easy building of a one-group cross-section of libraries and burnup/decay chain matrices by selecting the nuclides of interest and the ORIGEN formatted library for each reactor type (e.g., pressurized water reactors (PWRs) and FBRs). The basic nuclear data in this code were taken from ORLIBJ40 [19], which is the nuclear data library for ORIGEN and based on JENDL4.0. Since the fission yield contained in ORLIBJ40 is given for only 7 nuclides, NMB4.0 employs JENDL FPY [20] for fission yield and GEFY 7.2 [21] for spontaneous fission yield. This framework supports incorporation of original core information designed using nuclear core design codes such as MARBLE and SRAC.

Fig. 2 Flowchart of vector and matrix creation for depletion calculations.

2.4 Front-end module (mining, enrichment, and fuel fabrication)

The enrichment of UO₂ fuel and the fissionable Pu content of MOX fuel are determined using the infinite multiplication factor (k_∞) set for each reactor and equation (1). This makes it possible to consider the needs of natural uranium and enrichment of UO₂ and MOX fuel depending on the core, fuel, and stored Pu isotopic composition in the scenario being analyzed. The cross-section used in the k_∞ calculation can take into account burnup dependencies.

NMB4.0 provides conditions that can be configured for user-defined fabrication of MOX fuel. These are as follows.

• Setting the used material (multiple materials can be mixed)

• Setting the order of materials used (e.g., Last In First Out (LIFO), First In First Out (FIFO), uniformly)

• Setting an upper limit on Pu content in MOX fuel

• Setting either the beginning of cycle (BOC) or end of cycle (EOC) as the target value of k_∞ in fuel production. BOC and EOC target values are provided since the nuclear reaction behavior differs depending on the Pu composition in the case of recycling of Pu. The calculation cost is relatively high owing to iteration of the depletion calculation when calculating the Pu content in MOX fuel by EOC.

$$k_{\infty }=\frac{{\sum }_i^{n}v{\sum }_{f,i}{w}_i}{{\sum }_i^n\left({\sum }_{c,i}+{\sum }_{f,i}{w}_i\right)}\text{,}$$(1)

2.5 Reactor module

2.5.1 Reactor/fuel

NMB4.0 can calculate multiple reactors in parallel. Moreover, various combinations of reactor operating conditions such as burnup and batch length can be considered based on user definitions. The standard nuclear reactor cores in NMB4.0 are LWRs (boiling water reactors (BWRs) and PWRs), CANDU reactors, gas-cooled reactors, sodium-cooled FBRs, molten salt reactors, and ADSs. These cores can be combined with multiple fuels such as UO₂, MOX, UO₂ blanket, and nitride fuel. If the user wants to calculate another design, including original designs, the user can calculate the core data using NMB4.0 by using the framework discussed in Section 2.3. The depletion calculation uses 6 nuclear reaction cross sections to ground state ((n, γ), (n, fission), (n, 2n), (n, 3n), (n, p), (n, α)), and 2 cross sections to metastable state ((n, γ), (n, 2n)). In addition, the cross-section used in the depletion calculation can take into account burnup dependencies. The combination of the reactor core and fuel can be flexibly considered by user definition. However, the feasibility of the combinations was not investigated, except for common combinations such as LWRs and UO₂.

2.5.2 Decay and depletion calculations

Nuclide decay and depletion calculations in NMB4.0 are obtained by solving the following equation:$$\frac{d{\displaystyle \overrightarrow{N}}}{\mathrm{}\mathrm{dt}}=A{\displaystyle \overrightarrow{N}}\text{,}$$(2) $$A=D+\frac{P}{{\sum }_i{P}_{F,i}{F}_i{N}_i}\left(C_{fission}+C_{capture}+C_{n,2n\cdots }\right)\text{,}$$(3)where ${\displaystyle \overrightarrow{N}}=\left(\begin{array}{c}\hfill N_1\hfill \\ \hfill ⋮\hfill \\ \hfill N_E\hfill \end{array}\right)$, nuclide quantity vector; E, Number of nuclides; P, reactor power; ${\displaystyle \stackrel{}{P_F}}$, vector of a fission energy; D, decay chain matrix; ${\displaystyle \overrightarrow{F}}$, vector of fission cross-section per neutron; C_x, burnup chain matrix for reaction x; A, transition matrix.

In the case of nuclide decay, the reactor power P = 0 and A=D independent of the reactor. For example, if the decay matrix D_1year, which represents the decay matrix after one year, is calculated in advance, the nuclide quantities after one year of decay can be calculated using the following equation:$${\displaystyle \overrightarrow{N}}\left(T+1year\right)=D_{1year}{\displaystyle \overrightarrow{N}}\left(T\right)\text{.}$$(4)

NMB4.0 provides multiple decay matrices with different decay times and the decay calculation is performed using these matrices. Advance preparation of these decay matrices contributes greatly to reducing the calculation time.

Depletion calculations in NMB3.0 were performed using the matrix exponential method in which the nuclide composition in fuel is updated every time step Δt using the following equation:$${\displaystyle \overrightarrow{N}}\left(T+\mathrm{\Delta }t\right)=\left(I+A\mathrm{\Delta }t+\frac{1}{2}{\left(A\mathrm{\Delta }t\right)}^2+\frac{1}{3}{\left(A\mathrm{\Delta }t\right)}^3+\cdots \right){\displaystyle \overrightarrow{N}}\left(T\right)\text{.}$$(5)

When the nuclide i has a short half-life, the i-th diagonal element of the burnup matrix A is approximately the negative value of the decay constant −λ_i. In the first-order approximation of matrix exponential method, the nuclide vector ${\displaystyle \overrightarrow{N}}$ is expressed as ${\displaystyle \overrightarrow{N}}\left(T+\mathrm{\Delta }t\right)\cong \left(I+A\mathrm{\Delta }t\right){\displaystyle \overrightarrow{N}}\left(T\right)$. It is approximated as $\left(1-{\lambda }_i\mathrm{\Delta }t\right){\displaystyle \overrightarrow{N_i}}\left(T\right)$ for short half-life nuclide i. When λ_iΔt is more than 1, the nuclide number density becomes negative, and the calculation diverges. To mitigate the divergence of the calculation, it is necessary to use the higher order approximation and/or sufficiently small Δt in the matrix exponential method. However, both require large computational costs. In addition, depletion calculations account for most of the calculation time in a dynamic fuel cycle simulation because depletion calculations are performed several thousand times per scenario. Therefore, it was necessary to find a method to solve burnup calculations involving short half-life nuclides with as much accuracy as possible and with low computational cost. This led to the development of OEM.

In OEM, the corrected time step $\mathrm{\Delta }{\displaystyle \tilde{t}}$ is given so that ${\lambda }_i\mathrm{\Delta }{\displaystyle \widetilde{t_i}}<1$ is valid even when the diagonal elements of burnup matrix take large negative values, by giving a time step for each nuclide to take the same number of densities as the analytical solution. The corrected time step $\mathrm{\Delta }{\displaystyle \tilde{t}}$ is represented in following equation:$$\mathrm{\Delta }{\displaystyle \widetilde{t_i}}=\frac{1-\exp \left(-{\lambda }_i\mathrm{\Delta }t\right)}{{\lambda }_i}\text{.}$$(6)

The corrected time step $\mathrm{\Delta }{\displaystyle \widetilde{t_i}}$ is given so that it has the same number density as the analytical solution for each nuclide. With the definition of a vector $\mathrm{\Delta }{\displaystyle \overrightarrow{t}}$ which consists of $\mathrm{\Delta }{\displaystyle \widetilde{t_i}}$, the normal and transposed OEM are defined as equations (7) and (8), respectively:$${\displaystyle \overrightarrow{N}}\left(T+\mathrm{\Delta }t\right)={\displaystyle \overrightarrow{N}}\left(T\right)+A\left({\displaystyle \overrightarrow{N}}\left(T\right)\circ \mathrm{\Delta }{\displaystyle \overrightarrow{t}}\right)\text{,}$$(7) $$N_i\left(T+\mathrm{\Delta }t\right)=N_i\left(T\right)+\mathrm{\Delta }{\displaystyle \overrightarrow{t}}°\left(A{\displaystyle \overrightarrow{N}}\left(t\right)\right)\text{,}$$(8)where a ° b takes element-wise product of a and b. The features of OEM are summarized below.

• OEM can reduce the calculation time compared to the matrix exponential method since OEM is a first order approximation and can use a longer time step.

• About 10 to 10³ time steps, depending on the type of OEM, are required to obtain the sufficient accuracy for the purpose of nuclear fuel cycle simulation.

• Transposed OEM has better accuracy than normal OEM when the daughter nuclides have a larger decay constant than their parent nuclides.

• The computational cost considering the sparse matrix in OEM is O (E²) ⋅ β when the calculated number of nuclides is E and the percentage of the number of elements in the burnup matrix is β (≈ 15%), which is fundamentally smaller than CRAM (O (E³)).

Further details of OEM such as the derivation and the comparison of computational performance with other depletion calculation methods were reported in [22].

2.6 Back-end module (reprocessing, stabilization, storage, and geological disposal)

The back-end module was the biggest update in NMB4.0. Figure 3 shows the structure of the back-end module. The back-end module consists of reprocessing, partitioning, stabilization, and geological disposal sections. The general features of the back-end module are as follows.

• Performance of waste partitioning: waste partitioning within the range of nuclides, including those selected for back-end simulation handled in NMB4.0 can be assumed by user definition.

• Flexibility of back-end simulation: although there are some limitations, reprocessing, waste partitioning, stabilization, storage, and geological disposal conditions can be set flexibly by user definition.

• Agreement of accuracy of the back-end simulation with static simulation using ORIGEN and COMSOL: details are described in Section 3.

Fig. 3 Back-end module structure.

2.6.1 Reprocessing

The role of the reprocessing section in this code is to separate the elements that are recycled to the front-end fuel cycle from the elements that become waste. This means that U, Pu, minor actinides (MAs), and waste are separated in the reprocessing section, with U, Pu, and MAs recycled and everything else going to waste management, which is handled by the waste partitioning section.

The reprocessing specification (including recovery ratio), type of spent fuel to be reprocessed, and schedule can be defined by the user. Not only actinides but also FPs can be reprocessed in the fuel cycle. For example, it is very difficult to separate MAs from rare earths (REs) by current MA separation technology. NMB guarantees that separation conditions can be configured that take into account conditions such as REs accompanying MAs. The fuel fabrication is then calculated by taking into account impurities by using the structure discussed in Section 2.4.

The reprocessing calculation realizes consideration of mixed reprocessing of different types of spent fuels. For example, it is possible to mix spent fuels of BWRs and PWRs, or to mix UO₂ and MOX. Waste generated from mixed reprocessing can also be assumed. The ratios of mixed reprocessing of different fuels are user-defined.

2.6.2 Waste partitioning

As shown in Figure 3, NMB4.0 divides waste into separate types of waste and sets different process conditions for each type such as stabilization, predisposal storage and geological disposal. For example, waste generated from reprocessing is assumed to be separated into rare gas, iodine, and high-level waste (HLW) in the partitioning section according to the typical PUREX process. In advanced technologies, waste also consists of Sr-Cs and Mo waste, and so on. The separation ratios are defined by the user.

2.6.3 Stabilization

Depending on the specifications of the reprocessing process, multiple types of waste with different nuclide compositions are generated. In addition, the corresponding stabilization methods depend on the chemical and physical properties of the individually separated waste. Therefore, the stabilization methods, size, waste content in solid, storage period before disposal, and restrictions on stabilization methods are defined by the user. Basic stabilization methods that are provided include vitrification, sodalite, and cement.

In terms of the restrictions on stabilization methods, it is said that the upper limit on the waste content in solid is constrained by limits including heat generation and the content of inhibitory elements for stabilization. For example, although the amount of vitrified waste can be reduced by increasing the waste content, the waste content is determined by limits on the heat generation of the solid in place to prevent glass crystallization by high heat generation. Therefore, NMB4.0 allows user-defined restrictions on conditions such as the waste content and heat generation, the content of inhibitory elements, and the temperature of the buffer material at the time of geological disposal. This makes it possible to perform the back-end simulation depending on the characteristics or properties of the solid.

2.6.4 Geological disposal

2.6.4.1 Geological repository modeling

The waste generated by the nuclear fuel cycle is stabilized and stored for several decades and will eventually be disposed of deep underground. Analysis of the size of the geological repository is required for nuclear fuel cycle simulation. The footprint of the geological repository is estimated from the product of the amount of waste by the area needed to dispose of one waste package (occupied area). In the Japanese geological program [23,24], the occupied area is determined from the view point of disposal tunnel stability and the engineered barriered system (EBS) temperature limit. In particular, the temperature limit of the buffer material, 100 °C, is considered to be the most critical condition in the EBS. This temperature limit is set to prevent illitization of bentonite, the main material of the buffer material. When the stability of the disposal tunnel is guaranteed, the occupied area is designed so that the temperature of buffer material is below the limit. Therefore, to calculate the occupied area, it is required to analyze the temperature of the EBS, especially the buffer material, when the waste package is disposed deep underground.

The basic concept of the thermal analysis model of the EBS in NMB4.0 is presented. The geological disposal system is designed based on the results of coupled thermal-hydro-mechanical-chemical simulation (THMC) in the near field (nearby the EBS) [23,25]. However, according to the results of coupling analysis of thermal and hydro on the geological repository by Taniguchi et al., the movement speed of groundwater in the relatively early stage after geological disposal is very small, and the influence on the thermal conductivity of near-field components is limited [26]. In another report, coupling analysis of thermal, hydro, and mechanical was used to examine the characteristics of heat transfer, infiltration, and swelling of bentonite, which is the main raw material of the buffer material [27]. As a result, it was found that conservative results can be obtained by thermal analysis using the thermophysical properties in the hydrous state when EBS is stationed, although the thermophysical properties of the buffer material fluctuate owing to swelling by groundwater. Therefore, in the NMB code, the occupied area is calculated by performing a near-field thermal analysis using the thermophysical properties in the hydrous state without considering hydro, mechanical, and chemical effects.

Figure 4 shows the thermal analysis model. This model and the thermophysical properties were taken from the Japanese geological disposal program [23,24]. Thermal analysis at the geological repository can be obtained by solving the following equation:$$\sigma \rho \frac{d}{\mathrm{}\mathrm{dt}}u\left(r,t\right)=\kappa {\nabla }^{2}u\left(r,t\right)+Q\left(r,t\right)\text{,}$$(9)where u(r, t), temperature at location r and time t (K); σ, specific heat (J/K/kg); ρ, density (kg/m³); κ, thermal conductivity (W/K/m); Q (r, t), heat generation density (W/m³).

The heat generation density is described as the sum of the decay heat of nuclides,$$Q\left(r,t\right)={\sum }_i{Q}_i\left(r,t\right).$$(10)

The following equation is obtained by standardizing Q_i:$$Q_i\left(r,t\right)=A_i{q}_i\left(r,t\right)\text{,}$$(11)

q_i(r, t) normally has a uniform spatial distribution in the waste and its time evolution follows an exponential function when the nuclides released from the waste are ignored. The spatial and time distribution of the heat generation is therefore unique when a unit amount of a single nuclide exists in the waste. The same solution as in equation (9) is obtained by solving u_i for each nuclide using equation (12) below and taking the linear sum as in equation (13):$$\sigma \rho \frac{d}{\mathrm{}\mathrm{dt}}{u}_i\left(r,t\right)=\kappa {\nabla }^2{u}_i\left(r,t\right)+q_i\left(r,t\right)\text{,}$$(12) $$u={\sum }_i{A}_i{u}_i.$$(13)

The time distribution of temperature u_i (r_buffer, t) at the position where the temperature of the buffer material becomes the highest (r_buffer) is created in a database by solving equation (12) for each nuclide that affects the heat generation of waste after disposal using COMSOL, which is a general finite element method code. This database can be used to calculate the time course of the temperature at r_buffer corresponding to any nuclide composition. The database was created for the four waste disposal methods (vertical, horizontal, PEM, compacted) used in the Japanese geological disposal program by taking into considering the mechanical stability of the disposal tunnel [23,24]. The thermophysical properties in the hydrous state for the thermal analysis model are taken from reports [24]. By using this database, calculation with an accuracy the same as COMSOL can be conducted with a greatly reduced calculation cost. The temperature limit of buffer material can be set by user defined.

Fig. 4 Model of thermal analysis for geological repository.

2.6.4.2 Nuclide migration

Long-term safety assessment of geological disposal systems is conducted through analysis of nuclide migration considering the characteristics of the geological environment, waste properties, and design options related to geological repository. The migration behavior of tracked nuclides released into the environment from waste packages placed in the deep geological repository and the radiation effects by assuming a suitable nuclide migration model for the actual environment of the geological repository can be evaluated.

Analysis of nuclide migration of HLW in Japan was conducted using the single waste condition calculated by ORIGEN [23,24]. However, it is predicted that the waste condition changes depending on the upstream conditions of the nuclear fuel cycle. Therefore, the analysis function of nuclide migration is implemented in the NMB code to successively evaluate the nuclide migration for a variety of nuclear fuels. The nuclide migration is calculated using the database prepared in the code. As shown in equation (14), the release rate D_i(r, t) of nuclide i at location r and time t is evaluated by solving the nuclide migration model with a nuclide composition in waste x_i in the basic analysis. However, it is difficult to solve the complex model including THMC because of the calculation cost. In order to solve the nuclide migration simply, equation (14) was approximated by equation (15) with D_i for discrete x_j, r, and t prepared in the code by solving equation (15) using GoldSim^®, which is commonly used for nuclide migration analysis [28]. The nuclide migration model was taken from the H12 report [29]. This method cannot collectively evaluate the interrelationships between some nuclides, such as accurate decay chains and solubility limiting concentrations of nuclides. However, nuclide migration analysis with reasonable accuracy and high speed is realized in NMB4.0 by introducing the database. It is thought that this function is a major advancement in the field of fuel cycle simulation.$${\displaystyle \overrightarrow{D}}\left(r,t;x_1,x_2,\cdots x_n\right)=\left(\begin{array}{c}\hfill D_1\left(r,t;x_1,x_2,\cdots x_n\right)\hfill \\ \hfill D_2\left(r,t;x_1,x_2,\cdots x_n\right)\hfill \\ \hfill ⋮\hfill \\ \hfill D_n\left(r,t;x_1,x_2,\cdots x_n\right)\hfill \end{array}\right)$$(14) $$D_i\approx {\sum }_j{D}_i\left(r,t;{\delta }_{jj}\cdot {{\displaystyle \overrightarrow{x}}}^t\right)$$ $${\delta }_{lm}=\{\begin{array}{c}\hfill 1,\left(l=m=j\right)\hfill \\ \hfill 0,\left(l\ne jorm\ne j\right).\hfill \end{array}$$(15)

3 Verification of NMB4.0 performance

3.1 Verification against NEA benchmark study

The computational performance of NMB4.0 was verified against other codes by referring to the results of the international benchmark study for nuclear fuel cycle simulation conducted at the NEA in 2012. The codes for comparison were FAMILY-21, COSI, EVOLCODE, VISION, and DESAE, which were summarized in the report, as well as ANICCA, which published the benchmark results referring to the same publication. There are three scenarios involving depletion calculations for PWR UO₂, PWR MOX, and fast reactors. Scenario 1 assumes only operation of PWRs with UO₂ fuel. Scenario 2 is a sort of mono-recycling scenario assuming the reprocessing of UO₂ fuel, and MOX fuel is burned in PWRs. Scenario 3 is a case of transition from PWRs to FBRs. Pu and MAs are assumed to be recycled. The benchmark study in all scenarios were conducted by NMB4.0, and the results are shown below. Detailed calculation conditions for each scenario are not described here.

3.1.1 Scenario 1

For scenario 1, shown in Table 2, the calculation results from NMB4.0 agree with the other codes in all evaluation indexes. The results are slightly different between NMB4.0 and DESAE2.2. According to the report, it was considered that this difference is considered by difference of the assumption of the fuels loaded into reactor core is different.

3.1.2 Scenario 2

The calculation results for scenario 2 are shown in Table 3 and Figure 5. The reprocessing amount shown in Table 3 was the average value of the reprocessing amount for 120 years. In scenario 1, although the calculation results from NMB4.0 agreed with the other codes in respect of front-end indexes, different values were observed in the average annual reprocessing amount. This is because different interpretations exist on the loading and exchange of fuel in each code. Regarding the amount of transuranic waste (TRU), NMB4.0 can assume both cumulation and decay of tracked nuclides in all processes the same as COSI. Thus, the results calculated under both conditions are shown in Figure 5. The trends in NMB4.0 agreed with the other codes. The NMB4.0 result values were close to COSI except for the amount of Am, for which there was a difference of about 15 tons between NMB4.0 and COSI. It is thought that this is because the amount of Pu-241, which is the parent nuclide of Am-241, is different in the depletion calculation owing to differences in nuclear data, particularly cross-sections.

Fig. 5 Results of scenario 2.

3.1.3 Scenario 3

Figure 6 shows the annual production of fuels and the amount of TRU. In the annual production of FBR MOX fuel and blankets, similar trends of results were obtained by NMB4.0 and the other codes. In particular, the results from NMB4.0 agreed with FAMILY-21 and DESAE 2.1. However, the output of 0 every 8 years as seen in COSI did not occur in NMB4.0. According to the report, the cause of this trend is thought to be that COSI cannot assume different burnup batches for fuels loaded into a single reactor because the fast reactor is calculated as a single fission reactor in this version, whereas the oscillation in COSI did not occur in NMB4.0 since the MOX fuel and blankets could be calculated individually. The trend in the amount of TRU calculated by NMB4.0, except for the amount of Am, agreed with the other codes. It is thought that the cause of the difference in Am is not due to differences between the scenarios because the slope of the time course of amount of Am during mono-recycling was the same in scenarios 2 and 3. Therefore, the difference in Am is presumed to be differences in the amount of Pu-241 caused by cross-sections.

Fig. 6 Results of scenario 3.

3.2 Verification with ORIGEN and COMSOL

Verification of the performance of the back-end module, which is one of the main features of NMB4.0, is shown here. NMB4.0 was compared with the general static method that is commonly used for back-end simulation [23,24]. Depletion and decay calculations were conducted using ORIGEN, and thermal analysis of the geological repository was performed using COMSOL. Since the static method simulates equilibrium conditions, it was compared with the simulation result in the equilibrium period even in NMB4.0. As shown in Table 4, three benchmark conditions were prepared. S1 is a condition close to the reference case in the Japanese geological disposal program. S2 is a condition that extends the cooling period of spent fuel. S3 is a condition assuming high burnup. In the depletion calculation by ORIGEN, ORLIBJ40 was used as the cross-section library, and the specific power was 37.5 MW/tHM. The initial fuel enrichment was 4.1 wt.% for S1 and S2 and 4.7 wt.% for S3. The reprocessing conditions were assumed to be the typical PUREX process. The evaluated indexes are the number of vitrified wastes generated, heat generation at the time of vitrification and disposal, radioactivity at the time of vitrification and disposal, Mo content in vitrified waste, and maximum temperature of buffer material. Each scenario has a different restriction on the stabilization. S1 has waste content (=10 wt.%), S2 has Mo content (≤1.5 wt.%), and S3 has initial heat generation (≤2.3 kW).

The calculation results are shown in Table 5. Differences in the calculation are within 1% for all evaluation indexes. The difference in the number of vitrified wastes, which is the amount of waste, was about 0.2% for S1 and S2 and 0.4% for S3, assuming high burnup. In the verification of the thermal analysis method of the geological repository shown in Section 2.6.4.1, it was found that the difference in the maximum temperature of the buffer material was about 0.4%, which is sufficiently good performance. The time evolution of the temperature of the buffer material is shown in Figure 7 for cases S1 and S2 calculated by both codes. In S1, the peak at around 20 years after disposal was reproduced. In S2, the peak at around 300 years after disposal was also reproduced. These time variations in temperature of buffer material tend to be the same as the analysis results performed in a previous study [2]. According to the studies, the effects of Cs and Sr occur in a relatively early period after disposal like in S1, and the peak in S2 occurs due to Am-241 from beta-decay of Pu-241. Therefore, the performance of thermal analysis method in NMB4.0 was clarified to consider the decay heat. In addition, it can be seen that NMB4.0 agrees with the static calculation method of ORIGEN and COMSOL considering that the depletion calculation method and the burnup chain are different.

Finally, the performance of nuclide migration analysis, that is, the time course of the release rate at the EBS, is shown as an example case in Figure 8. The nuclide composition in waste generated from S1 was assumed. The release rate of each nuclide was expressed corresponding to the nuclide composition of the input. The general trends are consistent, although different trends occurred in some actinide nuclides compared with reference [27]. Such calculations can be performed on a time scale of up to 10⁸ at three locations, namely, the EBS, host rock, and biosphere, based on the model from the H12 report.

Fig. 7 Time course of temperature of buffer material on S1 and S2.

Fig. 8 Example case for performance of nuclide migration analysis in NMB4.0 (Time course of the release rate of each nuclide at the EBS in the case of S1).

4 Conclusion

NMB4.0, which enables integrated analysis of the entire nuclear fuel cycle from the front- to back-end, was developed and introduced in this paper. NMB4.0 was verified against existing nuclear fuel cycle simulation codes. The features of NMB4.0 are summarized below.

• 179 nuclides (26 actinides and 153 FPs) selected from the viewpoint of burnup chain and nuclides required for mass balance of waste management are tracked throughout the nuclear fuel cycle.

• The nuclear data library obtained from the core analysis can be utilized with the nuclear data library editing system.

• Light water reactors, CANDU reactors, gas reactors, sodium-cooled FBR, and ADS are available. Each core can be combined with fuels such as UO₂, MOX and nitride fuel.

• OEM is implemented in depletion calculation to reduce the calculation time while maintaining the accuracy of short half-life nuclides.

• A module is implemented for flexible back-end simulation. The accuracy of the back-end simulation differs from ORIGEN and COMSOL by about 1%.

• Nuclide migration in geological disposal can be calculated using a prepared database

As a result of the development and verification of these functions, it was shown that NMB4.0 is a code that can perform integrated analysis of a wide range of sections of the nuclear fuel cycle from mining to nuclide migration.

Our next step is to prepare a database of indexes related to the nuclear fuel cycle such as economic efficiency, environmental impact, and risk, and to expand the calculation scope of the code. In addition, improvement of the framework on nuclide migration study will be conducted through collaborating with geological disposal experts. Finally, we will continue development aiming for an integrated nuclear fuel cycle code platform.

Funding

This work was performed through the “Development of evaluation method of mass-balance in back-end process of nuclear fuel cycle and implementation to NMB code-1,2,3” project, a collaborative research effort between TokyoTech&JAEA funded by the Laboratory for Advanced Nuclear Energy, Tokyo Institute of Technology (2019-11, 2020-14, 2021-10). In addition, this work was supported by a Grant-in-Aid for JSPS Research Fellows (20J14956).

Author contribution statement

All authors contributed to the development of this code as a project member. Mr. Tomohiro Okamura was primarily responsible for its development and writing. Dr. Kenji Nishihara directed this work, and Mr. Ryota Katano contributed to the derivation of the OEM.

Acknowledgements

We are grateful to Dr. Ivan MERINO for providing the verification results discussed in Section 3.

References

All Tables

All Figures

	Fig. 1 Overview of NMB4.0.
In the text

	Fig. 2 Flowchart of vector and matrix creation for depletion calculations.
In the text

	Fig. 3 Back-end module structure.
In the text

	Fig. 4 Model of thermal analysis for geological repository.
In the text

	Fig. 5 Results of scenario 2.
In the text

	Fig. 6 Results of scenario 3.
In the text

	Fig. 7 Time course of temperature of buffer material on S1 and S2.
In the text

	Fig. 8 Example case for performance of nuclide migration analysis in NMB4.0 (Time course of the release rate of each nuclide at the EBS in the case of S1).
In the text