Next Article in Journal
Discussing Landscape Compositional Scenarios Generated with Maximization of Non-Expected Utility Decision Models Based on Weighted Entropies
Next Article in Special Issue
A LiBr-H2O Absorption Refrigerator Incorporating a Thermally Activated Solution Pumping Mechanism
Previous Article in Journal
Complex and Fractional Dynamics
Previous Article in Special Issue
Scaling Relations of Lognormal Type Growth Process with an Extremal Principle of Entropy
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Response Surface Methodology Control Rod Position Optimization of a Pressurized Water Reactor Core Considering Both High Safety and Low Energy Dissipation

School of Energy Science and Engineering, Harbin Institute of Technology, Harbin 150001, China
*
Author to whom correspondence should be addressed.
Entropy 2017, 19(2), 63; https://doi.org/10.3390/e19020063
Submission received: 30 November 2016 / Revised: 15 January 2017 / Accepted: 6 February 2017 / Published: 10 February 2017
(This article belongs to the Special Issue Advances in Applied Thermodynamics II)

Abstract

:
Response Surface Methodology (RSM) is introduced to optimize the control rod positions in a pressurized water reactor (PWR) core. The widely used 3D-IAEA benchmark problem is selected as the typical PWR core and the neutron flux field is solved. Besides, some additional thermal parameters are assumed to obtain the temperature distribution. Then the total and local entropy production is calculated to evaluate the energy dissipation. Using RSM, three directions of optimization are taken, which aim to determine the minimum of power peak factor Pmax, peak temperature Tmax and total entropy production Stot. These parameters reflect the safety and energy dissipation in the core. Finally, an optimization scheme was obtained, which reduced Pmax, Tmax and Stot by 23%, 8.7% and 16%, respectively. The optimization results are satisfactory.

Graphical Abstract

1. Introduction

In China, energy consumption has grown rapidly driven by the improvement of the economic level. In order to balance economic development and environmental protection, China’s government has established a green development strategy. An important measure is to adjust the energy structure and promote the consumption of clean energy sources such as nuclear energy. Thus, the use of nuclear energy increased from 44.19 billion kWh in 2000 to 290.75 billion kWh in 2013 [1], and 26 nuclear power units (28,528 MWe) are currently under construction in China [2,3].
The design of a nuclear reactor system includes shielding design, thermodynamics, fluid flow and heat transfer, fuel element design, radioactivity releases, etc. [4]. In nuclear physics calculations, the neutron flux density and its distribution are usually the primary target. The Neutron Transport Equation (NTE) is used to describe the variation of neutron flux density in a physical field. For simplification, the Neutron Diffusion Equation (NDE) is also used. These equations can be solved by the Finite Different Method (FDM) [5,6,7], the Finite Volume Method (FVM) [8,9,10,11], the Finite Element Method (FEM) [10,12,13,14], etc. Normally, the neutron flux density distribution matches the fission heat source distribution.
Heat transfer is also an important physical process in nuclear reactor core. It is an irreversible process which increases the entropy of the whole system. The differences in heat transfer in different regions causes a distribution of local entropy production. In 1979, Bejan [15] presented the fundamental equations describing entropy production due to finite temperature gradients and fluid flow. Then Bejan [16,17,18] discussed the irreversibility in thermodynamic systems by entropy production analysis. Since then, many researchers have started to use the entropy production method to evaluate heat transfer processes in various structures. For example, Ibáñez et al. [19] investigated the entropy production minimization of a solid slab with uniform internal heating and asymmetric convective cooling. Makinde and Aziz [20] analyzed the inherent irreversibility and thermal stability in the model of a long hollow cylinder with temperature-dependent internal heating and asymmetric convective cooling boundaries. Aziz and Khan [21] analyzed the entropy production of steady conduction in the model of a hollow sphere with temperature-dependent internal heating and asymmetric convective cooling boundaries. Malvandi et al. [22] studied the entropy production of steady two-dimensional boundary layer flow of nano-fluids over a flat plate. Torabi and Zhang [23] used classical entropy production analysis to investigate heat transfer in cooled homogenous and functionally graded material slabs with variation of internal heat generation with temperature, and convective-radiative boundaries.
For coupling the neutronic and thermal-hydraulic calculation, the traditional technique is the Operator Splitting (OS) method [24,25]. By this method, different physical processes are solved separately. The output of neutronic calculation is the input of thermal-hydraulic part. Some other advanced multi-physics algorithms have also been developed. One of the most famous algorithms is the Jacobian-Free Newton-Krylov (JFNK) method [26,27]. By the JFNK method, several physical processes are treated together as an entire unit. The biggest advantage of the OS method is its ease of implementation, while JFNK is more suitable for nonlinear problems and has higher accuracy.
Response Surface Methodology (RSM) is a collection of mathematical and statistical techniques used to explore the relationships between several independent variables and one or more response variables. This methodology was proposed by Box and Wilson in the 1950s [28,29]. Since then this methodology has been widely used in analytical chemistry [28], bioprocessing [29,30], structural reliability [31,32], food chemistry [33,34,35], etc. The biggest advantage of this methodology is that the RSM model is easy and convenient to establish, even when little information about the process can be obtained.
This paper aimed to optimize the insertion positions of the control rods in a Pressurized Water Reactor (PWR). The neutron diffusion and heat transfer processes are coupled and both high safety and low energy dissipation are taken into consideration. The chosen safety factors are the power peak factor and maximum temperature in the core, and the energy dissipation is calculated by the entropy production method. RSM is utilized to find a more improved scheme.

