- Email: [email protected]
The axial power distribution validation of the SCWR fuel assembly with coupled neutronics-thermal hydraulics method
Nuclear Engineering and Design 258 (2013) 157–163
Contents lists available at SciVerse ScienceDirect
Nuclear Engineering and Design journal homepage: www.elsevier.com/locate/nucengdes
The axial power distribution validation of the SCWR fuel assembly with coupled neutronics-thermal hydraulics method Xi Xi, Zejun Xiao ∗ , Xiao Yan, Yongliang Li, Yanping Huang CNNC Key Laboratory on Nuclear Reactor Thermal Hydraulics Technology, Nuclear Power Institute of China, Chengdu 610041, China
h i g h l i g h t s CFX and MCNP codes are suitable to calculate the axial power profile of the FA. The partition method in the calculation will affect the final result. The density feedback has little effect on the axial power profile of CSR1000 FA.
a r t i c l e
i n f o
Article history: Received 6 September 2012 Received in revised form 13 December 2012 Accepted 8 January 2013
a b s t r a c t SCWR (super critical water reactor) is one of the IV generation nuclear reactors in the world. In a typical SCWR the water enters the reactor from the cold leg with a temperature of 280 ◦ C and then leaves the core with a temperature of 500 ◦ C. Due to the sharp change in temperature, there is a huge density change of the water along the axial direction of the fuel assembly (FA), which will affect the moderating power of the water. So the axial power distribution of the SCWR FA could be different from the traditional PWR FA.In this paper, it is the first time that the thermal hydraulics code CFX and neutronics code MCNP are used to analyze the axial power distribution of the SCWR FA. First, the factors in the coupled method which could affect the result are analyzed such as the initialization value or the partition method especially in the MCNP code. Then the axial power distribution of the Europe HPLWR FA is obtained by the coupled method with the two codes and the result is compared with that obtained by Waata and Reiss. There is a good agreement among the three kinds of results. At last, this method is used to calculate the axial power distribution of the Chinese SCWR (CSR1000) FA. It is found the axial power profile of the CSR1000 FA is not so sensitive to the change of the moderator density. © 2013 Elsevier B.V. All rights reserved.
1. Introduction SCWR is one kind of the IV generation nuclear reactors which was recommended by GIF in 2002. Due to the high thermal efficiency, compact system as well as the abundant operation and manufacture experience of the fossil supercritical power plant, SCWR is the most promising one in the future (US DOE, 2002). But before the commercial deployment of the SCWR, a lot of key techniques should be solved, one of which is the neutronics-T/H characteristics of the FA. In a typical SCWR, when the coolant flows into the core through the cold legs, it will be divided into three parts, one part flows towards the downcomer and the other parts flow towards the upper plenum and then into the water rod and the gap between the neighboring FAs as the moderator. Then the three different water will mix in the lower plenum and flow together into
∗ Corresponding author. Tel.: +86 028 85906199; fax: +86 028 85908889. E-mail address: fabulous [email protected] (Z. Xiao). 0029-5493/$ – see front matter © 2013 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.nucengdes.2013.01.031
the core again as the coolant shown in Fig. 1 (Bitterman et al., 2003). During this process, not only the density of the coolant changes a lot along the axial direction but also the density of the moderator will change simultaneously, which could result in the heterogeneous moderating power distribution along the axial direction of the FA. So the axial power distribution of the SCWR FA could be different from the traditional PWR FA which follows a cosine shape axial power profile. In 2005, a coupled neutronics-thermal hydraulics calculation was carried out by Waata to analyze the axial power distribution of the Europe HPLWR FA as well as the factors which could affect it such as the burn-up, density feedback effect, Doppler effect, etc. (Waata et al., 2005). In the calculation, neutronics code MCNP and subchannel T-H code STAFAS are used and two inactive parts of the FA are also modeled. The result shows that the axial power profile of the HPLWR FA does not follow the cosine shape when the burn-up is not taken into account. In fact there are two peaks, the larger one is close to the lower part and the smaller one is close to the upper part of the FA. The asymmetric power profile is mainly caused by the density feedback effect and has little relationship with the Doppler
158
X. Xi et al. / Nuclear Engineering and Design 258 (2013) 157–163
Nomenclature Cp E H N ˙ m q Q S ST Su Sv Sw T u u
v w
constant pressure specific heat (J/K kg) fission energy of the FA specific enthalpy (J/kg) partition number of the FA in the axial direction Mass flow rate (kg/s) average heat flux of the FA power of the FA surface of the fuel rods source term in energy equation (kg K/m3 s) source term in x momentum equations (Pa/m) source term in y momentum equations (Pa/m) source term in z momentum equations (Pa/m) temperature (K) vector of the velocity (m/s) velocity in x direction (m/s) velocity in y direction (m/s) velocity in z direction (m/s)
Greek letters source term in thermal conductivity equation (W/m3 ) thermal conductivity (W/m K) viscosity (Pa s) density of the fluent (kg/m3 ) Subscript C coolant region D downcomer G assembly gap region i the ith section of the FA in inlet out outlet t total W water rod region
although the axial power profile obtained by Reiss agrees well with Waata’s, there is only one power peak which is close to the bottom of the FA. Furthermore, both their results indicate that not only the burn-up but also the enrichment has a great influence on the power profile of the Europe HLWR FA. In 2009, the T/H code ATHAS and MCNP code were coupled to analyze the T/H characteristics as well as the axial power distribution of the SCWR-CANDU FA by Shan et al. (2008). The result shows that the density change of the water has little effect on the axial power distribution of the FA. The power profile follows the cosine shape approximately. Except for the research mentioned above, there are also some other coupled neutronics-thermal hydraulics analysis of the SCWR, which focuses on the power distribution of the whole core but they will not be discussed in detail here anymore (Cao et al., 2009; Hu et al., 2009). In this paper, the CFX 13.0 code and MCNP4C code are coupled for the axial power distribution analysis of the SCWR FA. It should be mentioned here that to avoid hot spots caused by a nonuniform power profile or by uncertainties and allowances for operation, a three passes SCWR is designed by Schulenberg and Starflinger (2007) to replace the single pass SCWR recommended by Bitterman et al. (2003). But the axial power profile of the single pass Europe SCWR is still selected for the code validation. In the calculation, only the density change effect of the water is taken into account and the Doppler effect is neglected. The following factors such as initialization value for the calculation and zone partition method in the MCNP code which could affect the final coupled result are discussed in detail. The burn-up and enrichment of the FA are not taken into consideration here. Finally, the axial power distribution of the CSR1000 FA is obtained with the same method. But it should be noted here that the CSR1000 is a two pass core design and for each pass the T/H boundary conditions such as mass flow rate, inlet temperature are different, so two types of axial power distributions will be attained finally.
2. Geometrical models and boundary conditions 2.1. Geometrical models
effect. In 2008, another coupled neutronics-thermal hydraulics analysis was carried out by Reiss with a self-developed T/H code and MCNP code (Reiss et al., 2008). In Reiss’s model, only the active part of the FA is taken into consideration and the result shows that
Both the Europe SCWR and CSR1000 FA are 1/8th symmetrical, so in the coupled numerical analysis only 1/8th of the FA is taken into account. Table 1 shows the main geometrical parameters of the
Fig. 1. Water route of the Europe and Chinese SCWR.
X. Xi et al. / Nuclear Engineering and Design 258 (2013) 157–163
159
Table 1 Parameters of the fuel assemblies. Parameters
Active height Fuel rod diameter Cladding thickness Fuel pin diameter Pitch/diameter ratio Moderator box lengthen Gap between FA Gap between rod and box wall Moderator box wall thickness FA box thickness Gaseous core diameter
Version Europe
Chinese
4200 8 0.5 6.7 1.150 26 5 1 0.4 1 –
4200 9.5 0.57 8.19 1.105 50.9 10 1 1 2 1.5
Europe and Chinese FAs. The Europe SCWR is the same as which analyzed by Waata and Reiss. Although the configurations of the two fuel assemblies are nearly the same, the CSR1000 FA is larger than the Europe FA, which has 56 fuel rods. Furthermore, there is a gaseous core in the center of the CSR1000 fuel pin, filled with the helium. And the control rod is not taken into account here. To find how the partition method in the MCNP code could affect the final axial power profile of the SCWR FA, four kinds of partition methods are used in this paper, which are P1, P2, P3 and P4 respectively. Figs. 2–4 show the zone partition methods of the cross-section of the Europe FA, which correspond to P1, P2 and P3, the partition method of P4 is similar to the P2, the only difference is that the FA is divided into 21 zones axially in P1–P3, while in P4, there are 42 axial zones. Fig. 5 shows the partition method of the CSR1000 FA which will be used for the coupled analysis. It can be found that the water rod, coolant region and the gap between the FAs are not further divided for the MCNP code as shown in Fig. 2. While in Fig. 3, the coolant region is divided into 7 zones. Moreover, not only the water rod but also the FA gap is further divided. Considering that the gradient of the water density is vertical to the water rod box and FA box, the water in these two regions are divided into three parts from the inner to the outer respectively. In Fig. 4, the partition method of the coolant region is similar to Waata’s, which is a 9 subzones division and the partition method of the water rod region and gap region is the same as that in Fig. 3. With the change of the partition method for the MCNP code, the mesh type for the CFX code should be changed
Fig. 2. Partition method of the Europe FA (P1) 101 assembly gap, 102 coolant, 103 water rod, 201–207 fuel rods, 301–303 cladding and box walls, 401 helium.
Fig. 3. Partition method of the Europe FA (P2) 101–103 water rod, 104–110 coolant, 111–113 assembly gap, 201–207 fuel rods, 301–303 cladding and box walls, 401 helium.
simultaneously. For example, the mesh type of the coolant region in P3 is triangular prism, whereas it is hexahedral in P1 and P2. But the mesh types of the water rod and gap regions are the same from P2 to P4, which are hexahedral. The mesh numbers for the Europe fuel assembly are 3,200,000 and 3,000,000 respectively with triangular prism and hexahedral mesh. It is 3,500,000 for Chinese fuel assembly with triangular prism mesh. 2.2. Boundary conditions The boundary conditions mentioned here not only include the temperature, mass flow rate, etc., but also the material properties such as the density and the enrichment of the fuel pin. Table 2 shows the thermal hydraulics boundary conditions for CFX code and Table 3 is the neutronics boundary conditions for MCNP code. It should be noted here again that the CSR1000 is a typical two flow pass reactor. The symbol (*) in Table 2 means that the inlet temperature of the coolant could not be input directly because it is related to the outlet
Fig. 4. Partition method of the Europe FA (P3) 101–103 water rod, 104–1120 coolant, 113–115 assembly gap, 201–207 fuel rods, 301–303 cladding and box walls, 401 helium.
160
X. Xi et al. / Nuclear Engineering and Design 258 (2013) 157–163
and energy conservation equations but also the 3D heat conductivity equation should be solved which are as shown below. But it should be mentioned here that the conductivity in the fuel rod is not taken into account here because the Doppler effect of the fuel rod is neglected.
∇ · (u) = 0
(2)
∂p ∇ · (uu ) = ∇ · (∇ u) − + Su ∂x
(3)
∂p + Sv ∂y
(4)
∇ · (vu ) = ∇ · (∇ v) −
∇ · (wu ) = ∇ · (∇ w) −
∇ · (T u ) = ∇ · c Fig. 5. Partition method of the CSR1000 FA 101 assembly gap, 102 coolant, 103 water rod, 201–202 fuel rods, 301–303 cladding and box walls, 401 helium. Table 2 Boundary conditions for T/H analysis of 1/8th FA. Parameters
System pressure Water properties Tin of the water rod Tin of the assembly gap Tin of the coolant region Min of the water rod Min of the assembly gap Min of the coolant region FA power
Version Europe
Chinese I
Chinese II
25 MPa IAPWS IF97 553 K 553 K * 0.0139 kg/s 0.0278 kg/s 0.167 kg/s 327.5 kW
25 MPa IAPWS IF97 553 K 553 K 553 K 0.0389 kg/s 0.0706 kg/s 0.2549 kg/s 199 kW
25 MPa IAPWS IF97 553k 553 K 657 K 0.0517 kg/s 0.0706 kg/s 0.3697 kg/s 199 kW
temperature of the water rod and the assembly gap. The following correlation can be used to obtain HinC in the coolant region: HinC =
˙ inD · HinD + m ˙ inW · HoutW + m ˙ inG · HoutG m ˙t m
(1)
∂2 T ∂2 T ∂2 T + + 2 2 ∂x ∂y ∂z 2
(5)
+ ST
(6)
=
Cp
(7)
The difference scheme of the convection terms in the conservation equations is second order upwind scheme, and the SST model is selected as the turbulence model in the calculation, which is suitable for the pseudocritical conditions. When the residual is less than 10−3 , the T/H solution is thought to be converged and the water density distribution in the FA will be obtained. The neutron transport equation in the FA will be solved by MCNP which is a general-purpose, continuous-energy, generalizedgeometry, time-dependent, coupled neutron/photon/electron Monte Carlo transport code. It can be used in several transport modes: neutron only, photon only, electron only, combined neutron/photon transport where the photons are produced by neutron interactions, neutron/photon/electron, photon/electron, or electron/photon (LANL, 2003). With the MCNP code, not only the keff but also the power distribution of the FA will be obtained. The power of the fuel rod here will be calculated with F7 tally. The particle number in one cycle is 50,000, and inactive cycle is 500, and total cycle is 1000. With the MCNP calculation, the mass averaged fission energy can be obtained for each section of the FA, and the average heat flux of the fuel rods in this section can be attained with the following equation: qi =
3. Physical models and coupled method
∇T Cp
∂p + Sw ∂z
Ei · Q
N E i=1 i
(8) · Si
3.1. Physical models 3.2. Coupled method The physical models here refer to the T/H model and neutronics model in the coupled calculations. The T/H phenomenon in the SCWR FA includes the heat conductivity in the water box wall and FA wall, the convection heat transfer in the water rod, coolant, assembly gap, etc. So not only the steady state mass, momentum Table 3 Boundary conditions for neutronics analysis of 1/8th FA. Parameters
Corner fuel rod enrichment Other fuel rod enrichment Fuel rod density Material of the cladding Material of the water box Material of the assembly box Density of the cladding water box and assembly box Material of the gas Cross section library
Version Europe
Chinese
4% 5% 10,600 kg/m3 Alloy 316 Alloy 316 Alloy 316 7450 kg/m3
4.3% 5.7% 10,600 kg/m3 Alloy 310S Alloy 310S Alloy 310S 7450 kg/m3
Helium ENDF6
Helium ENDF6
The external coupled method is used here and the flowchart is shown in Fig. 6. If the T/H calculation will be carried out first, an assumed power profile will be given as the boundary condition for the CFX code and the water density distribution will be obtained when the solution is converged. Then the density will be input to the MCNP code to obtain the new power profile. If the relative error of the new profile and the old one is small enough, the coupled numerical simulation will be stopped and it is assumed that a converged solution is attained. The calculation is run in a parallel computer which has 24 CPUs, so for one coupled iteration it will take only 30 min. To decrease the iteration steps and get a converged result, a relaxation factor 0.8 is used in the coupled calculation. 4. Results and analysis 4.1. Influence of the partition method As a numerical simulation approach, not only the mesh in the CFX code but also the partition method in the MCNP code could
X. Xi et al. / Nuclear Engineering and Design 258 (2013) 157–163
161
Fig. 6. Flowchart of the coupled method. Fig. 8. Comparison of the average axial heat flux with different partition methods in MCNP code.
affect the precision and convergence process of the solution. The mesh influence of the CFD code has already been discussed by so many researchers when the T/H analysis of the FA was carried out and there will be no further discussion here on this (Gu et al., 2008; Laltu Chandra et al., 2009). But the influence of the partition method of the MCNP code is not well analyzed especially for the coupled neutronics-T/H calculation, so how it will affect the coupled solution is discussed in detail here and four different partition methods are selected which correspond to P1–P4. During the coupled numerical simulation, it is found that the solution is too difficult to converge with the partition method P1. After 10 iterations, there is still a periodic oscillation of the axial power profile which is shown in Fig. 7. To make the comparison with different partition methods possible, the result of P1 is the algebraic average value from iteration 10 to iteration 15. The comparison in Fig. 8 shows that although different partition methods are used in the coupled calculation, there is a good agreement among these four axial heat flux profiles. With the increase of the partition zones in the transverse coolant region, the peak value of heat flux will decrease slightly. On the other hand, with the increase of the partition zones in the axial region, there is a 4% increase of the peak value. So it can be concluded that the partition method in the MCNP code will affect the average axial heat flux profile of the FA. But compared to the axial partition method, the transverse partition method plays a minor influence on the axial power distribution.
Fig. 7. Average axial heat flux of different iteration steps with partition method P1.
4.2. Influence of the initialization value In a numerical simulation, the initialization value will probably affect not only the convergence process but also the final result. To analyze the influence of the initialization value in the coupled calculation, two different heat flux profiles are used as initial condition in the CFX code, which are homogeneous profile and standard cosine profile respectively. The partition method in the calculation is P3 and the obtained final coupled axial heat flux profile is shown in Fig. 9. The comparison shows that the average heat flux profiles obtained with different initialization values are nearly the same, which means it has no effect on the final result. But when the homogeneous heat flux is selected as the initialization value, it costs 8 iteration steps for the solution to be converged, while the cosine shape heat flux is selected, 11 iteration steps are needed to obtain a converged solution. 4.3. Results comparison To validate whether the coupled method in this paper is appropriate for the coupled analysis, the axial power distribution of the Europe SCWR FA obtained with CFX and MCNP codes is compared with that obtained by Waata and Reiss. The partition method here is P3 and Fig. 10 shows the comparison of the average axial heat flux distributions of the Europe FA.
Fig. 9. Comparison of the average axial heat flux with different initialization values in MCNP code.
162
X. Xi et al. / Nuclear Engineering and Design 258 (2013) 157–163
Fig. 10. Comparison of the average axial heat flux.
Fig. 11. Flow directions in two passes of Chinese FA.
From the comparison it can be found that the average axial heat flux profile obtained by Reiss and Xi is nearly the same. There is a power peak at the close of the coolant inlet and then decreases gradually. The maximum heat flux in Xi’s result is larger than Reiss’s with a relative error of 9%. But the heat flux profile obtained by Waata is a little different from Xi’s, which has two power peaks. This is thought to be caused by the two inactive parts in Waata’s model, especially the one close to the coolant outlet where the water density in the water rod and FA gap is larger, which results in a stronger moderating power. So the upper part heat flux of the FA in Watta’s result is larger than Reiss’s and Xi’s. The result comparison demonstrates that CFX is an appropriate code which could be used for the coupled netronics-T/H analysis of the FA. And the results obtained by CFX code and subchannel code are so close to each other if the same geometrical model are used. 4.4. Axial power distribution of the CSR1000 FA Fig. 12. Average axial heat flux of Chinese FA.
The CSR1000 is the Chinese SCWR designed by Nuclear Power Institute of China. It is a typical two pass SCWR with an inlet temperature of 280 ◦ C and an outlet temperature of 500 ◦ C. Due to the special geometrical configuration of the CSR1000 FA as well as the two flow passes, the axial power distribution of the CSR1000 could be different from the Europe SCWR FA. So the CFX and MCNP codes are coupled here to calculate the axial power profile of the CSR1000 FA. It has been demonstrated in this paper that the transverse partition method in the MCNP code has minor effect on the average axial power distribution. When the coupled neutronics-T/H analysis of the CSR1000 FA is carried out, the simplest transverse partition method is used, which is shown in Fig. 5. The flow direction in two passes and the two types average heat flux profiles of the FA are shown in Figs. 11 and 12. Figs. 13 and 14 show the volume average water density and temperature in pass 1 and pass 2 respectively. It is found that the changes of the water density and temperature along the FA are totally different in the two passes. For example in pass 1, the water density increases with the increase of the axial height, while in pass 2, it will decrease first and then increase to 600 kg/m3 . Furthermore, the change of the water density and temperature in pass 1 is much larger than that in pass 2. But in Fig. 12, it is found that the average heat flux profiles of FA in different passes are nearly the same, which is a cosine profile. So it is assumed that the power distribution is not sensitive to the change of the water density. To validate the assumption mentioned above and find out the influence of the density feedback effect on the power distribution, another neutronics analysis is carried out with MCNP code. In this
calculation, there is no axial zone partition for the FA and the average density is used. The density of the water in the coolant region is supposed to be 300 kg/m3 and will not change anymore, then the density of the water in the water rod and assembly gap will decrease from 800 kg/m3 to 300 kg/m3 and for each 100 kg/m3 a calculation will be run to obtain the keff . Table 4 shows the change of the keff with the decrease of the moderator density.
Fig. 13. Volume average density of the water in two passes of Chinese FA.
X. Xi et al. / Nuclear Engineering and Design 258 (2013) 157–163
163
Table 4 keff of the FA with different moderator densities. Water rod density (kg/m3 )
Coolant density (kg/m3 )
Assembly gap density (kg/m3 )
keff of the Europe
keff of the CSR1000
800 700 600 500 400 300
300 300 300 300 300 300
800 700 600 500 400 300
1.15960 1.14801 1.13400 1.11530 1.09188 1.06047
1.27265 1.27350 1.27068 1.26054 1.23942 1.20447
appropriate initialization is applied, the solution will converge more quickly. (2) CFX is an appropriate T/H code which could be used for the coupled neutronics-T/H analysis of the SCWR FA. The average axial heat flux profile of the Europe FA obtained with CFX code and MCNP code in this paper is similar to those obtained by Reiss, where only one peak exists. But due to the distinct geometrical model used in Watta’s analysis, there are two power peaks; that means the geometrical model has a great influence on the power distribution of the FA.
Fig. 14. Volume average temperature of the water in two passes of Chinese FA.
It can be found from Table 4 that with the decrease of the water density in water rod and assembly gap, the keff of the Europe FA will decrease more and faster. On the contrary, the keff of CSR1000 FA does not change a lot especially when the water density decreases from 800 kg/m3 to 700 kg/m3 ; there is even a slight increase of the keff , until the water density is lower than 500 kg/m3 , and the decrease of the keff will be faster. So it can be concluded that under the normal operation condition, the density feedback effect plays a minor role on the axial power distribution of the CSR1000 FA. And from the comparison of Europe and Chinese FAs, it is thought that the ratio of water volume to fuel volume is the main factor which results in the cosine shape axial power profile. In the Europe FA, the ratio is about 2.583 whereas it is 3.279 in Chinese FA. The larger water volume to fuel volume ratio will lead to the over-moderation effect; thus, the density feedback is not so sensitive to the change of the water density. 5. Conclusions In this paper, the coupled neutronics-T/H numerical simulation is carried out to obtain the axial power distribution of the FA. The following conclusions can be drawn from the comparison of the different results: (1) The partition method in the MCNP code will influence not only the final axial power distribution but also on the convergence process. Due to the large density change along the axial direction of the FA, it is suggested that there should be more partition zones in that direction. While the initialization value will not affect the final axial power distribution, but if an
(3) The coupled neutronics-T/H analysis of the CSR1000 FA shows that the power distribution is not so sensitive to the change of the water density, especially when the density of the moderator is more than 500 kg/m3 . This is thought to be caused by the lager water volume to fuel volume ratio of the Chinese FA. The cosine shape axial power profile could result in an excessive high cladding temperature. To flatten the axial power distribution and decrease the maximum cladding temperature, there should be a further improvement of the FA. Acknowledgement This work is supported by CNNC Key Laboratory on Nuclear Reactor Thermal Hydraulics Technology. References Bitterman, D., Squarer, D., Schulenberg, T., Oka, Y., Aksan, N., Maracy, C., KyrkiRajmarki, R., Soiyri, A., Dumaz, P., Steuwe, D., 2003. A potential plant characteristics of high performance light water reactor. In: ICAPP, Cordoba. Cao, L.Z., Oka, Y., Ishiwatar, Y., Ikejiri, S., 2009. Three-dimensional core analysis on a super fast reactor with negative local void reactivity. J. Nucl. Eng. Des. 239, 408–417. Gu, H.Y., Cheng, X., Yang, Y.H., 2008. CFD analysis of thermal-hydraulic behavior of supercritical water in sub-channels. J. Nucl. Eng. Des., 1–11. Hu, P., Wilson, P.P.H., 2009. Supercritical water reactors steady state, burnup and transient analysis with extended PARCS/RELAP5. In: 4th International Symposium on Supercritical Water-Cooled Reactors, Heidelberg, Germany. Laltu Chandra, Lycklama à Nijeholt, J.-A., Visser, D.C., Roelofs, F., 2009. CFD analyses of the influence of wire wraps spacers on heat transfer at supercritical conditions. In: 4th International Symposium on Supercritical Water-Cooled Reactors, Heidelberg, Germany. Los Alamos National Laboratory, 2003. MCNP5 Manual. Reiss, T., Fehér, S., Czifrus, Sz., 2008. Coupled neutronics and thermohydraulics calculations with burn-up for HPLWRs. J. Prog. Nucl. Energy, 1–10. Schulenberg, T., Starflinger, J., 2007. Core design concepts for high performance light water reactors. Nucl. Eng. Technol. 39 (4). Shan, J.Q., Chen, W., Leung, L.K.H., 2008. Coupled neutronics/thermal-hydraulics analysis of CANDU-SCWR fuel assembly. In: 4th International Symposium on Supercritical Water-Cooled Reactors, Heidelberg, Germany. US DOE Nuclear Energy Research Advisory Committee, 2002. A technology road map for generation IV nuclear energy systems, technical report. In: Generation IV International Forum. Waata, C., Schulenberg, T., Cheng, X., Starfinger, J., Lauruen, E., 2005. Coupling of MCNP with a sub-channel code for analysis of a HPLWR fuel. In: NURETH-11, Avigon, France (paper 021).
Contents lists available at SciVerse ScienceDirect
Nuclear Engineering and Design journal homepage: www.elsevier.com/locate/nucengdes
The axial power distribution validation of the SCWR fuel assembly with coupled neutronics-thermal hydraulics method Xi Xi, Zejun Xiao ∗ , Xiao Yan, Yongliang Li, Yanping Huang CNNC Key Laboratory on Nuclear Reactor Thermal Hydraulics Technology, Nuclear Power Institute of China, Chengdu 610041, China
h i g h l i g h t s CFX and MCNP codes are suitable to calculate the axial power profile of the FA. The partition method in the calculation will affect the final result. The density feedback has little effect on the axial power profile of CSR1000 FA.
a r t i c l e
i n f o
Article history: Received 6 September 2012 Received in revised form 13 December 2012 Accepted 8 January 2013
a b s t r a c t SCWR (super critical water reactor) is one of the IV generation nuclear reactors in the world. In a typical SCWR the water enters the reactor from the cold leg with a temperature of 280 ◦ C and then leaves the core with a temperature of 500 ◦ C. Due to the sharp change in temperature, there is a huge density change of the water along the axial direction of the fuel assembly (FA), which will affect the moderating power of the water. So the axial power distribution of the SCWR FA could be different from the traditional PWR FA.In this paper, it is the first time that the thermal hydraulics code CFX and neutronics code MCNP are used to analyze the axial power distribution of the SCWR FA. First, the factors in the coupled method which could affect the result are analyzed such as the initialization value or the partition method especially in the MCNP code. Then the axial power distribution of the Europe HPLWR FA is obtained by the coupled method with the two codes and the result is compared with that obtained by Waata and Reiss. There is a good agreement among the three kinds of results. At last, this method is used to calculate the axial power distribution of the Chinese SCWR (CSR1000) FA. It is found the axial power profile of the CSR1000 FA is not so sensitive to the change of the moderator density. © 2013 Elsevier B.V. All rights reserved.
1. Introduction SCWR is one kind of the IV generation nuclear reactors which was recommended by GIF in 2002. Due to the high thermal efficiency, compact system as well as the abundant operation and manufacture experience of the fossil supercritical power plant, SCWR is the most promising one in the future (US DOE, 2002). But before the commercial deployment of the SCWR, a lot of key techniques should be solved, one of which is the neutronics-T/H characteristics of the FA. In a typical SCWR, when the coolant flows into the core through the cold legs, it will be divided into three parts, one part flows towards the downcomer and the other parts flow towards the upper plenum and then into the water rod and the gap between the neighboring FAs as the moderator. Then the three different water will mix in the lower plenum and flow together into
∗ Corresponding author. Tel.: +86 028 85906199; fax: +86 028 85908889. E-mail address: fabulous [email protected] (Z. Xiao). 0029-5493/$ – see front matter © 2013 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.nucengdes.2013.01.031
the core again as the coolant shown in Fig. 1 (Bitterman et al., 2003). During this process, not only the density of the coolant changes a lot along the axial direction but also the density of the moderator will change simultaneously, which could result in the heterogeneous moderating power distribution along the axial direction of the FA. So the axial power distribution of the SCWR FA could be different from the traditional PWR FA which follows a cosine shape axial power profile. In 2005, a coupled neutronics-thermal hydraulics calculation was carried out by Waata to analyze the axial power distribution of the Europe HPLWR FA as well as the factors which could affect it such as the burn-up, density feedback effect, Doppler effect, etc. (Waata et al., 2005). In the calculation, neutronics code MCNP and subchannel T-H code STAFAS are used and two inactive parts of the FA are also modeled. The result shows that the axial power profile of the HPLWR FA does not follow the cosine shape when the burn-up is not taken into account. In fact there are two peaks, the larger one is close to the lower part and the smaller one is close to the upper part of the FA. The asymmetric power profile is mainly caused by the density feedback effect and has little relationship with the Doppler
158
X. Xi et al. / Nuclear Engineering and Design 258 (2013) 157–163
Nomenclature Cp E H N ˙ m q Q S ST Su Sv Sw T u u
v w
constant pressure specific heat (J/K kg) fission energy of the FA specific enthalpy (J/kg) partition number of the FA in the axial direction Mass flow rate (kg/s) average heat flux of the FA power of the FA surface of the fuel rods source term in energy equation (kg K/m3 s) source term in x momentum equations (Pa/m) source term in y momentum equations (Pa/m) source term in z momentum equations (Pa/m) temperature (K) vector of the velocity (m/s) velocity in x direction (m/s) velocity in y direction (m/s) velocity in z direction (m/s)
Greek letters source term in thermal conductivity equation (W/m3 ) thermal conductivity (W/m K) viscosity (Pa s) density of the fluent (kg/m3 ) Subscript C coolant region D downcomer G assembly gap region i the ith section of the FA in inlet out outlet t total W water rod region
although the axial power profile obtained by Reiss agrees well with Waata’s, there is only one power peak which is close to the bottom of the FA. Furthermore, both their results indicate that not only the burn-up but also the enrichment has a great influence on the power profile of the Europe HLWR FA. In 2009, the T/H code ATHAS and MCNP code were coupled to analyze the T/H characteristics as well as the axial power distribution of the SCWR-CANDU FA by Shan et al. (2008). The result shows that the density change of the water has little effect on the axial power distribution of the FA. The power profile follows the cosine shape approximately. Except for the research mentioned above, there are also some other coupled neutronics-thermal hydraulics analysis of the SCWR, which focuses on the power distribution of the whole core but they will not be discussed in detail here anymore (Cao et al., 2009; Hu et al., 2009). In this paper, the CFX 13.0 code and MCNP4C code are coupled for the axial power distribution analysis of the SCWR FA. It should be mentioned here that to avoid hot spots caused by a nonuniform power profile or by uncertainties and allowances for operation, a three passes SCWR is designed by Schulenberg and Starflinger (2007) to replace the single pass SCWR recommended by Bitterman et al. (2003). But the axial power profile of the single pass Europe SCWR is still selected for the code validation. In the calculation, only the density change effect of the water is taken into account and the Doppler effect is neglected. The following factors such as initialization value for the calculation and zone partition method in the MCNP code which could affect the final coupled result are discussed in detail. The burn-up and enrichment of the FA are not taken into consideration here. Finally, the axial power distribution of the CSR1000 FA is obtained with the same method. But it should be noted here that the CSR1000 is a two pass core design and for each pass the T/H boundary conditions such as mass flow rate, inlet temperature are different, so two types of axial power distributions will be attained finally.
2. Geometrical models and boundary conditions 2.1. Geometrical models
effect. In 2008, another coupled neutronics-thermal hydraulics analysis was carried out by Reiss with a self-developed T/H code and MCNP code (Reiss et al., 2008). In Reiss’s model, only the active part of the FA is taken into consideration and the result shows that
Both the Europe SCWR and CSR1000 FA are 1/8th symmetrical, so in the coupled numerical analysis only 1/8th of the FA is taken into account. Table 1 shows the main geometrical parameters of the
Fig. 1. Water route of the Europe and Chinese SCWR.
X. Xi et al. / Nuclear Engineering and Design 258 (2013) 157–163
159
Table 1 Parameters of the fuel assemblies. Parameters
Active height Fuel rod diameter Cladding thickness Fuel pin diameter Pitch/diameter ratio Moderator box lengthen Gap between FA Gap between rod and box wall Moderator box wall thickness FA box thickness Gaseous core diameter
Version Europe
Chinese
4200 8 0.5 6.7 1.150 26 5 1 0.4 1 –
4200 9.5 0.57 8.19 1.105 50.9 10 1 1 2 1.5
Europe and Chinese FAs. The Europe SCWR is the same as which analyzed by Waata and Reiss. Although the configurations of the two fuel assemblies are nearly the same, the CSR1000 FA is larger than the Europe FA, which has 56 fuel rods. Furthermore, there is a gaseous core in the center of the CSR1000 fuel pin, filled with the helium. And the control rod is not taken into account here. To find how the partition method in the MCNP code could affect the final axial power profile of the SCWR FA, four kinds of partition methods are used in this paper, which are P1, P2, P3 and P4 respectively. Figs. 2–4 show the zone partition methods of the cross-section of the Europe FA, which correspond to P1, P2 and P3, the partition method of P4 is similar to the P2, the only difference is that the FA is divided into 21 zones axially in P1–P3, while in P4, there are 42 axial zones. Fig. 5 shows the partition method of the CSR1000 FA which will be used for the coupled analysis. It can be found that the water rod, coolant region and the gap between the FAs are not further divided for the MCNP code as shown in Fig. 2. While in Fig. 3, the coolant region is divided into 7 zones. Moreover, not only the water rod but also the FA gap is further divided. Considering that the gradient of the water density is vertical to the water rod box and FA box, the water in these two regions are divided into three parts from the inner to the outer respectively. In Fig. 4, the partition method of the coolant region is similar to Waata’s, which is a 9 subzones division and the partition method of the water rod region and gap region is the same as that in Fig. 3. With the change of the partition method for the MCNP code, the mesh type for the CFX code should be changed
Fig. 2. Partition method of the Europe FA (P1) 101 assembly gap, 102 coolant, 103 water rod, 201–207 fuel rods, 301–303 cladding and box walls, 401 helium.
Fig. 3. Partition method of the Europe FA (P2) 101–103 water rod, 104–110 coolant, 111–113 assembly gap, 201–207 fuel rods, 301–303 cladding and box walls, 401 helium.
simultaneously. For example, the mesh type of the coolant region in P3 is triangular prism, whereas it is hexahedral in P1 and P2. But the mesh types of the water rod and gap regions are the same from P2 to P4, which are hexahedral. The mesh numbers for the Europe fuel assembly are 3,200,000 and 3,000,000 respectively with triangular prism and hexahedral mesh. It is 3,500,000 for Chinese fuel assembly with triangular prism mesh. 2.2. Boundary conditions The boundary conditions mentioned here not only include the temperature, mass flow rate, etc., but also the material properties such as the density and the enrichment of the fuel pin. Table 2 shows the thermal hydraulics boundary conditions for CFX code and Table 3 is the neutronics boundary conditions for MCNP code. It should be noted here again that the CSR1000 is a typical two flow pass reactor. The symbol (*) in Table 2 means that the inlet temperature of the coolant could not be input directly because it is related to the outlet
Fig. 4. Partition method of the Europe FA (P3) 101–103 water rod, 104–1120 coolant, 113–115 assembly gap, 201–207 fuel rods, 301–303 cladding and box walls, 401 helium.
160
X. Xi et al. / Nuclear Engineering and Design 258 (2013) 157–163
and energy conservation equations but also the 3D heat conductivity equation should be solved which are as shown below. But it should be mentioned here that the conductivity in the fuel rod is not taken into account here because the Doppler effect of the fuel rod is neglected.
∇ · (u) = 0
(2)
∂p ∇ · (uu ) = ∇ · (∇ u) − + Su ∂x
(3)
∂p + Sv ∂y
(4)
∇ · (vu ) = ∇ · (∇ v) −
∇ · (wu ) = ∇ · (∇ w) −
∇ · (T u ) = ∇ · c Fig. 5. Partition method of the CSR1000 FA 101 assembly gap, 102 coolant, 103 water rod, 201–202 fuel rods, 301–303 cladding and box walls, 401 helium. Table 2 Boundary conditions for T/H analysis of 1/8th FA. Parameters
System pressure Water properties Tin of the water rod Tin of the assembly gap Tin of the coolant region Min of the water rod Min of the assembly gap Min of the coolant region FA power
Version Europe
Chinese I
Chinese II
25 MPa IAPWS IF97 553 K 553 K * 0.0139 kg/s 0.0278 kg/s 0.167 kg/s 327.5 kW
25 MPa IAPWS IF97 553 K 553 K 553 K 0.0389 kg/s 0.0706 kg/s 0.2549 kg/s 199 kW
25 MPa IAPWS IF97 553k 553 K 657 K 0.0517 kg/s 0.0706 kg/s 0.3697 kg/s 199 kW
temperature of the water rod and the assembly gap. The following correlation can be used to obtain HinC in the coolant region: HinC =
˙ inD · HinD + m ˙ inW · HoutW + m ˙ inG · HoutG m ˙t m
(1)
∂2 T ∂2 T ∂2 T + + 2 2 ∂x ∂y ∂z 2
(5)
+ ST
(6)
=
Cp
(7)
The difference scheme of the convection terms in the conservation equations is second order upwind scheme, and the SST model is selected as the turbulence model in the calculation, which is suitable for the pseudocritical conditions. When the residual is less than 10−3 , the T/H solution is thought to be converged and the water density distribution in the FA will be obtained. The neutron transport equation in the FA will be solved by MCNP which is a general-purpose, continuous-energy, generalizedgeometry, time-dependent, coupled neutron/photon/electron Monte Carlo transport code. It can be used in several transport modes: neutron only, photon only, electron only, combined neutron/photon transport where the photons are produced by neutron interactions, neutron/photon/electron, photon/electron, or electron/photon (LANL, 2003). With the MCNP code, not only the keff but also the power distribution of the FA will be obtained. The power of the fuel rod here will be calculated with F7 tally. The particle number in one cycle is 50,000, and inactive cycle is 500, and total cycle is 1000. With the MCNP calculation, the mass averaged fission energy can be obtained for each section of the FA, and the average heat flux of the fuel rods in this section can be attained with the following equation: qi =
3. Physical models and coupled method
∇T Cp
∂p + Sw ∂z
Ei · Q
N E i=1 i
(8) · Si
3.1. Physical models 3.2. Coupled method The physical models here refer to the T/H model and neutronics model in the coupled calculations. The T/H phenomenon in the SCWR FA includes the heat conductivity in the water box wall and FA wall, the convection heat transfer in the water rod, coolant, assembly gap, etc. So not only the steady state mass, momentum Table 3 Boundary conditions for neutronics analysis of 1/8th FA. Parameters
Corner fuel rod enrichment Other fuel rod enrichment Fuel rod density Material of the cladding Material of the water box Material of the assembly box Density of the cladding water box and assembly box Material of the gas Cross section library
Version Europe
Chinese
4% 5% 10,600 kg/m3 Alloy 316 Alloy 316 Alloy 316 7450 kg/m3
4.3% 5.7% 10,600 kg/m3 Alloy 310S Alloy 310S Alloy 310S 7450 kg/m3
Helium ENDF6
Helium ENDF6
The external coupled method is used here and the flowchart is shown in Fig. 6. If the T/H calculation will be carried out first, an assumed power profile will be given as the boundary condition for the CFX code and the water density distribution will be obtained when the solution is converged. Then the density will be input to the MCNP code to obtain the new power profile. If the relative error of the new profile and the old one is small enough, the coupled numerical simulation will be stopped and it is assumed that a converged solution is attained. The calculation is run in a parallel computer which has 24 CPUs, so for one coupled iteration it will take only 30 min. To decrease the iteration steps and get a converged result, a relaxation factor 0.8 is used in the coupled calculation. 4. Results and analysis 4.1. Influence of the partition method As a numerical simulation approach, not only the mesh in the CFX code but also the partition method in the MCNP code could
X. Xi et al. / Nuclear Engineering and Design 258 (2013) 157–163
161
Fig. 6. Flowchart of the coupled method. Fig. 8. Comparison of the average axial heat flux with different partition methods in MCNP code.
affect the precision and convergence process of the solution. The mesh influence of the CFD code has already been discussed by so many researchers when the T/H analysis of the FA was carried out and there will be no further discussion here on this (Gu et al., 2008; Laltu Chandra et al., 2009). But the influence of the partition method of the MCNP code is not well analyzed especially for the coupled neutronics-T/H calculation, so how it will affect the coupled solution is discussed in detail here and four different partition methods are selected which correspond to P1–P4. During the coupled numerical simulation, it is found that the solution is too difficult to converge with the partition method P1. After 10 iterations, there is still a periodic oscillation of the axial power profile which is shown in Fig. 7. To make the comparison with different partition methods possible, the result of P1 is the algebraic average value from iteration 10 to iteration 15. The comparison in Fig. 8 shows that although different partition methods are used in the coupled calculation, there is a good agreement among these four axial heat flux profiles. With the increase of the partition zones in the transverse coolant region, the peak value of heat flux will decrease slightly. On the other hand, with the increase of the partition zones in the axial region, there is a 4% increase of the peak value. So it can be concluded that the partition method in the MCNP code will affect the average axial heat flux profile of the FA. But compared to the axial partition method, the transverse partition method plays a minor influence on the axial power distribution.
Fig. 7. Average axial heat flux of different iteration steps with partition method P1.
4.2. Influence of the initialization value In a numerical simulation, the initialization value will probably affect not only the convergence process but also the final result. To analyze the influence of the initialization value in the coupled calculation, two different heat flux profiles are used as initial condition in the CFX code, which are homogeneous profile and standard cosine profile respectively. The partition method in the calculation is P3 and the obtained final coupled axial heat flux profile is shown in Fig. 9. The comparison shows that the average heat flux profiles obtained with different initialization values are nearly the same, which means it has no effect on the final result. But when the homogeneous heat flux is selected as the initialization value, it costs 8 iteration steps for the solution to be converged, while the cosine shape heat flux is selected, 11 iteration steps are needed to obtain a converged solution. 4.3. Results comparison To validate whether the coupled method in this paper is appropriate for the coupled analysis, the axial power distribution of the Europe SCWR FA obtained with CFX and MCNP codes is compared with that obtained by Waata and Reiss. The partition method here is P3 and Fig. 10 shows the comparison of the average axial heat flux distributions of the Europe FA.
Fig. 9. Comparison of the average axial heat flux with different initialization values in MCNP code.
162
X. Xi et al. / Nuclear Engineering and Design 258 (2013) 157–163
Fig. 10. Comparison of the average axial heat flux.
Fig. 11. Flow directions in two passes of Chinese FA.
From the comparison it can be found that the average axial heat flux profile obtained by Reiss and Xi is nearly the same. There is a power peak at the close of the coolant inlet and then decreases gradually. The maximum heat flux in Xi’s result is larger than Reiss’s with a relative error of 9%. But the heat flux profile obtained by Waata is a little different from Xi’s, which has two power peaks. This is thought to be caused by the two inactive parts in Waata’s model, especially the one close to the coolant outlet where the water density in the water rod and FA gap is larger, which results in a stronger moderating power. So the upper part heat flux of the FA in Watta’s result is larger than Reiss’s and Xi’s. The result comparison demonstrates that CFX is an appropriate code which could be used for the coupled netronics-T/H analysis of the FA. And the results obtained by CFX code and subchannel code are so close to each other if the same geometrical model are used. 4.4. Axial power distribution of the CSR1000 FA Fig. 12. Average axial heat flux of Chinese FA.
The CSR1000 is the Chinese SCWR designed by Nuclear Power Institute of China. It is a typical two pass SCWR with an inlet temperature of 280 ◦ C and an outlet temperature of 500 ◦ C. Due to the special geometrical configuration of the CSR1000 FA as well as the two flow passes, the axial power distribution of the CSR1000 could be different from the Europe SCWR FA. So the CFX and MCNP codes are coupled here to calculate the axial power profile of the CSR1000 FA. It has been demonstrated in this paper that the transverse partition method in the MCNP code has minor effect on the average axial power distribution. When the coupled neutronics-T/H analysis of the CSR1000 FA is carried out, the simplest transverse partition method is used, which is shown in Fig. 5. The flow direction in two passes and the two types average heat flux profiles of the FA are shown in Figs. 11 and 12. Figs. 13 and 14 show the volume average water density and temperature in pass 1 and pass 2 respectively. It is found that the changes of the water density and temperature along the FA are totally different in the two passes. For example in pass 1, the water density increases with the increase of the axial height, while in pass 2, it will decrease first and then increase to 600 kg/m3 . Furthermore, the change of the water density and temperature in pass 1 is much larger than that in pass 2. But in Fig. 12, it is found that the average heat flux profiles of FA in different passes are nearly the same, which is a cosine profile. So it is assumed that the power distribution is not sensitive to the change of the water density. To validate the assumption mentioned above and find out the influence of the density feedback effect on the power distribution, another neutronics analysis is carried out with MCNP code. In this
calculation, there is no axial zone partition for the FA and the average density is used. The density of the water in the coolant region is supposed to be 300 kg/m3 and will not change anymore, then the density of the water in the water rod and assembly gap will decrease from 800 kg/m3 to 300 kg/m3 and for each 100 kg/m3 a calculation will be run to obtain the keff . Table 4 shows the change of the keff with the decrease of the moderator density.
Fig. 13. Volume average density of the water in two passes of Chinese FA.
X. Xi et al. / Nuclear Engineering and Design 258 (2013) 157–163
163
Table 4 keff of the FA with different moderator densities. Water rod density (kg/m3 )
Coolant density (kg/m3 )
Assembly gap density (kg/m3 )
keff of the Europe
keff of the CSR1000
800 700 600 500 400 300
300 300 300 300 300 300
800 700 600 500 400 300
1.15960 1.14801 1.13400 1.11530 1.09188 1.06047
1.27265 1.27350 1.27068 1.26054 1.23942 1.20447
appropriate initialization is applied, the solution will converge more quickly. (2) CFX is an appropriate T/H code which could be used for the coupled neutronics-T/H analysis of the SCWR FA. The average axial heat flux profile of the Europe FA obtained with CFX code and MCNP code in this paper is similar to those obtained by Reiss, where only one peak exists. But due to the distinct geometrical model used in Watta’s analysis, there are two power peaks; that means the geometrical model has a great influence on the power distribution of the FA.
Fig. 14. Volume average temperature of the water in two passes of Chinese FA.
It can be found from Table 4 that with the decrease of the water density in water rod and assembly gap, the keff of the Europe FA will decrease more and faster. On the contrary, the keff of CSR1000 FA does not change a lot especially when the water density decreases from 800 kg/m3 to 700 kg/m3 ; there is even a slight increase of the keff , until the water density is lower than 500 kg/m3 , and the decrease of the keff will be faster. So it can be concluded that under the normal operation condition, the density feedback effect plays a minor role on the axial power distribution of the CSR1000 FA. And from the comparison of Europe and Chinese FAs, it is thought that the ratio of water volume to fuel volume is the main factor which results in the cosine shape axial power profile. In the Europe FA, the ratio is about 2.583 whereas it is 3.279 in Chinese FA. The larger water volume to fuel volume ratio will lead to the over-moderation effect; thus, the density feedback is not so sensitive to the change of the water density. 5. Conclusions In this paper, the coupled neutronics-T/H numerical simulation is carried out to obtain the axial power distribution of the FA. The following conclusions can be drawn from the comparison of the different results: (1) The partition method in the MCNP code will influence not only the final axial power distribution but also on the convergence process. Due to the large density change along the axial direction of the FA, it is suggested that there should be more partition zones in that direction. While the initialization value will not affect the final axial power distribution, but if an
(3) The coupled neutronics-T/H analysis of the CSR1000 FA shows that the power distribution is not so sensitive to the change of the water density, especially when the density of the moderator is more than 500 kg/m3 . This is thought to be caused by the lager water volume to fuel volume ratio of the Chinese FA. The cosine shape axial power profile could result in an excessive high cladding temperature. To flatten the axial power distribution and decrease the maximum cladding temperature, there should be a further improvement of the FA. Acknowledgement This work is supported by CNNC Key Laboratory on Nuclear Reactor Thermal Hydraulics Technology. References Bitterman, D., Squarer, D., Schulenberg, T., Oka, Y., Aksan, N., Maracy, C., KyrkiRajmarki, R., Soiyri, A., Dumaz, P., Steuwe, D., 2003. A potential plant characteristics of high performance light water reactor. In: ICAPP, Cordoba. Cao, L.Z., Oka, Y., Ishiwatar, Y., Ikejiri, S., 2009. Three-dimensional core analysis on a super fast reactor with negative local void reactivity. J. Nucl. Eng. Des. 239, 408–417. Gu, H.Y., Cheng, X., Yang, Y.H., 2008. CFD analysis of thermal-hydraulic behavior of supercritical water in sub-channels. J. Nucl. Eng. Des., 1–11. Hu, P., Wilson, P.P.H., 2009. Supercritical water reactors steady state, burnup and transient analysis with extended PARCS/RELAP5. In: 4th International Symposium on Supercritical Water-Cooled Reactors, Heidelberg, Germany. Laltu Chandra, Lycklama à Nijeholt, J.-A., Visser, D.C., Roelofs, F., 2009. CFD analyses of the influence of wire wraps spacers on heat transfer at supercritical conditions. In: 4th International Symposium on Supercritical Water-Cooled Reactors, Heidelberg, Germany. Los Alamos National Laboratory, 2003. MCNP5 Manual. Reiss, T., Fehér, S., Czifrus, Sz., 2008. Coupled neutronics and thermohydraulics calculations with burn-up for HPLWRs. J. Prog. Nucl. Energy, 1–10. Schulenberg, T., Starflinger, J., 2007. Core design concepts for high performance light water reactors. Nucl. Eng. Technol. 39 (4). Shan, J.Q., Chen, W., Leung, L.K.H., 2008. Coupled neutronics/thermal-hydraulics analysis of CANDU-SCWR fuel assembly. In: 4th International Symposium on Supercritical Water-Cooled Reactors, Heidelberg, Germany. US DOE Nuclear Energy Research Advisory Committee, 2002. A technology road map for generation IV nuclear energy systems, technical report. In: Generation IV International Forum. Waata, C., Schulenberg, T., Cheng, X., Starfinger, J., Lauruen, E., 2005. Coupling of MCNP with a sub-channel code for analysis of a HPLWR fuel. In: NURETH-11, Avigon, France (paper 021).