Preliminary analysis of reactivity reconstruction capability based on inverse kinetics method under different initial reactivity states

Preliminary analysis of reactivity reconstruction capability based on inverse kinetics method under different initial reactivity states 1 2 first-author 1 1 1 1 corresp E-mail address of corresponding author: qi.yang@fds.org.cn E-mail address of corresponding author: qi.yang@fds.org.cn Key Laboratory of Neutronics and Radiation Safety, Institute of Nuclear Energy Safety Technology, Chinese Academy of Sciences, Hefei 230031, China; Key Laboratory of Neutronics and Radiation Safety, Institute of Nuclear Energy Safety Technology, Chinese Academy of Sciences, Hefei 230031, China; University of Science and Technology of China, Hefei 230027, China University of Science and Technology of China, Hefei 230027, China On-line reactivity monitoring plays an important role in operation and safety analyses of fission reactor systems. The inverse kinetics method, which is based on a point kinetics model, is the most widely used method for reactivity reconstruction of critical water reactors. However, this method is seldom applied to the reactivity reconstruction of subcritical reactors. In this study, an inverse kinetics method was employed for the reactivity reconstruction of a lead-based reactor under different initial reactivity states ( ρ 0 = 0, −2786, −5486, −8367, and −12371 pcm). The results showed that the deviation in the reactivity of the lead-based subcritical reactor was greater when ρ 0 became smaller. The reactivity reconstructed using the inverse kinetics method was globally underestimated. At a given reactivity perturbation, the relative and absolute errors increased with the decrease in the initial reactivity. At a given initial reactivity, with the increase in the reactivity perturbation, the absolute error increased, whereas the relative error remained the same. This deviation is due to the variation in the external neutron source, spatial-spectral effects, and sub-diffusive effects, which require further study. Reactivity reconstruction Preliminary analysis Inverse kinetics method Lead-based reactor 1 Introduction 1 On-line reactivity monitoring is fundamental to the safe operation of nuclear power systems. For critical reactors in normal operation, the reactivity is maintained at approximately zero using multiple control strategies; moreover, it is always lower than the effective delayed neutron fraction to avoid prompt criticality. For subcritical systems with an extraneous neutron source, the reactivity is maintained within a particular subcriticality range. Under no circumstances should the variation in the reactivity lead to a delay criticality, which is used to characterize the inherent safety of source-driven subcritical systems [ 1 1 ]. The following dual operation modes are implemented in lead-based reactors: critical and subcritical. Owing to their excellent nuclear, thermal-hydraulic, and inherent safety properties, lead-based reactors are considered the most promising candidates for fourth-generation nuclear energy systems [ 2 2 , 3 3 ], accelerator-driven sub-critical systems (ADSs) [ 4 4 , 5 5 ], and fusion-driven sub-critical systems (FDSs) [ 6 6 - 9 9 6-9 ]. The Institute of Nuclear Energy Safety Technology·FDS Team, Chinese Academy of Sciences, has made considerable progress in lead-based (pure lead or lead-based alloy coolant) reactors, including various designs of liquid-lead-lithium-cooled fusion reactors [ 10 10 - 13 13 10-13 ] and China LEAd-based reactor (CLEAR) [ 3 3 - 5 5 3-5 ]. As lead-based reactors are designed for both critical and different subcritical states, it is crucial to employ a robust method to monitor the reactivity under different initial states. The inverse kinetics method can respond to delayed neutrons, which is considered to be a general method for reconstructing the reactivity of critical water reactors[ 14 14 ]. However, the inverse kinetics method is seldom applied to subcritical reactivity reconstruction [ 15 15 ]. Moreover, whether the inverse kinetics method can be applied to subcritical reactivity reconstruction under different initial reactivity states is unclear. Hence, in this study, the inverse kinetics method is employed for reconstructing the reactivity of a lead-based reactor under different initial reactivity states. In the following sections, the principle behind the inverse kinetics method is first introduced. Second, for a lead-based critical reactor, the evaluation is conducted using power data obtained via the point kinetic calculation. For a lead-based subcritical reactor, the method is analyzed using power data obtained via the space–time kinetics calculation. The reactivity transient under different initial reactivity states is employed to quantitatively verify the applicability of the proposed method. 2 Principle behind the method 1 The exact point-kinetics equation along with six delayed neutron precursor terms is given as follows [ 16 16 ]: <math display="block" id="M1"><mrow><mfrac><mrow><mi>d</mi><mi>p</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow></mrow><mrow><mi>d</mi><mi>t</mi></mrow></mfrac><mo>=</mo><mfrac><mrow><mi>ρ</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow><mo>−</mo><mi>β</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow></mrow><mrow><mi>Λ</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow></mrow></mfrac><mi>p</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow><mo>+</mo><mstyle displaystyle="true"><munderover><mo>∑</mo><mrow><mi>k</mi><mo>=</mo><mn>1</mn></mrow><mn>6</mn></munderover><mrow><msub><mi>λ</mi><mi>k</mi></msub><msub><mi>ξ</mi><mi>k</mi></msub></mrow></mstyle><mo>+</mo><mi>s</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow><mtext> </mtext><mo>,</mo></mrow></math> <math display="block" id="M2"><mrow><mfrac><mrow><mi>d</mi><msub><mi>ξ</mi><mi>k</mi></msub><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow></mrow><mrow><mi>d</mi><mi>t</mi></mrow></mfrac><mo>=</mo><mfrac><mrow><msub><mi>β</mi><mi>k</mi></msub><mi>p</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow></mrow><mrow><mi>Λ</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow></mrow></mfrac><mo>−</mo><msub><mi>λ</mi><mi>k</mi></msub><msub><mi>ξ</mi><mi>k</mi></msub><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow><mo>.</mo></mrow></math> Here, p ( t ) can be considered the neutron density (neutrons·com -3 ) or the integrated power <math id="M3"><mrow><mrow><mo>(</mo><mi>W</mi><mo>)</mo></mrow></mrow></math> . ρ ( t ), β ( t ), and Λ ( t ), are the reactivity <math id="M4"><mrow><mrow><mo>(</mo><mrow><mtext>pcm</mtext></mrow><mo>)</mo></mrow></mrow></math> , effective delayed neutron fraction <math id="M5"><mrow><mrow><mo>(</mo><mrow><mtext>pcm</mtext></mrow><mo>)</mo></mrow></mrow></math> , and prompt neutron generation time <math id="M6"><mrow><mrow><mo>(</mo><mi>s</mi><mo>)</mo></mrow></mrow></math> , respectively. ξ k ( t ) is the concentration of the delayed neutron precursor (neutrons·com -3 ), and λ k is the decay constant of the precursors. s ( t ) denotes the <math id="M7"><mrow><mrow><mo>(</mo><mrow><mi>γ</mi><mo>,</mo><mo> </mo><mi>n</mi></mrow><mo>)</mo></mrow></mrow></math> or <math id="M8"><mrow><mrow><mo>(</mo><mrow><mi>α</mi><mo>,</mo><mo> </mo><mi>n</mi></mrow><mo>)</mo></mrow></mrow></math> source when the critical reactor starts or when an extraneous high intensity neutron source is detected for a subcritical reactor. The above equation can be inversely transformed into Eq. (3), which is the basis for the reactivity measurement [ 16 16 ]. <math display="block" id="M9"><mrow><mi>ρ</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow><mo>=</mo><mi>β</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow><mo>+</mo><mfrac><mrow><mi>Λ</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow></mrow><mrow><mi>p</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow></mrow></mfrac><mfrac><mrow><mi>d</mi><mi>p</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow></mrow><mrow><mi>d</mi><mi>t</mi></mrow></mfrac><mo>−</mo><mfrac><mrow><mi>Λ</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow></mrow><mrow><mi>p</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow></mrow></mfrac><munder><mstyle mathsize="140%" displaystyle="true"><mo>∑</mo></mstyle><mi>k</mi></munder><msub><mi>λ</mi><mi>k</mi></msub><msub><mi>ξ</mi><mi>k</mi></msub><mrow><mo>(</mo><mrow><msub><mi>t</mi><mn>0</mn></msub></mrow><mo>)</mo></mrow><msup><mi>e</mi><mrow><mo>−</mo><msub><mi>λ</mi><mi>k</mi></msub><mrow><mo>(</mo><mrow><mi>t</mi><mo>−</mo><msub><mi>t</mi><mn>0</mn></msub></mrow><mo>)</mo></mrow></mrow></msup><mo>−</mo><mfrac><mn>1</mn><mrow><mi>p</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow></mrow></mfrac><munder><mstyle mathsize="140%" displaystyle="true"><mo>∑</mo></mstyle><mi>k</mi></munder><msub><mi>λ</mi><mi>k</mi></msub><msub><mi>β</mi><mi>k</mi></msub><munderover><mstyle mathsize="140%" displaystyle="true"><mo>∫</mo></mstyle><mrow><msub><mi>t</mi><mn>0</mn></msub></mrow><mi>t</mi></munderover><mi>p</mi><mrow><mo>(</mo><mi>τ</mi><mo>)</mo></mrow><msup><mi>e</mi><mrow><mo>−</mo><msub><mi>λ</mi><mi>k</mi></msub><mrow><mo>(</mo><mrow><mi>t</mi><mo>−</mo><mi>τ</mi></mrow><mo>)</mo></mrow></mrow></msup><mi>d</mi><mi>τ</mi><mo>−</mo><mfrac><mrow><mi>s</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow><mi>Λ</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow></mrow><mrow><mi>p</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow></mrow></mfrac><mo>.</mo></mrow></math> It is difficult to determine the power history, i.e., the integral of the fourth term in the right-hand side of Eq. (3). To derive this term from the equation, some studies have used the Laplace transform and some others have used numerical procedures such as the Lagrange method [ 17 17 ], hamming method [ 18 18 ], and the method given elsewhere [ 19 19 ]. In this study, the power history term is written as follows. <math display="block" id="M10"><mrow><msub><mi>I</mi><mrow><mi>k</mi><mo>,</mo><mi>j</mi></mrow></msub><mo>=</mo><munderover><mstyle mathsize="140%" displaystyle="true"><mo>∫</mo></mstyle><mrow><msub><mi>t</mi><mn>0</mn></msub></mrow><mrow><msub><mi>t</mi><mi>j</mi></msub></mrow></munderover><mi>p</mi><mrow><mo>(</mo><mi>τ</mi><mo>)</mo></mrow><msup><mi>e</mi><mrow><mo>−</mo><msub><mi>λ</mi><mi>k</mi></msub><mrow><mo>(</mo><mrow><msub><mi>t</mi><mi>j</mi></msub><mo>−</mo><mi>τ</mi></mrow><mo>)</mo></mrow></mrow></msup><mi>d</mi><mi>τ</mi><mo>,</mo></mrow></math> <math display="block" id="M11"><mrow><msub><mi>I</mi><mrow><mi>k</mi><mo>,</mo><mi>j</mi><mo>+</mo><mn>1</mn></mrow></msub><mo>=</mo><munderover><mstyle mathsize="140%" displaystyle="true"><mo>∫</mo></mstyle><mrow><msub><mi>t</mi><mn>0</mn></msub></mrow><mrow><msub><mi>t</mi><mrow><mi>j</mi><mo>+</mo><mn>1</mn></mrow></msub></mrow></munderover><mi>p</mi><mrow><mo>(</mo><mi>τ</mi><mo>)</mo></mrow><msup><mi>e</mi><mrow><mo>−</mo><msub><mi>λ</mi><mi>k</mi></msub><mrow><mo>(</mo><mrow><msub><mi>t</mi><mrow><mi>j</mi><mo>+</mo><mn>1</mn></mrow></msub><mo>−</mo><mi>τ</mi></mrow><mo>)</mo></mrow></mrow></msup><mi>d</mi><mi>τ</mi><mo>,</mo></mrow></math> where, t j +1 = t j +Δ t , Δ t is the time interval, and t j denotes the discrete time instants for the j th point. After integrating by parts and omitting the higher-order terms of integration and differentiation, the following discrete recursive relationship, i.e., Eq. (6), is obtained. <math display="block" id="M12"><mrow><msub><mi>I</mi><mrow><mi>k</mi><mo>,</mo><mi>j</mi><mo>+</mo><mn>1</mn></mrow></msub><mo>=</mo><msub><mi>I</mi><mrow><mi>k</mi><mo>,</mo><mi>j</mi></mrow></msub><msup><mi>e</mi><mrow><mo>−</mo><msub><mi>λ</mi><mi>k</mi></msub><mi>Δ</mi><mi>t</mi></mrow></msup><mo>+</mo><mfrac><mrow><mi>n</mi><mrow/><mrow><mi>j</mi><mo>+</mo><mn>1</mn></mrow></mrow><mrow><msub><mi>λ</mi><mi>k</mi></msub></mrow></mfrac><mrow><mo>[</mo><mrow><mn>1</mn><mo>−</mo><mfrac><mn>1</mn><mrow><msub><mi>λ</mi><mi>k</mi></msub><mi>Δ</mi><mi>t</mi></mrow></mfrac><mrow><mo>(</mo><mrow><mn>1</mn><mo>−</mo><msup><mi>e</mi><mrow><mo>−</mo><msub><mi>λ</mi><mi>k</mi></msub><mi>Δ</mi><mi>t</mi></mrow></msup></mrow><mo>)</mo></mrow></mrow><mo>]</mo></mrow><mo>−</mo><mfrac><mrow><mi>n</mi><mrow/><mi>j</mi></mrow><mrow><msub><mi>λ</mi><mi>k</mi></msub></mrow></mfrac><mrow><mo>[</mo><mrow><msup><mi>e</mi><mrow><mo>−</mo><msub><mi>λ</mi><mi>k</mi></msub><mi>Δ</mi><mi>t</mi></mrow></msup><mo>−</mo><mfrac><mn>1</mn><mrow><msub><mi>λ</mi><mi>k</mi></msub><mi>Δ</mi><mi>t</mi></mrow></mfrac><mrow><mo>(</mo><mrow><mn>1</mn><mo>−</mo><msup><mi>e</mi><mrow><mo>−</mo><msub><mi>λ</mi><mi>k</mi></msub><mi>Δ</mi><mi>t</mi></mrow></msup></mrow><mo>)</mo></mrow></mrow><mo>]</mo></mrow><mo>.</mo></mrow></math> In the process of transient reactivity reconstruction, the kinetic parameters, such as the effective delayed neutron fraction and prompt neutron generation time, can be assumed time-independent except in the event of a major core damage accident. An optimal position for the neutron detector should be determined such that the effect of reactivity on the shape of the time-dependent flux is negligible. Based on the aforementioned assumptions, Eq. (3) can be transformed into Eq. (7), which is discrete. <math display="block" id="M13"><mrow><msub><mi>ρ</mi><mrow><mi>j</mi><mo>+</mo><mn>1</mn></mrow></msub><mo>=</mo><mi>β</mi><mo>+</mo><mfrac><mi>Λ</mi><mrow><msub><mi>p</mi><mrow><mi>j</mi><mo>+</mo><mn>1</mn></mrow></msub></mrow></mfrac><mfrac><mrow><msub><mi>p</mi><mrow><mi>j</mi><mo>+</mo><mn>1</mn></mrow></msub><mo>−</mo><msub><mi>p</mi><mi>j</mi></msub></mrow><mrow><mi>Δ</mi><mi>t</mi></mrow></mfrac><mo>−</mo><mfrac><mn>1</mn><mrow><msub><mi>p</mi><mrow><mi>j</mi><mo>+</mo><mn>1</mn></mrow></msub></mrow></mfrac><munder><mstyle mathsize="140%" displaystyle="true"><mo>∑</mo></mstyle><mi>k</mi></munder><msub><mi>β</mi><mi>k</mi></msub><mi>p</mi><mrow><mo>(</mo><mn>0</mn><mo>)</mo></mrow><msup><mi>e</mi><mrow><mo>−</mo><msub><mi>λ</mi><mi>k</mi></msub><mi>j</mi><mi>Δ</mi><mi>t</mi></mrow></msup><mo>−</mo><mfrac><mn>1</mn><mrow><msub><mi>p</mi><mrow><mi>j</mi><mo>+</mo><mn>1</mn></mrow></msub></mrow></mfrac><munder><mstyle mathsize="140%" displaystyle="true"><mo>∑</mo></mstyle><mi>k</mi></munder><msub><mi>λ</mi><mi>k</mi></msub><msub><mi>β</mi><mi>k</mi></msub><msub><mi>I</mi><mrow><mi>k</mi><mo>,</mo><mi>j</mi><mo>+</mo><mn>1</mn></mrow></msub><mo>−</mo><mfrac><mrow><msub><mi>s</mi><mrow><mi>j</mi><mo>+</mo><mn>1</mn></mrow></msub><mi>Λ</mi></mrow><mrow><msub><mi>p</mi><mrow><mi>j</mi><mo>+</mo><mn>1</mn></mrow></msub></mrow></mfrac><mo>.</mo></mrow></math> Thus, based on Eq. (7), an algorithm is coded with s j +1 = 0 for the critical reactor; and s j +1 = s j = s 0 , which can be calculated at the initial subcritical state, for the subcritical reactor, as shown in Eq. (8). For practical application, the initial reactivity ρ (0) could be measured using a pulsed neutron method at a lower power or using a source-jerk method at a higher power. The initial neutron density p (0) could be measured directly by monitoring the neutron flux at the optimal detector position. <math display="block" id="M14"><mrow><msub><mi>s</mi><mrow><mi>j</mi><mo>+</mo><mn>1</mn></mrow></msub><mo>=</mo><msub><mi>s</mi><mn>0</mn></msub><mo>=</mo><mfrac><mrow><mi>ρ</mi><mrow><mo>(</mo><mn>0</mn><mo>)</mo></mrow><mi>p</mi><mrow><mo>(</mo><mn>0</mn><mo>)</mo></mrow></mrow><mi>Λ</mi></mfrac><mo>.</mo></mrow></math> 3 Reactivity reconstruction for the lead-based reactor 1 To verify the applicability of the inverse kinetics method for the reactivity reconstruction of the lead-based reactor, a lead-based critical reactor with natural circulation [ 20 20 ] and a lead-based subcritical reactor, designed per the preliminary design studies of eXperimental accelerator-driven subcritical system (PDS-XADS) project [ 21 21 ], were selected as reference reactor models. 3.1 Lead-based critical reactor 2 A lead-based critical reactor with a capacity of 10 MWth [ 20 20 ] is designed to be a pool-type reactor driven via natural circulation. UO 2 with a mass enrichment of approximately 19.75% is selected as the fuel. At the beginning of the cycle, the effective multiplication factor k eff is 1.0142. The effective delayed neutron fractions from the k th group of precursors β eff, k are 0.00025, 0.00081, 0.00236, 0.00208, 0.00168, and 0.00026, corresponding to delayed neutron decay constants λ k of 0.0129, 0.0311, 0.134, 0.331, 1.26, and 3.21, respectively. The prompt neutron generation time is 0.521 μs. The reactivity is reconstructed using the inverse kinetics method with a time interval of 0.1 s for various power forms. The power data can be obtained using a numerical solution with a point kinetics code using a high-order end floating method [ 22 22 ] and via an analytical solution as the reference. Subsequently, the reconstruction results of the reactivity are compared to two numerical benchmark solutions and three analytical benchmark solutions, listed in Table 1 Table 1 . The absolute and relative errors are defined in Eqs. (9) and (10), respectively. 1 1 Reference solution 1 1 Power type 1 1 Reactivity type 1 1 Numerical solution 1 1 Point kinetics solution with high-order end floating method 1 1 ρ ( t )=20, 50, 100, 200 1 1   1 1   1 1 ρ ( t )=50sin(π t /2.5) 1 1 Analytical solution 1 1 p ( t )=e wt ( w =0.1599) 1 1 <math id="M15"><mrow><mi>ρ</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow><mo>=</mo><mi>β</mi><mo>+</mo><mi>Λ</mi><mi>w</mi><mo>−</mo><mstyle displaystyle="true"><munderover><mo>∑</mo><mrow><mi>k</mi><mo>=</mo><mn>1</mn></mrow><mn>6</mn></munderover><mrow><mrow><mo>[</mo><mrow><mfrac><mrow><msub><mi>λ</mi><mi>i</mi></msub><msub><mi>β</mi><mi>i</mi></msub><mo>+</mo><msub><mi>β</mi><mi>i</mi></msub><mi>w</mi><msup><mtext>e</mtext><mrow><mo>−</mo><mrow><mo>(</mo><mrow><msub><mi>λ</mi><mi>i</mi></msub><mo>+</mo><mi>w</mi></mrow><mo>)</mo></mrow><mi>t</mi></mrow></msup></mrow><mrow><msub><mi>λ</mi><mi>i</mi></msub><mo>+</mo><mi>w</mi></mrow></mfrac></mrow><mo>]</mo></mrow></mrow></mstyle></mrow></math> 1 1   1 1 p ( t )=1+ wt ( w =1) 1 1 <math id="M16"><mrow><mi>ρ</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow><mo>=</mo><mi>β</mi><mo>+</mo><mfrac><mrow><mi>Λ</mi><mi>w</mi><mo>−</mo><msubsup><mstyle mathsize="140%" displaystyle="true"><mo>∑</mo></mstyle><mrow><mi>k</mi><mo>=</mo><mn>1</mn></mrow><mn>6</mn></msubsup><msub><mi>β</mi><mi>i</mi></msub><mrow><mo>[</mo><mrow><msup><mtext>e</mtext><mrow><mo>−</mo><msub><mi>λ</mi><mi>i</mi></msub><mi>t</mi></mrow></msup><mo>+</mo><mfrac><mrow><mo stretchy="false">(</mo><msub><mi>λ</mi><mi>i</mi></msub><mo>−</mo><mi>w</mi><mo stretchy="false">)</mo><mrow><mo>(</mo><mrow><mn>1</mn><mo>−</mo><msup><mtext>e</mtext><mrow><mo>−</mo><msub><mi>λ</mi><mi>i</mi></msub><mi>t</mi></mrow></msup></mrow><mo>)</mo></mrow></mrow><mrow><msub><mi>λ</mi><mi>i</mi></msub></mrow></mfrac></mrow><mo>]</mo></mrow><msub><mi>λ</mi><mi>i</mi></msub><msub><mi>β</mi><mi>i</mi></msub><mo>+</mo><mi>w</mi><mi>t</mi></mrow><mrow><mn>1</mn><mo>+</mo><mi>w</mi><mi>t</mi></mrow></mfrac><mo>−</mo></mrow></math> 1 1   1 1 p ( t )=1+ A sin( wt ) ( A =1, w =0.3142) 1 1 <math id="M17"><mrow><mi>ρ</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow><mo>=</mo><mi>β</mi><mo>+</mo><mfrac><mrow><mrow><mo>{</mo><mrow><mi>A</mi><mi>w</mi><mi>Λ</mi><mtext>cos</mtext><mrow><mo>(</mo><mrow><mi>w</mi><mi>t</mi></mrow><mo>)</mo></mrow><mo>−</mo><msubsup><mstyle mathsize="140%" displaystyle="true"><mo>∑</mo></mstyle><mrow><mi>k</mi><mo>=</mo><mn>1</mn></mrow><mn>6</mn></msubsup><mrow><mo>[</mo><mrow><msub><mi>β</mi><mi>i</mi></msub><mo>+</mo><mfrac><mrow><mi>A</mi><msub><mi>β</mi><mi>i</mi></msub><msub><mi>λ</mi><mi>i</mi></msub><mrow><mo>[</mo><mrow><msub><mi>λ</mi><mi>i</mi></msub><mtext>sin</mtext><mrow><mo>(</mo><mrow><mi>w</mi><mi>t</mi></mrow><mo>)</mo></mrow><mo>−</mo><mi>w</mi><mtext>cos</mtext><mrow><mo>(</mo><mrow><mi>w</mi><mi>t</mi></mrow><mo>)</mo></mrow><mo>+</mo><mi>w</mi><msup><mtext>e</mtext><mrow><mo>−</mo><msub><mi>λ</mi><mi>i</mi></msub><mi>t</mi></mrow></msup></mrow><mo>]</mo></mrow></mrow><mrow><msub><mi>λ</mi><mi>i</mi></msub><msup><mrow/><mn>2</mn></msup><mo>+</mo><msup><mi>w</mi><mn>2</mn></msup></mrow></mfrac></mrow><mo>]</mo></mrow></mrow><mo>}</mo></mrow></mrow><mrow><mn>1</mn><mo>+</mo><mi>A</mi><mtext>sin</mtext><mrow><mo>(</mo><mrow><mi>w</mi><mi>t</mi></mrow><mo>)</mo></mrow></mrow></mfrac></mrow></math> <math display="block" id="M18"><mrow><mi>A</mi><mi>b</mi><mi>s</mi><mi>o</mi><mi>l</mi><mi>u</mi><mi>t</mi><mi>e</mi><mtext> </mtext><mtext> </mtext><mtext> </mtext><mi>e</mi><mi>r</mi><mi>r</mi><mi>o</mi><mi>r</mi><mo>=</mo><mrow><mo>|</mo><mrow><msub><mi>ρ</mi><mrow><mi>I</mi><mi>K</mi></mrow></msub><mo>−</mo><msub><mi>ρ</mi><mrow><mi>r</mi><mi>e</mi><mi>f</mi><mi>e</mi><mi>r</mi><mi>e</mi><mi>n</mi><mi>c</mi><mi>e</mi></mrow></msub></mrow><mo>|</mo></mrow><mo>,</mo></mrow></math> <math display="block" id="M19"><mrow><mi>R</mi><mi>e</mi><mi>l</mi><mi>a</mi><mi>t</mi><mi>i</mi><mi>v</mi><mi>e</mi><mtext> </mtext><mtext> </mtext><mtext> </mtext><mi>e</mi><mi>r</mi><mi>r</mi><mi>o</mi><mi>r</mi><mo>=</mo><mrow><mo>|</mo><mrow><mfrac><mrow><msub><mi>ρ</mi><mrow><mi>I</mi><mi>K</mi></mrow></msub><mo>−</mo><msub><mi>ρ</mi><mrow><mi>r</mi><mi>e</mi><mi>f</mi><mi>e</mi><mi>r</mi><mi>e</mi><mi>n</mi><mi>c</mi><mi>e</mi></mrow></msub></mrow><mrow><msub><mi>ρ</mi><mrow><mi>r</mi><mi>e</mi><mi>f</mi><mi>e</mi><mi>r</mi><mi>e</mi><mi>n</mi><mi>c</mi><mi>e</mi></mrow></msub></mrow></mfrac></mrow><mo>|</mo></mrow><mo>.</mo></mrow></math> Table 2 Table 2 and Figs. 1(a) and 1(b) Figs. 1(a) and 1(b) present the results of the inferred reactivity for the numerical benchmark cases. For positive step insertions of 20, 50, 100, and 200 pcm, the final reactivity was in good agreement with the induced reactivity, with the largest relative error of 1.89%, as listed in Table 2 Table 2 . A relative error of 4% was selected as the index to consider the real-time capability. The real-time reconstruction capability of the algorithm was good because it took only 1.6 s to obtain a relative error as low as 4%, as shown in Fig. 1(a) Fig. 1(a) . For a more complicated reactivity perturbation such as a sinusoidal insertion, this algorithm could effectively follow the reactivity variation, as shown in Fig. 1(b) Fig. 1(b) . 1 1 Step reactivity (pcm) 1 1 IKfinal reactivity (pcm) 1 1 Time for 4% relative error (s) 1 1 RE (%) 1 1 20 1 1 20.29 1 1 1.6 1 1 1.47 1 1 50 1 1 50.76 1 1 1.6 1 1 1.52 1 1 100 1 1 101.63 1 1 1.6 1 1 1.63 1 1 200 1 1 203.78 1 1 1.6 1 1 1.89 In the analytical benchmark cases, the following three power functions were considered: exponential power p ( t ) = e wt with w = 0.1599, ramp power p ( t ) = 1+ wt with w = 1, and sinusoidal power p ( t ) = 1+ A sin( wt ) with A = 1 and w = 0.3142. Figure 2(a) Figure 2(a) shows the calculated reactivity determined using the inverse kinetics method with respect to the reference solution, on the basis of the power function p ( t ) = e 0.1599 t . The relative error decreases with time after a maximum relative error of 4.12%, which was obtained at 0.4 s. The calculated reactivity is in good agreement with the reference solution, though the absolute error increases with time. The maximum absolute error is 6.88 pcm at t = 15 s when the reference reactivity is 300.05 pcm. Figure 2(b) Figure 2(b) shows the reactivity reconstruction for the power function p ( t ) = 1+ t . The maximum absolute error found in this case is 11.70 pcm at t = 1.3 s when the reference reactivity is 327.03 pcm. The relative error decreases with time as the maximum relative error of 4.13% is reached at t = 0.4 s. In the next example, a sinusoidal power function p ( t ) = 1+sin(0.3142 t ) is assumed, and the corresponding reactivity calculated is benchmarked against the reference solution listed in Table 1 Table 1 . The maximum absolute error found is 6.13 pcm at t = 1.9 s when the reference reactivity is 190.38 pcm. The relative error is maintained at approximately 4% for 7.5 s, as shown in Fig. 2(c) Fig. 2(c) . (Analytical solution as the reference) The comparison results of the numerical and analytical benchmark cases show that the inverse kinetics algorithm can be used in monitoring the reactivity in real time in lead-based critical reactors by following the reactivity variation with some uncertainty. The reactivity reconstructed using the inverse kinetics method is slightly overestimated compared to the reference solution. This can be improved in the future; nevertheless, it does not affect the use of the method applied to subcritical conditions. 3.2 Lead-based subcritical reactor 2 The PDS-XADS [ 21 21 ] is a lead/bismuth cooled subcritical system with a thermal power of 80 MW. The proton beam has an energy of 600 MeV with a maximum beam current of 6 mA. MOX fuel with a PU enrichment of 23% is considered. The effective multiplication factor at the beginning of the cycle is 0.97. There are 19 assembly regions in the center for accommodating the spallation target. There are 120 fuel assemblies in the subcritical core. Each fuel assembly has 90 fuel rods. The length of the flat-to-flat fuel assembly is 137.6 mm, with an assembly gap of 4 mm. There are 165 dummy assemblies encapsulating the core for reflecting the neutrons and 12 peripheral absorb assemblies for refueling. Figure 3 Figure 3 and Table 3 Table 3 present the cross plane of the PDS-XADS and its related design parameters, respectively. 1 1 Parameter 1 1 Design value 1 1 Power (MWth) 1 1 80 1 1 Fuel 1 1 MOX 1 1 Effective multiplication factor (BOC) 1 1 0.97 1 1 PU mass enrichment 1 1 23% 1 1 Effective delayed neutron fraction (pcm) 1 1 321 (8, 67, 59, 117, 53, 17) 1 1 Prompt neutron generation time (μs) 1 1 1.0142 1 1 Decay constant for precursors (in group) (s -1 ) 1 1 0.0129, 0.0313, 0.1346, 0.3443, 1.3764, 3.7425 The neutronic transient in the subcritical system can be subdivided into two categories: reactivity transient and extraneous neutron source transient. In this study, the reactivity transient was considered with a constant external neutron source, which can be determined using the initial reactivity and neutron flux or the density under the subcritical steady state, as shown in Eq. (8). In experimental nuclear reactor physics, the initial reactivity can be determined using a pulsed neutron source method or the ADS source-jerk method [ 23 23 ]. The steady state neutron flux can be easily measured using neutron monitors in the core. The super Monte Carlo program for nuclear and radiation simulation (SuperMC) [ 24 24 - 27 27 24-27 ] is used to obtain the subcritical neutron data, which serve as the reference for testing the inverse kinetics method. The initial reactivity was −2786 pcm at the beginning of the cycle. To test the reconstruction capability of the subcritical inverse kinetics method, different initial reactivities, such as −5486, −8367, and −12371 pcm, were selected by adjusting the plutonium mass enrichment. Four types of reactivity step insertions were considered: $0.1 (1 dollar is equal to the total delayed neutron fraction β ), $0.3, $0.5, and $0.9. The induced reactivity sequence comprises four phases: the first phase with a step insertion χ at 90 s, the second phase with a step insertion χ at 100 s, the third phase with a step insertion χ at 110 s, and the last phase with a linear insertion ρ = (1−0.9 t ) χ at 110 s until 135 s. Figures 4(a) and 4(b) Figures 4(a) and 4(b) show the reactivity reconstruction results for an initial reactivity of −2786 pcm, corresponding to reactivity insertion values of $0.5 and $0.9, respectively. The variation in the reactivity was detected in 1 s, thereby demonstrating the monitoring performance of the proposed method. The trend in the reactivity remains largely the same after the step perturbation; hence, the average reactivity from 90 s to 100 s was considered for the comparative analysis. Tables 4 Tables 4 and 5 5 list the test results of the reconstructed reactivity and the corresponding relative errors under different insertion values and initial reactivity states, respectively. The reactivity calculated using the inverse kinetics method was compared with the SuperMC results. 1 1 Inverse kinetics reactivity (reference reactivity) (pcm) 4 1 Initial reactivity ρ 0 (pcm) 1 1   1 1 Δ ρ = $0.1 1 1 Δ ρ =$0.3 1 1 Δ ρ = $0.5 1 1 Δ ρ = $0.9 1 1 −12371 1 1 31.96 (36.04) 1 1 93.27 (107.74) 1 1 154.61 (179.33) 1 1 277.34 (322.26) 1 1 −8367 1 1 31.73 (34.85) 1 1 93.00 (104.20) 1 1 154.31 (173.25) 1 1 277.05 (311.70) 1 1 −5486 1 1 31.72 (33.87) 1 1 92.81 (101.30) 1 1 153.89 (168.63) 1 1 276.19 (303.02) 1 1 −2786 1 1 30.80 (32.94) 1 1 91.66 (98.55) 1 1 152.61 (164.05) 1 1 274.86 (294.72) 1 1 Initial reactivity ρ 0 (pcm) 5 1 Relative error (%) 1 1   1 1 Δ ρ = $0.1 1 1 Δ ρ = $0.3 1 1 Δ ρ = $0.5 1 1 Δ ρ = $0.9 1 1 Average R. E. 1 1 −12371 1 1 11.32 1 1 13.43 1 1 13.78 1 1 13.94 1 1 13.12 1 1 −8367 1 1 8.95 1 1 10.75 1 1 10.93 1 1 11.12 1 1 10.44 1 1 −5486 1 1 6.11 1 1 8.38 1 1 8.74 1 1 8.85 1 1 8.02 1 1 −2786 1 1 6.53 1 1 6.99 1 1 6.97 1 1 6.74 1 1 6.81 4 Discussion 1 The reactivity reconstruction results obtained using the inverse kinetics method under different initial reactivity states show that this method is capable of real-time reactivity reconstruction of both critical and subcritical reactors. In the lead-based critical reactor, the reactivity reconstructed using the inverse kinetics method was compared with those reconstructed using the numerical and analytical methods. In terms of accuracy, the largest relative error was below 1.89% for an insertion of 200 pcm in the numerical benchmark cases. In the analytical benchmark cases, the largest relative error was approximately 4%. This demonstrated that the method was effective in reconstructing the reactivity of the lead-based critical reactor. For the lead-based subcritical reactor, as in the case of the PDS-XADS, the reactivity reconstruction capability is not as high as that in the critical reactor. The absolute error depends on the initial reactivity and the corresponding induced reactivity. For a given initial reactivity of −2786 pcm, the absolute errors were 2.14, 6.89, 11.44, and 19.86 pcm, corresponding to reactivity perturbations of $0.1, $0.3, $0.5, and $0.9, respectively. The absolute error tends to increase for stronger reactivity perturbations. However, the relative errors were 6.53, 6.99, 6.97, and 6.74%, which are approximately the same. Hence, the average relative error depends on the initial reactivity. The same trend was observed for other initial reactivities. For the same positive reactivity insertion, such as $0.5, the absolute errors were 11.44, 14.74, 18.94, and 24.72 pcm, corresponding to initial reactivities of −2786, −5486, −8367, and −12371 pcm, respectively. The relative errors tend to increase as follows: 6.97, 8.74, 10.93, and 13.78%. For any other given reactivity perturbation, the same trend was observed. This shows that both the absolute and relative errors increase with the decrease in the initial reactivity. If the relative error is controlled to be below 13.12%, the initial reactivity obtained using the inverse kinetics method can reach up to −12371 pcm. From all the cases pertaining to the subcritical reactor, the reactivity obtained using the inverse kinetics method was globally underestimated. The inverse kinetics method is more reliable and robust for reconstructing the reactivity of critical reactors. In source-driven subcritical systems, even if the beam current remains unchanged, the effective source would change with respect to the reactivity perturbation for a given initial reactivity. This scenario indicates that the assumption of a constant external neutron source is questionable. Moreover, the inverse kinetics method is based on a simplified point-kinetics model, wherein the variation in the shape function of the neutron flux and the complex neutron energy spectrum were neglected. Further, recent advancements show that the classical point-kinetics method has limitations in the application of the neutron diffusion theory [ 28 28 , 29 29 ]. For example, the sub-diffusive effects are neglected in this method [ 30 30 ]. Nevertheless, this can be explained using a time-fractional point kinetics model [ 31 31 , 32 32 ], which has some applications [ 33 33 - 35 35 33-35 ]. In the future, we will employ the inverse equation in the fractional point-neutron kinetics model for subcritical reactivity reconstruction. 5 Conclusion 1 In this study, a preliminary analysis was performed on the inverse kinetics method used for reactivity reconstruction under various initial reactivity states. The results show that the reactivity of a lead-based critical reactor was effectively reconstructed using the inverse kinetics method. However, for an initial reactivity below zero, as observed in lead-based subcritical reactors, the reactivity reconstructed using this method is globally underestimated. For a selected reactivity perturbation, both the relative and absolute errors increase with the decrease in the initial reactivity. For a particular initial reactivity, the absolute error increases with the increase in the reactivity perturbation, whereas the relative error remains the same. This deviation is due to the variation in the external neutron source, spatial-spectral effects, and sub-diffusive effects, which need to be investigated in the future. References W. Johanna , Dissertation , Chalmers University of Technology , 2005 . Dissertation A. Alessandro , S. Valery , F.S. Craig et al ., Overview of lead-cooled fast reactor activities . Prog. Nucl. Energy 77 , 300 - 307 ( 2014 ). doi: 10.1016/j.pnucene.2013.11.011 http://doi.org/10.1016/j.pnucene.2013.11.011 Overview of lead-cooled fast reactor activities Y.C. Wu , M.H. Wang , Q.Y. Huang et al ., Development status and prospects of lead-based reactors . Nucl. Sci. Eng . 35 ( 2 ), 213 - 221 ( 2015 ). doi: 10.3969/j.issn.0258-0918.2015.02.004 http://doi.org/10.3969/j.issn.0258-0918.2015.02.004 . (In Chinese) Development status and prospects of lead-based reactors Y.C. Wu , Y.Q. Bai , Y. Song et al ., Development strategy and conceptual design of China lead-based research reactor . Ann. Nucl. Energy 87 , 511 - 516 ( 2016 ). doi: 10.1016/j.anucene.2015.08.015 http://doi.org/10.1016/j.anucene.2015.08.015 Development strategy and conceptual design of China lead-based research reactor Y.C. Wu . Design and R&D progress of China lead-based reactor for ADS research facility . Eng . 2 ( 1 ), 124 - 131 ( 2016 ). doi: 10.1016/j.eng.2016.01.023 http://doi.org/10.1016/j.eng.2016.01.023 Design and R&D progress of China lead-based reactor for ADS research facility Y.C. Wu , J.P. Qian , J.N. Yu , The fusion-driven hybrid system and its material selection . J. Nucl. Mater. S 307-311 ( 2 ), 1629 - 1636 ( 2002 ). doi: 10.1016/S0022-3115(02)01272-2 http://doi.org/10.1016/S0022-3115(02)01272-2 The fusion-driven hybrid system and its material selection Y.C. Wu , FDS Team , Conceptual design activities of FDS series fusion power plants in China . Fusion. Eng. Des . 81 ( 23-24 ), 2713 - 2718 ( 2006 ). doi: 10.1016/j.fusengdes.2006.07.068 http://doi.org/10.1016/j.fusengdes.2006.07.068 Conceptual design activities of FDS series fusion power plants in China Y.C. Wu , FDS Team , Conceptual design of the China fusion power plant FDS-II . Fusion. Eng. Des . 83 ( 10 ), 1683 - 1689 ( 2008 ). doi: 10.1016/j.fusengdes.2008.06.048 http://doi.org/10.1016/j.fusengdes.2008.06.048 Conceptual design of the China fusion power plant FDS-II Y.C. Wu , J.Q. Jiang , M.H. Wang et al ., A fusion-driven subcritical system concept based on viable technologies . Nucl. Fusion 51 ( 10 ), 532 - 542 ( 2011 ). doi: 10.1088/0029-5515/51/10/103036 http://doi.org/10.1088/0029-5515/51/10/103036 A fusion-driven subcritical system concept based on viable technologies L.J. Qiu , Y.C. Wu , B.J. Xiao et al ., A low aspect ratio tokamak transmutation system . Nucl. Fusion 40 , 629 - 633 ( 2000 ). doi: 10.1088/0029/5515/40/3Y/325 http://doi.org/10.1088/0029/5515/40/3Y/325 A low aspect ratio tokamak transmutation system Y.C. Wu , FDS Team , Design status and development strategy of China liquid lithium-lead blankets and related material technology . J. Nucl. Mater . 367-370 ( 4 ), 1410 - 1415 ( 2007 ). doi: 10.1016/j.jnucmat.2007.04.031 http://doi.org/10.1016/j.jnucmat.2007.04.031 Design status and development strategy of China liquid lithium-lead blankets and related material technology Y.C. Wu , FDS Team , Design analysis of the China dual-functional lithium lead (DFLL) test blanket module in ITER . Fusion Eng. Des . 82 ( 15-24 ), 1893 - 1903 ( 2007 ). doi: 10.1016/j.fusengdes.2007.08.012 http://doi.org/10.1016/j.fusengdes.2007.08.012 Design analysis of the China dual-functional lithium lead (DFLL) test blanket module in ITER Y.C. Wu , FDS Team , Conceptual design and testing strategy of a dual functional lithium-lead test blanket module in ITER and EAST . Nucl. Fusion 47 ( 11 ), 1533 - 1539 ( 2007 ). doi: 10.1088/0029-5515/47/11/015 http://doi.org/10.1088/0029-5515/47/11/015 Conceptual design and testing strategy of a dual functional lithium-lead test blanket module in ITER and EAST H. Khalafi , S.H. Mosavi , S.M. Mirvakili , Design & construction of a digital real time reactivity meter for Tehran research reactor . Prog. Nucl. Energy 53 ( 1 ), 100 - 105 ( 2011 ). doi: 10.1016/j.pnucene.2010.08.001 http://doi.org/10.1016/j.pnucene.2010.08.001 Design & construction of a digital real time reactivity meter for Tehran research reactor S. Dulla , M. Nervo , P. Ravetto , A method for the continuous monitoring of reactivity in subcritical source-driven systems . Ann. Nucl. Energy 87 , 1 - 11 ( 2016 ). doi: 10.1016/j.anucene.2014.12.001 http://doi.org/10.1016/j.anucene.2014.12.001 A method for the continuous monitoring of reactivity in subcritical source-driven systems K.O. Ott , R.J. Neuhold , Introductory nuclear reactor dynamics . ( American Nuclear Society , La Grange Park, Illinois , 1985 ), pp. 212 - 217 . Introductory nuclear reactor dynamics H. Malmir , N. Vosoughi , On-line reactivity calculation using Lagrange method . Ann. Nucl. Energy 62 ( 12 ), 463 - 467 ( 2013 ). doi: 10.1016/j.anucene.2013.07.006 http://doi.org/10.1016/j.anucene.2013.07.006 On-line reactivity calculation using Lagrange method D.S. Diaz , J.F.F. Ospina , J.A.R. Sarasty , Hamming method for solving the delayed precursor concentration for reactivity calculation . Ann. Nucl. Energy 42 , 47 - 49 ( 2012 ). doi: 10.1016/j.anucene.2011.12.019 http://doi.org/10.1016/j.anucene.2011.12.019 Hamming method for solving the delayed precursor concentration for reactivity calculation D.S. Diaz , A.S. Martinez , F.C. Da Silva , Formulation for the calculation of reactivity without nuclear power history . J. Nucl. Sci. Technol . 44 ( 9 ), 1149 - 1155 ( 2007 ). doi: 10.3327/Jnst.44.1149 http://doi.org/10.3327/Jnst.44.1149 Formulation for the calculation of reactivity without nuclear power history Z.X. Gu , G. Wang , Z. Wang et al ., Transient analyses on loss of heat sink and overpower transient of natural circulation LBE-cooled fast reactor . Prog. Nucl. Energy 81 , 60 - 66 ( 2015 ). doi: 10.1016/j.pnucene.2015.01.010 http://doi.org/10.1016/j.pnucene.2015.01.010 Transient analyses on loss of heat sink and overpower transient of natural circulation LBE-cooled fast reactor Ansaldo Nucleare , Core configuration specification of the LBE-Cooled XADS (2001) , http://nuklear-server.nuklear.kit.edu/ADOPT/ADOPT_Transfer_FZK/Files/deliverables/DEL03010.pdf http://nuklear-server.nuklear.kit.edu/ADOPT/ADOPT_Transfer_FZK/Files/deliverables/DEL03010.pdf . Accessed 28 Feb 2017 Core configuration specification of the LBE-Cooled XADS (2001) H.Q. Yuan , D.P. Hu , High-order end floating method for solving point reactor neutron kinetics equations . Nucl. Power. Eng . 16 ( 2 ), 24 - 128 ( 1995 ). (In Chinese) High-order end floating method for solving point reactor neutron kinetics equations F. Mellier , The MUSE experiments for subcritical neutronics validation, FIKW-CT-2000-00063 , 2005 . The MUSE experiments for subcritical neutronics validation, FIKW-CT-2000-00063 Y.C. Wu , FDS Team , CAD-based interface programs for fusion neutron transport simulation . Fusion. Eng. Des . 84 ( 7-11 ), 1987 - 1992 ( 2009 ). doi: 10.1016/j.fusengdes.2008.12.041 http://doi.org/10.1016/j.fusengdes.2008.12.041 CAD-based interface programs for fusion neutron transport simulation Y.C. Wu , J. Song , H.Q. Zheng , et al ., CAD-based monte carlo program for integrated simulation of nuclear system SuperMC . Ann. Nucl. Energy . 82 , 161 - 168 ( 2014 ). doi: 10.1016/j.anucene.2014.08.058 http://doi.org/10.1016/j.anucene.2014.08.058 CAD-based monte carlo program for integrated simulation of nuclear system SuperMC Y.C. Wu , Z.S. Xie , U. Fischer , A discrete ordinates nodal method for one-dimensional neutron transport calculation in curvilinear geometries . Nucl. Sci. Eng . 133 ( 3 ), 350 - 357 ( 1999 ). A discrete ordinates nodal method for one-dimensional neutron transport calculation in curvilinear geometries Y.C. Wu , Z.B. Chen , L.Q. Hu et al ., Identification of safety gaps for fusion demonstration reactors . Nature Energy ( 2016 ). doi: 10.1038/NENERGY.2016.154 http://doi.org/10.1038/NENERGY.2016.154 Identification of safety gaps for fusion demonstration reactors G. Espinosa-Paredes , J.B. Morales-Sandoval , R. Vázquez-Rodríguez et al ., Constitutive laws for the neutron density current . Ann. Nucl. Energy . 35 ( 10 ), 1963 - 1967 ( 2008 ). doi: 10.1016/j.anucene.2008.05.002 http://doi.org/10.1016/j.anucene.2008.05.002 Constitutive laws for the neutron density current G. Espinosa-Paredes , R. Vázquez-Rodríguez , E. del Valle Gallegos et al ., Fractional-space law for the neutron current density . Ann. Nucl. Energy . 55 ( 55 ), 120 - 125 ( 2013 ). doi: 10.1016/j.anucene.2016.08.007 http://doi.org/10.1016/j.anucene.2016.08.007 Fractional-space law for the neutron current density G. Espinosa-Paredes , M.A. Polo-Labarrios , J. Alvarez-Ramirez , Anomalous diffusion processes in nuclear reactors . Ann. Nucl. Energy . 54 , 227 - 232 ( 2013 ). doi: 10.1016/j.anucene.2012.11.024 http://doi.org/10.1016/j.anucene.2012.11.024 Anomalous diffusion processes in nuclear reactors G. Espinosa-Paredes , M.A. Polo-Labarrios , E.G. Espinosa-Martinez et al ., Fractional neutron point kinetics equations for nuclear reactor dynamics . Ann. Nucl. Energy . 38 ( 2 ), 307 - 330 ( 2011 ). doi: 10.1016/j.anucene.2010.10.012 http://doi.org/10.1016/j.anucene.2010.10.012 Fractional neutron point kinetics equations for nuclear reactor dynamics G. Espinosa-Paredes , Fractional-space neutron point kinetics (F-SNPK) equations for nuclear reactor dynamics . Ann. Nucl. Energy . ( 2016 ). doi: 10.1016/j.anucene.2016.08.007 http://doi.org/10.1016/j.anucene.2016.08.007 Fractional-space neutron point kinetics (F-SNPK) equations for nuclear reactor dynamics M.A. Polo-Labarrios , G. Espinosa-Paredes , Numerical analysis of startup PWR with fractional neutron point kinetic equation . Prog. Nucl. Energy . 60 ( 2071 ), 38 - 46 ( 2012 ). doi: 10.1016/j.pnucene.2012.05.003 http://doi.org/10.1016/j.pnucene.2012.05.003 Numerical analysis of startup PWR with fractional neutron point kinetic equation T.K. Nowak , K. Duzinkiewicz , R. Piotrowski , Numerical solution analysis of fractional point kinetics and heat exchange in nuclear reactor . Nucl. Eng. Des . 281 , 121 - 130 ( 2015 ). doi: 10.1016/j.nucengdes.2014.11.028 http://doi.org/10.1016/j.nucengdes.2014.11.028 Numerical solution analysis of fractional point kinetics and heat exchange in nuclear reactor S.S. Ray , A. Patra , On the solution of the nonlinear fractional neutron point-kinetics equation with newtonian temperature feedback reactivity . Nucl. Technol . 189 ( 1 ), 103 - 109 ( 2015 ). doi: 10.13182/NT13-148 http://doi.org/10.13182/NT13-148 On the solution of the nonlinear fractional neutron point-kinetics equation with newtonian temperature feedback reactivity The authors would also like to thank the other members of the FDS team. 10.1007/s41365-018-0367-3 This work was supported by the Strategic Priority Science & Technology Program of the Chinese Academy of Sciences (No. XDA03040000), the National Natural Science Foundation of China (NSFC) (Nos. 11305205, 11305203, and 11405204), the Special Program for Informatization of the Chinese Academy of Sciences (No. XXH12504-1-09), the Anhui Provincial Special project for High Technology Industry, and the Special Project of Youth Innovation Promotion Association of Chinese Academy of Sciences, the Industrialization Fund. 2017-01-16 2017-03-08 2017-03-19 2018-02-05 2018-02-01 2021-04-11 © China Science Publishing & Media Ltd. (Science Press), Shanghai Institute of Applied Physics, the Chinese Academy of Sciences, ChineseNuclear Society and Springer Nature Singapore Pte Ltd. 2018 This is a post-peer-review, pre-copyedit version of an article published in Nuclear Science and Techniques. The final authenticated version is available online at: http://dx.doi.org/10.1007/s41365-018-0367-3 http://dx.doi.org/10.1007/s41365-018-0367-3 . 2018 Nuclear Science and Techniques 02 29