2. Description of Problem

The International Atomic Energy Agency (IAEA) has published a three-dimensional PWR benchmark problem named the 3D-IAEA [36]. It is a typical PWR problem which is widely used in verification of neutron diffusion calculation codes. This PWR is taken as a study object in this paper to demonstrate the RSM optimization of control rod insertion positions. The nuclear data given in [36] are used to solve the neutron flux distribution. Some essential thermal parameters are also assumed to calculate the temperature field and local and total entropy production.
There are 177 groups of Fuel Assemblies (FAs) arranged in the reactor core, and 13 groups of these FAs are have inserted control rods. Five kinds of materials are applied in this core. A quarter of horizontal and one half of vertical cross section of the 3D-IAEA problem are illustrated in Figure 1 and Figure 2, respectively. The multi-group neutron diffusion approximation that most widely used for commercial reactors is the two-energy groups approximation [37,38], which is also applied in 3D-IAEA benchmark problem. The neutron energy spectrum is divided into two groups, which are called fast and thermal neutrons, respectively. The two-group nuclear data are given in Table 1. Some additional thermal parameters are given in Table 2.

3. Numerical Methodology

3.1. Neutron Diffusion

The NTE based on transport theory can describe the neutron behaviors in more details, but to balance of amount of calculation and perduring accuracy, in the engineering area the diffusion theory and NDE are always used to solve the neutron flux, which can be written as:
[ 1 v g t · D g ( r , t ) + Σ t , g ( r , t ) ] φ g ( r , t ) = g = 1 G s , g g ( r , t ) φ g ( r , t ) + 1 k eff χ g g = 1 G υ g Σ f , g ( r , t ) φ g ( r , t ) + S e , g ( r , t )
where g is the index of energy group, g = 1,2,…,G; G is the total number which the energy spectrum is divided into; v is the velocity of neutrons, cm·s−1; t is the time, s; D is the diffusion constant, cm; Σt is the total cross section, cm−1; Σs,g’→g is the scattering cross section from g’ energy group to g energy group, cm−1; Σf is the fission cross section, cm−1; r is the position, cm; φ is the neutron scalar flux density, neutrons·cm−2·s−1; υ is the average number of neutrons that are emitted from each fission process; χ is the fission spectrum function; keff is the multiplication factor; Se is the extra neutron source, neutrons·cm−3·s−1.
When using in time independent problem with no extra neutron source and no upscattering neutrons, Equation (1) can be simplified as:
D g ( r ) φ g ( r ) + Σ t , g ( r ) φ g ( r ) = g = 1 g 1 s , g g ( r ) φ g ( r ) + 1 k eff χ g g = 1 G υ g Σ f , g ( r ) φ g ( r )
Thus, for the two-group diffusion problem, the governing equations can be written as:
{ D 1 ( r ) φ 1 ( r ) + [ a , 1 ( r ) + s , 1 2 ( r ) ] φ 1 ( r ) = 1 k eff υ 1 Σ f , 2 ( r ) φ 2 ( r ) D 2 ( r ) φ 2 ( r ) + a , 2 ( r ) = s , 1 2 ( r ) φ 1 ( r )
where Σa is the absorb cross section, cm−1, and t , 1 ( r ) = a , 1 ( r ) + s , 1 2 ( r ) , t , 2 ( r ) = a , 2 ( r ) .
The local fission reaction rate per unit volume is:
R f ( r ) = g = 1 G ν f , g ( r ) φ g ( r )
where Rf is the fission reaction rate, reaction times·cm3·s−1.
It can be assumed that the total heat is generated at the position where the fission reaction occurs, therefore the power of local fission reaction can be defined as [39]:
P ( r ) = E f R f ( r )
where Ef is the total heat generation in one fission reaction.
According to the thermal parameters shown in Table 2, the reactor core is designed at the normal thermal power Pcore, which means that under the normal conditions, the summary of local powers is Pcore:
P core = V core P ( r )   d V
Average subassembly power Pk is defined as:
P k = 1 V k V k P ( r )   d V
where k is the designated number of FA which can be seen in Figure 1.
In order to compare our results with Ref. [36], a type of normalized neutron fluxes φ ¯ g ( r ) is introduced. In the normalized form, the volumetrically weighted average of released fission neutron fluxes in the reactor core is set as 1, which can be written as:
1 V core V core g = 1 G ν f , g ( r ) φ ¯ g ( r )   d V = 1
The aim of normalization is to eliminate the impact of initial conditions on the neutron flux distribution.

3.2. Heat Transfer and Energy Dissipation

The governing equation of conductive heat transfer can be written as:
ρ ( r , t ) c P ( r , t ) T ( r , t ) t = λ ( r , t ) T ( r , t ) + P ( r )
where ρ is the density, kg·m−3; cp is the specific heat capacity, J·kg−1·K−1; λ is the thermal conductivity, W·m−1·K−1; T is the temperature, K; Pcore is the total thermal power of the reactor core, W.
For the time independent problem, Equation (7) can be simplified as:
λ ( r ) T ( r ) + P ( r ) = 0
The irreversible energy dissipation can be calculated by the entropy production. Normally, the local entropy production of conductive heat transfer can be written as [17,18]:
S gen ( r ) = λ ( r ) [ T ( r ) ] 2 2 T ( r )
The total entropy production is:
S tot = V core S gen ( r ) d V

3.3. Response Surface Methodology

By RSM, if all the independent variables can be measured, the response surface can be expressed as:
Y = F ( X 1 , X 2 , ... , X n )
where Xi represents the independent variable, Y is the response variable.
It is assumed that the independent variables are continuous. The goal is to find a suitable approximate relationship between independent variables and the response variable. Usually, a second-order model is utilized [40,41], which can be written as:
Y = β 0 + i = 1 n β i X i + i = 1 n β i i X i 2 + i < j n β i j X i X j + ε
where β is the undetermined coefficient, ε is a random error.
Equation (14) can be written in matrix form:
Y = b X + ε
The solution of Equation (15) can be obtained by the matrix approach:
b = ( X T X ) 1 X T Y
In Equation (16), the superscript T represents matrix transposition, the superscript −1 represents matrix inversion.
There are several two-order designs for RSM. The Central Composite Design (CCD) and Box-Behnken Design (BBD) are two typical ones. In this paper, BBD is adopted to get the coefficients of Equation (14). Based on the symmetry of the reactor core along x and y axis, only one-eighth of the whole core need to be designated, which can be seen in Figure 1. Thus, the insertion positions of four control rods can be chosen as independent variables, namely Z1, Z2, Z3 and Z4, which represent the positions of FA.1, FA.18, FA.5 and FA.31, respectively.

4. Results and Discussion

4.1. Standard Problem Solution

The control rod insertion positions shown in Figure 2 were chosen as the standard and reference solution. This example solution is the basis of the optimization. The neutron diffusion equations and heat transfer equations are all solved by FVM. In each calculation, the considered domain is one quarter of the whole core with reflective and non-return external boundaries. The convergence criterion is set as maximum relative flux change on each inner iteration = 10−10, maximum keff change on outer iteration = 10−6. Four groups of regular hexahedral meshes are taken into consideration. The mesh quantities and calculation results are given in Table 3. The error in Table 3 represents the relative error between using a mesh group and its further refinement.
Considering both of the calculation accuracy and time cost, mesh No. 2 was selected to be applied for the next optimization scheme. This quantity of mesh is also used in [36], thus our neutron diffusion calculation results can be checked. The calculated keff is 1.02855 and the reference at the same mesh quantity is 1.02864, with the difference of 0.0087%. The results of fast and thermal neutron flux at the diagonal line on the x-y plane at the level of z = 195 cm are plotted in Figure 3 and Figure 4. A good match between calculated results and references can be seen.
The calculation local power, temperature and entropy production results are plotted in Figure 5, Figure 6 and Figure 7. It can be seen in Figure 5 that the control rods have an obvious effect on the neutron fission power. Both the vertical and horizontal cross sections show significant decreases of local power depending on the control rod positions. The comparison of horizontal cross section at different z levels shows that the local power near the center of the core is higher than at the upper level, as expected. Figure 6 shows that the regular temperature distribution is not as complex as the fission power distribution. It shows a continuous reduce from the core center to the edge. The control rod positions cannot be seen clearly in these temperature cloud pictures. The regular local entropy production distribution is far from the temperature distribution.
Figure 7 shows that the boundary centers of the reactor core are the peaks of the local entropy production, and the volume center of the reactor core is the valley of the local entropy production. The reason is that this core is symmetric at the vertical middle cross section, thus the middle cross section is heat insulated and the temperature gradient there is zero. The rapid decrease of temperature near the boundary center causes the entropy production peak.
In total, some lumped and feature parameters are selected to describe the thermal situation of the reactor core. These parameters are the power peak factor Pmax, the maximum temperature Tmax and the total entropy production Stot. These parameters’ values are shown in Table 4.

4.2. Response Surface Design

The BBD is introduced to apply the RSM. The insert positions of control rods Z1, Z2, Z3, Z4 are the input independent variables, and the feature parameters Pmax, Tmax, Stot are chosen as the response variables. The RSM can be written as follows:
Y = [ P max T max S tot ] = [ b 1 b 2 b 3 ] X + [ ε 1 ε 2 ε 3 ]
where:
X = [ 1 Z 1 Z 2 Z 3 Z 4 Z 1 2 Z 2 2 Z 3 2 Z 4 2 Z 1 Z 2 Z 1 Z 3 Z 1 Z 4 Z 2 Z 3 Z 2 Z 4 Z 3 Z 4 ]
In order to establish the coefficient matrix b, several experiments should be performed. Section 4.1 presents an example of a numerical experiment. The input variables are [Z1 Z2 Z3 Z4] = [202802020], and the response result Y is given in Table 4. Using BBD to design the RSM for this problem, 29 groups of numerical experiments must be executed. The input and response variables are shown in Table 5. The matrix b is obtained by Equation (17) and the fitting results are shown in Table 6.
The fitting results are usually evaluated by R-square. For this evaluation index, the closer to 1, the better result gained. The R-square value of the three response surface fitting is shown in Table 7.
It can be seen in Table 7 that, although Tmax is a local feature, depending on its simple regular distribution, the R-square of its response surface is the highest. The second is Stot since it is a total summary of local entropy production. It has the both global and local parameter features, so its sensitivity to the input variables is weaker than that of the local parameter. The R-square of the Pmax response surface is not as fine as the others, due to the complex local power distribution. As these response surfaces are introduced to find out the direction of optimization, rather than obtain the exact solution, this level of fitting can also be applied.

4.3. Rod Position Optimization

In this section, a simple code is programmed to traverse all four input variable ranges in the definition domain, which is set as [30 cm, 300 cm]. The coefficients in Table 6 and Equation (18) are used in this code to obtain the approximate response results. The interval step of the input variables is set as 10 cm. In every traversal calculation, one evaluation index (Pmax, Tmax or Stot) is selected as the optimization direction. The aim of optimization is to promote safety or reduce entropy dissipation, that is to say, Pmax, Tmax or Stot should be reduced.
Although these three indexes are key metrics for a reactor core, the importance rankings are not at the same level for each of them. For every nuclear energy system, safety is always a seriously concerned for the public [42], and its status is higher than economical efficiency. Tmax is the most obvious safety index, as the reactor materials may melt above a certain temperature limit, which would cause a serious accident. Pmax is also a type of safety index, which is associated with the level of power flattening. A more flattened power distribution means a smoother reactor operation is expected. The last is the economical index Stot, which represents the energy dissipation.
Therefore, our optimization procedure is: (1) use RSM to traverse all the control rod positions to find N1 positions approaching the lowest Tmax. The selected positions are put in set Ω1; (2) Use RSM again to traverse all the positions in the definition domain to find N2 positions approaching the lowest Pmax and name the selected positions set as Ω2; (3) Use FVM to obtain accurate indexes of position elements in Ω3, Ω3 = Ω1 ∩ Ω2, and sort them by Stot. The best solution is chosen as the final optimization scheme. The flow diagram of this optimization procedure is shown in Figure 8.
In this scheme, a smaller N1 or N2 means a narrower search range, but the indexes may become closer to the limit. In this case, they are set as N1 = N2 = 10. Finally, eight groups of control rod positions are selected, which are listed in Table 8.
It can be seen in Table 8 that the safety and entropy index will reach optimization when the central control rods are inserted deeply and the peripheral control rods are inserted shallowly. In the calculation series, it can be found that when the input rod positions are [30 30 300 300], the Pmax, Tmax and Stot reach the minimum, which means that both the power and temperature are flatter, and the energy dissipation reaches the minimum. In order to observe the physical field in depth, the distributions of local power, temperature and local entropy production are illustrated in Figure 9, Figure 10 and Figure 11.
The cloud pictures Figure 9 and Figure 11 share the same colorbar scale as Figure 5, Figure 6 and Figure 7, respectively. Thus the differences between these two groups can be clearly observed. Comparing Figure 9 with Figure 5, it can be found that the power flattening effect is better after the control rod adjustment. The central control rod has an obvious influence on the power distribution, which causes a power valley at the core center, so the temperature at the center line of x = y = 0, is not so high as the temperature in the standard problem, which can be seen in Figure 6 and Figure 10. The distribution of local entropy production is similar between Figure 7 and Figure 11, but the entropy production values are reduced after control rod adjustment.
In order to show some details of the optimization results, the local power, temperature and local entropy production at the diagonal line on midplane are shown in Figure 12. It can be seen in Figure 12a that the local power decreases near the core center and a power valley forms between the distance of 40 cm and 50 cm, as the control rods in FA.1 and FA.18 are inserted deeply. At the same time, with the control rods in FA.5 and FA.31 rise, the local power of the optimization scheme is higher than that of the standard problem in the distance range of [65 cm, 110 cm]. As a result, local power is flatter than in the standard problem and the power peak factor is reduced. For the temperature distribution, as shown in Figure 12b, the temperature and its gradient are reduced near the core center while a little increase is noted near the core boundary, which matches the local power distribution. The change of temperature gradient causes a lower local entropy production near the core center and the one higher near the core boundary, which is shown in Figure 12c.
As a result, the evaluation indexes of the standard problem and optimization scheme are listed in Table 9. An obvious decrease can be seen, which is satisfactory for the optimization.

5. Conclusions

In this paper, a typical PWR reactor core based on the 3D-IAEA problem is taken to be analysed. The neutron flux field, temperature and local entropy production distribution in this reactor core are calculated by FVM. The calculation results are illustrated in cloud pictures in order to observe them clearly. For evaluating the safety and energy dissipation, three characteristics are selected as evaluation indexes, which are the power peak factor Pmax, the maximum temperature Tmax, and the total entropy production Stot. Then the RSM is introduced to optimize the control rod insertion positions in order to get a lower Pmax, Tmax and Stot, which means higher safety and lower energy dissipation. Three directions of optimization are taken, and a final optimization scheme is obtained. The comparison of cloud pictures between the standard problem and optimization scheme shows that the central control rod has a great influence on the power distribution. In total, Pmax, Tmax and Stot are reduced by 23%, 8.7% and 16% after the adjustment of the control rod insert positions, which is satisfactory. RSM may be a quick and useful method in reactor optimization.

Acknowledgments

This work is supported by the Foundation for Innovative Research Groups of the National Natural Science Foundation of China (Grant No. 51421063) and the National Natural Science Foundation of China (Grant No. 51536001).

Author Contributions

Hao-Chun Zhang and Yi-Ning Zhang conceived and designed the research; Yi-Ning Zhang conducted the calculation; Yi-Ning Zhang and Chao Ma analyzed the data; Hai-Yan Yu contributed analysis tools; Yi-Ning Zhang wrote the paper, with the assistance of Hao-Chun Zhang.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Niu, S.; Liu, Y.; Ding, Y.; Qu, W. China’s energy systems transformation and emissions peak. Renew. Sustain. Energy Rev. 2016, 58, 782–795. [Google Scholar] [CrossRef]
  2. Guo, X.; Guo, X. Nuclear power development in China after the restart of new nuclear construction and approval: A system dynamics analysis. Renew. Sustain. Energy Rev. 2016, 57, 999–1007. [Google Scholar] [CrossRef]
  3. Shi, H.Y.; Song, H.T. Applying the real option approach on nuclear power project decision making. Energy Procedia 2013, 39, 193–198. [Google Scholar] [CrossRef]
  4. Weisman, J. Elements of Nuclear Reactor Design, 2nd ed.; Robert, E., Ed.; Krieger Publishing Co. Inc.: Malabar, FL, USA, 1983. [Google Scholar]
  5. Duderstadt, J.J.; Hamilton, L.J. Nuclear Reactor Analysis; John Wiley and Sons Inc.: New York, NY, USA, 1976. [Google Scholar]
  6. Nahla, A.A. Efficient computational system for transient neutron diffusion model via finite difference and theta methods. Ann. Nucl. Energy 2015, 89, 28–37. [Google Scholar] [CrossRef]
  7. Pazirandeh, A.; Jalili, B.M.H.; Abadi, P. Finite difference method for solving neutron diffusion equation in hexagonal geometry. In Proceedings of the International Conference Nuclear Energy for New Europe, Bled, Slovenia, 14–17 September 2009.
  8. Harvie, D.J.E. An implicit finite volume method for arbitrary transport equations. ANZIAM J. 2010, 52, C1126–C1145. [Google Scholar] [CrossRef]
  9. Bernal, A.; Roman, J.E.; Miro, R.; Verdú, G. Assembly discontinuity factors for the neutron diffusion equation discretized with the finite volume method application to BWR. Ann. Nucl. Energy 2016, 97, 76–85. [Google Scholar] [CrossRef]
  10. Theler, G. Unstructured grids and the multigroup neutron diffusion equation. Sci. Technol. Nucl. Instrall. 2013, 2013, 333–339. [Google Scholar] [CrossRef]
  11. Ge, J.; Zhang, D.; Tian, W.; Wang, K.; Qiu, S.; Su, G.H. Steady and transient solutions of neutronics problems based on finite volume method (FVM) with a CFD code. Prog. Nucl. Energy 2015, 85, 366–374. [Google Scholar] [CrossRef]
  12. Hosseini, S.A. Development of galerkin finite element method three-dimensional computational code for the on diffusion equation with unstructured tetrahedron elements. Nucl. Eng. Technol. 2015, 48, 43–54. [Google Scholar] [CrossRef]
  13. Eleodor, N.; Farzad, R. A heterogeneous finite element method in diffusion theory. Ann. Nucl. Energy 2003, 30, 317–347. [Google Scholar]
  14. Vidal-Ferrandiz, A.; Fayez, R.; Ginestar, D.; Verdú, G. Moving meshes to solve the time-dependent neutron diffusion equation in hexagonal geometry. J. Comput. Appl. Math. 2016, 291, 197–208. [Google Scholar] [CrossRef]
  15. Bejan, A. A study of entropy generation in fundamental convective heat transfer. J. Heat Transf. 1979, 101, 718–725. [Google Scholar] [CrossRef]
  16. Bejan, A. Entropy Generation through Heat and Fluid Flow. J. Appl. Mech. 1983, 50, 475. [Google Scholar] [CrossRef]
  17. Bejan, A. Advanced Engineering Thermodynamics, 2nd ed.; Wiley: New York, NY, USA, 1997. [Google Scholar]
  18. Bejan, A. Entropy Generation Minimization; CRC Press: New York, NY, USA, 1996. [Google Scholar]
  19. Ibáñez, G.; Cuevas, S.; Haro, M.L.D. Minimization of entropy generation by asymmetric convective cooling. Int. J. Heat Mass Transf. 2003, 46, 1321–1328. [Google Scholar] [CrossRef]
  20. Makinde, D.; Aziz, A. Analysis of entropy generation and thermal stability in a long hollow cylinder with asymmetry convective cooling. Heat Mass Transf. 2011, 47, 1407–1415. [Google Scholar] [CrossRef]
  21. Aziz, A.; Khan, W.A. Entropy generation in an asymmetrically cooled slab with temperature-dependent internal heat generation. Heat Trans. Res. 2012, 41, 260–271. [Google Scholar] [CrossRef]
  22. Malvandi, A.; Ganji, D.D.; Hedayati, F.; Yousefi Rad, E. An analytical study on entropy generation of nanofluids over a flat plate. Alex. Eng. J. 2013, 52, 595–604. [Google Scholar] [CrossRef]
  23. Torabi, M.; Zhang, K. Classical entropy generation analysis in cooled homogenous and functionally graded material slabs with variation of internal heat generation with temperature, and convective–radiative boundary conditions. Energy 2014, 65, 387–397. [Google Scholar] [CrossRef]
  24. Ivanov, K.; Avramova, M. Challenges in coupled thermal-hydraulics and neutronics simulations for LWR safety analysis. Energy 2007, 34, 501–513. [Google Scholar] [CrossRef]
  25. Sjenitzer, B.L.; Hoogenboom, J.E.; Escalante, J.J.; Espinoza, V.S. Coupling of dynamic Monte Carlo with thermal-hydraulic feedback. Ann. Nucl. Energy 2015, 76, 27–39. [Google Scholar] [CrossRef]
  26. Knoll, D.A.; Keyes, D.E. Jacobian-free Newton-Krylov methods: A survey of approaches and applications. J. Comput. Phys. 2004, 193, 357–397. [Google Scholar] [CrossRef]
  27. Song, H.; Luo, R.; Wan, J.; Li, S.; Zhao, F. Development of a novel 1D coupled neutronics/thermal-hydraulics code and its verification on PWR rod ejection accident benchmark. Prog. Nucl. Energy 2016, 92, 197–210. [Google Scholar] [CrossRef]
  28. Bezerra, M.A.; Santelli, R.E.; Oliveira, E.P.; Villar, L.S.; Escaleira, L.A. Response surface methodology (RSM) as a tool for optimization in analytical chemistry. Talanta 2008, 76, 965–977. [Google Scholar] [CrossRef] [PubMed]
  29. Gilmour, S.G. Response surface designs for experiments in bioprocessing. Biometrics 2006, 62, 323–331. [Google Scholar] [CrossRef] [PubMed]
  30. Kalil, S.J.; Maugeri, F.; Rodrigues, M.I. Response surface analysis and simulation as a tool for bioprocess design and optimization. Process Biochem. 2000, 36, 539–550. [Google Scholar] [CrossRef]
  31. Rajashekhar, M.R.; Ellingwood, B.R. A new look at the response surface approach for reliability analysis. Struct. Saf. 1993, 12, 205–220. [Google Scholar] [CrossRef]
  32. Bucher, C.G.; Bourgund, U. A fast and efficient response surface approach for structural reliability problems. Struct. Saf. 1990, 7, 57–66. [Google Scholar] [CrossRef]
  33. D’Archivio, A.A.; Maggi, M.A. Investigation by response surface methodology of the combined effect of pH and composition of water-methanol mixtures on the stability of curcuminoids. Food Chem. 2017, 219, 414–418. [Google Scholar] [CrossRef] [PubMed]
  34. Espadabellido, E.; Ferreirogonzalez, M.; Carrera, C.; Palma, M.; Barroso, C.G.; Barbero, G.F. Optimization of the ultrasound- assisted extraction of anthocyanins and total phenolic compounds in mulberry (Morus nigra) pulp. Food Chem. 2016, 219, 23–32. [Google Scholar] [CrossRef] [PubMed]
  35. Wang, Y.; Gao, Y.; Ding, H.; Liu, S.; Han, X.; Gui, J.; Liu, D. Subcritical ethanol extraction of flavonoids from Moringa oleifera leaf and evaluation of antioxidant activity. Food Chem. 2017, 218, 152–158. [Google Scholar] [CrossRef] [PubMed]
  36. Argonne Code Center: Benchmark Problem Book. ANL-7416 (Supplement 2). Available online: https://www.osti.gov/scitech/servlets/purl/5037820 (accessed on 8 February 2017).
  37. Stacey, W.M. Nuclear Reactor Physics; Wiley: New York, NY, USA, 2001. [Google Scholar]
  38. Bernal, A.; Roman, J.E.; Miró, R.; Ginestar, D. Development of a finite volume inter-cell polynomial expansion method for the neutron diffusion equation. J. Nucl. Sci. Technol. 2016, 53, 1212–1223. [Google Scholar] [CrossRef]
  39. Chen, W.; Yu, L.; Hao, J. Thermal-Hydraulics of Nuclear Power Facilities; Atomic Press: Beijing, China, 2012. [Google Scholar]
  40. Aslan, N.; Cebeci, Y. Application of Box-Behnken design and response surface methodology for modeling of some Turkish coals. Fuel 2007, 86, 90–97. [Google Scholar] [CrossRef]
  41. Kincl, M.; Turk, S.; Vrecer, F. Application of experimental design methodology in development and optimization of drug release method. Int. J. Pharm. 2005, 291, 39–49. [Google Scholar] [CrossRef] [PubMed]
  42. Lake, J.A. The fourth generation of nuclear power. Prog. Nucl. Energy 2002, 40, 301–307. [Google Scholar] [CrossRef]
Figure 1. Horizontal cross section of the 3D-IAEA problem.
Figure 1. Horizontal cross section of the 3D-IAEA problem.
Entropy 19 00063 g001
Figure 2. Vertical cross section of the 3D-IAEA problem.
Figure 2. Vertical cross section of the 3D-IAEA problem.
Entropy 19 00063 g002
Figure 3. Fast neutron flux at the diagonal line at the level of z = 195 cm.
Figure 3. Fast neutron flux at the diagonal line at the level of z = 195 cm.
Entropy 19 00063 g003
Figure 4. Thermal neutron flux at the diagonal line at the level of z = 195 cm.
Figure 4. Thermal neutron flux at the diagonal line at the level of z = 195 cm.
Entropy 19 00063 g004
Figure 5. Local power distribution of standard problem (left: vertical cross section cloud picture at y = 0; right-top: horizontal cross section cloud picture at z = 315 cm; right-bottom: horizontal cross section cloud picture at z = 195 cm. Pavg represents the average local power of standard problem).
Figure 5. Local power distribution of standard problem (left: vertical cross section cloud picture at y = 0; right-top: horizontal cross section cloud picture at z = 315 cm; right-bottom: horizontal cross section cloud picture at z = 195 cm. Pavg represents the average local power of standard problem).
Entropy 19 00063 g005
Figure 6. Temperature distribution of standard problem (left: vertical cross section cloud picture at y = 0; right-top: horizontal cross section cloud picture at z = 315 cm; right-bottom: horizontal cross section cloud picture at z = 195 cm).
Figure 6. Temperature distribution of standard problem (left: vertical cross section cloud picture at y = 0; right-top: horizontal cross section cloud picture at z = 315 cm; right-bottom: horizontal cross section cloud picture at z = 195 cm).
Entropy 19 00063 g006
Figure 7. Local entropy production distribution of standard problem (left: vertical cross section cloud picture at y = 0; right-top: horizontal cross section cloud picture at z = 315 cm; right-bottom: horizontal cross section cloud picture at z = 195 cm. S0 represents the average local entropy production of standard problem).
Figure 7. Local entropy production distribution of standard problem (left: vertical cross section cloud picture at y = 0; right-top: horizontal cross section cloud picture at z = 315 cm; right-bottom: horizontal cross section cloud picture at z = 195 cm. S0 represents the average local entropy production of standard problem).
Entropy 19 00063 g007
Figure 8. Flow diagram of optimization procedure.
Figure 8. Flow diagram of optimization procedure.
Entropy 19 00063 g008
Figure 9. Local power distribution of optimization scheme (left: vertical cross section cloud picture at y = 0; right-top: horizontal cross section cloud picture at z = 315 cm; right-bottom: horizontal cross section cloud picture at z = 195 cm. Pavg represents the average local power of standard problem).
Figure 9. Local power distribution of optimization scheme (left: vertical cross section cloud picture at y = 0; right-top: horizontal cross section cloud picture at z = 315 cm; right-bottom: horizontal cross section cloud picture at z = 195 cm. Pavg represents the average local power of standard problem).
Entropy 19 00063 g009
Figure 10. Temperature distribution of optimization scheme (left: vertical cross section cloud picture at y = 0; right-top: horizontal cross section cloud picture at z = 315 cm; right-bottom: horizontal cross section cloud picture at z = 195 cm).
Figure 10. Temperature distribution of optimization scheme (left: vertical cross section cloud picture at y = 0; right-top: horizontal cross section cloud picture at z = 315 cm; right-bottom: horizontal cross section cloud picture at z = 195 cm).
Entropy 19 00063 g010
Figure 11. Local entropy production distribution of optimization scheme (left: vertical cross section cloud picture at y = 0; right-top: horizontal cross section cloud picture at z = 315 cm; right-bottom: horizontal cross section cloud picture at z = 195 cm. S0 represents the average local entropy production of standard problem).
Figure 11. Local entropy production distribution of optimization scheme (left: vertical cross section cloud picture at y = 0; right-top: horizontal cross section cloud picture at z = 315 cm; right-bottom: horizontal cross section cloud picture at z = 195 cm. S0 represents the average local entropy production of standard problem).
Entropy 19 00063 g011
Figure 12. Local power, temperature and local entropy production of the standard problem and optimization scheme at the diagonal line on the midplane (the level of z = 195 cm).
Figure 12. Local power, temperature and local entropy production of the standard problem and optimization scheme at the diagonal line on the midplane (the level of z = 195 cm).
Entropy 19 00063 g012
Table 1. Two-group nuclear data.
Table 1. Two-group nuclear data.
RegionD1D2Σs,1→2Σa1Σa2υΣfMaterial
11.50.40.020.010.080.135Fuel 1
21.50.40.020.010.0850.135Fuel 1
31.50.40.020.010.130.135Fuel 2 + Rod
42.00.30.0400.010Reflector
52.00.30.0400.0550Refl. + Rod
Table 2. Additional thermal parameters.
Table 2. Additional thermal parameters.
ParameterPcore/MWλ/W·m−1·K−1Tboundary/K
λ1λ2λ3λ4λ5
Value1055.255.2555400
Table 3. Results on different mesh quantity.
Table 3. Results on different mesh quantity.
No.Mesh QuantitykeffError
117 × 17 × 19 (5491)1.02904-
234 × 34 × 38 (43,928)1.028550.048%
351 × 51 × 57 (148,257)1.028670.012%
468 × 68 × 76 (351,424)1.028730.006%
Table 4. Calculation results of standard problem.
Table 4. Calculation results of standard problem.
ParameterPmaxTmaxStot
Value2.5107657.40.8021
Table 5. Numerical experiments used to design response surface.
Table 5. Numerical experiments used to design response surface.
No.PositionsCalculation ResultsNo.PositionsCalculation Results
Z1Z2Z3Z4PmaxTmaxStotZ1Z2Z3Z4PmaxTmaxStot
1360190201903.258681.40.8159161901903603602.5946580.7794
21901901901903.24690.70.8615173601903601902.808673.10.8231
3201901903602.445641.40.760718190190360203.228688.60.8565
4360201901902.2116330.714119190190203602.889658.60.7531
519020190202.696649.80.749120190360190203.257694.40.8805
61901901901903.24690.70.86152119019020203.919703.40.8588
719020201902.052610.40.6448221903603601902.511663.80.8186
81901901901903.24690.70.861523201903601902.328639.50.7644
91903601903602.624663.70.7999241901901901903.24690.70.8615
10190201903602.071619.40.6795253601901903602.9386760.8196
1120190190202.919668.80.835526190360201902.89669.40.8045
12360190190203.649709.80.9074273603601901902.839680.20.8497
13190203601902.026622.40.7011281901901901903.24690.70.8615
14203601901902.332645.70.78422920201901902.041595.60.6406
1520190201902.486637.90.7365--------
Table 6. Coefficients in the response surface.
Table 6. Coefficients in the response surface.
Coefficientb1b2b3Coefficientb1b2b3
β01.825872573.34650.564135β23−3.10 × 10−6−1.50 × 10−4−3.60 × 10−7
β10.0057680.3749310.000685β24−7.00 × 10−8−2.20 × 10−6−9.50 × 10−8
β20.0100370.6288990.001546Β343.42 × 10−61.23 × 10−42.48 × 10−7
β30.0021260.2073160.000692β11−1.00 × 10−5−6.10 × 10−4−1.10 × 10−6
β4−0.00324−0.1222−0.00021β22−2.20 × 10−5−1.19 × 10−4−2.80 × 10−6
β122.91 × 10−6−2.50 × 10−5−6.90 × 10−8β33−7.10 × 10−6−5.10 × 10−4−1.50 × 10−6
β13−2.50 × 10−6−8.60 × 10−5−1.80 × 10−7β442.59 × 10−63.41 × 10−5−1.10 × 10−7
β14−2.10 × 10−6−5.60 × 10−5−1.10 × 10−7----
Table 7. R-square of response surface fitting.
Table 7. R-square of response surface fitting.
ResponsePmaxTmaxStot
R-square0.94950.99050.9892
Table 8. The selected control rod positions.
Table 8. The selected control rod positions.
No.PositionsStotPmaxTmax
1[30 30 300 300]0.67031.9408600.47
2[30 40 300 300]0.67071.9418600.58
3[50 30 300 300]0.67091.9433600.64
4[40 30 300 300]0.67101.9439600.66
5[30 30 290 300]0.67241.9602601.19
6[40 30 300 290]0.67311.9550601.25
7[30 30 300 290]0.67331.9561601.30
8[30 30 300 280]0.67521.9674601.86
Table 9. Comparison between calculation results.
Table 9. Comparison between calculation results.
ParameterStandard ProblemOptimization SchemeDecrease
Pmax2.51071.940823%
Tmax657.4600.478.7%
Stot0.80210.670316%

Share and Cite

MDPI and ACS Style

Zhang, Y.-N.; Zhang, H.-C.; Yu, H.-Y.; Ma, C. Response Surface Methodology Control Rod Position Optimization of a Pressurized Water Reactor Core Considering Both High Safety and Low Energy Dissipation. Entropy 2017, 19, 63. https://doi.org/10.3390/e19020063

AMA Style

Zhang Y-N, Zhang H-C, Yu H-Y, Ma C. Response Surface Methodology Control Rod Position Optimization of a Pressurized Water Reactor Core Considering Both High Safety and Low Energy Dissipation. Entropy. 2017; 19(2):63. https://doi.org/10.3390/e19020063

Chicago/Turabian Style

Zhang, Yi-Ning, Hao-Chun Zhang, Hai-Yan Yu, and Chao Ma. 2017. "Response Surface Methodology Control Rod Position Optimization of a Pressurized Water Reactor Core Considering Both High Safety and Low Energy Dissipation" Entropy 19, no. 2: 63. https://doi.org/10.3390/e19020063

APA Style

Zhang, Y. -N., Zhang, H. -C., Yu, H. -Y., & Ma, C. (2017). Response Surface Methodology Control Rod Position Optimization of a Pressurized Water Reactor Core Considering Both High Safety and Low Energy Dissipation. Entropy, 19(2), 63. https://doi.org/10.3390/e19020063

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop