A review of pebble flow study for pebble bed high temperature gas-cooled reactor

The pebble bed high temperature gas-cooled reactor is a promising generation-IV reactor, which uses large fuel pebbles and helium gas as coolant. The pebble bed flow is a fundamental issue for both academic investigation and engineering application, e.g., reactor core design and safety analysis. This work performed a review of recent progress on pebble flow study, focusing on the important issues like pebble flow, gas phase hydrodynamics, and inter-phase heat transfer (thermal hydraulics). Our group’s researches on pebble flow have also been reviewed through the aspects of phenomenological observation and measurement, voidage distribution, geometric and parameter optimization, pebble flow mechanisms, flow regime categorization, and fundamentals of modelings of pebble flow and radiation. Finally, the major problems or possible directions of research are concluded which would be some of our focuses on the pebble bed flow study.


Pebble bed high temperature gas-cooled reactor
As well known, the pebble bed reactors (PBR) and prismatic block reactors (PMR) are the two main types of high temperature gas-cooled reactors (HTGRs). The traditional PBR is a graphite-moderated, helium-cooled reactor. The reactor core consists of spherical fuel elements called pebbles. HTR-10 ( Fig. 1) is a test pebble bed real reactor built at the Institute of Nuclear and New Energy Technology (INET) of Tsinghua University (Hao et al., 2018). It uses 6-cm-diameter pebbles (Fig. 2), which are made of pyrolytic graphite (Kadak, 2005). Each fuel element contains approximate 12,000 microscopic tristructural-isotropic (TRISO) coated fuel particles of 0.92 mm diameter. The TRISO particles are scattered in a sphere of graphite matrix of a diameter f D = 5.0 cm to form the fuel zone, which is surrounded by a 5 mm thick fuel free zone also composed of graphite matrix. The UO 2 kernel of about k D = 0.5 mm is surrounded by a TRISO coating consisting of the buffer pyrocarbon layer, the inner high-density pyrocarbon layer (IPyC), the silicon carbon (SiC) layer, and the outer high-density pyrocarbon layer (OPyC). It was proved that, irradiated under the HTR-Module typical operation conditions, no additional particle failure was found when fuel elements were heated up to 1620 °C. Therefore, the temperature limitation of the HTR-PM fuel elements during design basis accident (DBA) is set at 1620 °C (Chen et al., 2017b). The high temperature gas-cooled reactor (HTGR) is featured by inherent safety with low power density and large graphite core design, in combination with environmental applicability, high efficiency, and industrial process heat applied in producing hydrogen. Therefore, it is generally recognized as a probable solution for the generation IV advanced reactors (Ryskamp et al., 2004).
Recently, a commercial-scale 200 MWe high temperature gas-cooled reactor pebble bed module project (HTR-PM, Fig. 3 (Sun et al., 2018)) is now under construction at Shandong Province in China. The HTR-PM nuclear power plant has two pebble-bed modular reactors, which are   (Kadak, 2005; reproduced with permission © Inderscience Enterprises Ltd. 2005). connected to one steam-turbine generator. The total power of them is 500 MWt, and the electrical generating efficiency is about 42%. The reactor core of HTR-PM is a randomly packed pebble bed of 3.0 m in diameter and 11 m in height on average, containing roughly 420,000 spherical fuel elements with a diameter of 6.0 cm . Thus, the pebble flow characteristics could be of great importance for the HTR-PM.

Pebble flow
In general, a large number of pebbles in reactor core of HTGR are randomly piled up in the core region. The pebble-bed HTGR runs in a recirculating way. When the fuel pebbles are drained out from discharge hole at the bottom of the core regularly and individually and loaded into the core from the top. When pebble bed reaches the equilibrium state, the number of pebbles in the core approximately remains constant (Yang et al., 2014b).
The basic physics of dense granular flow, to say precisely, the extremely slow pebble flow, is not fully understood yet. Different from other types of granular flow, the flow rate of pebbles in the reactor is so slow that the equivalent flow velocity is lower than most granular flows in several orders. In most applications, it can be viewed as a pebble system composed of equivalent stagnant pebbles. Therefore, we named it as a new flow regime of extremely slow pebble flow. For mechanical study of the dynamics of discrete particles, the local packing structure and equivalent stresstensor may be of the most importance. The packing structure is needed for either the design optimization or safety assessment, such as the packing density, radial distribution function, and coordination number distribution (Zhao et al., 2015). The packing factor is also important since it affects the mechanical behavior , void distribution and coolant flow drag and pressure drop, especially for multisized pebble bed (Reimann et al., 2013). For example, the mechanical behavior of mixed fusion pebble beds may have important features related to the stiffness between the softersolid system and stiffer-solids system, which is dependent on both solid material and the loading-unloading cycle numbers (see Fig. 4  ). On the other hand, the packing factors for pebbles in a piston-pressured container are fairly independent of bed height unless the height to diameter ratio becomes less than 10 (Reimann et al., 2013).
Besides, the mechanical and thermal properties of pebble beds have coupled effects, and the thermo-mechanical characteristics are critical for the reactor design especially for that uses a solid breeder (Lo Frano et al., 2016;Zheng et al., 2016). For example, the compression force and gas pressure were tested (Tehranian and Abdou, 1995), along with the thermal conductivity at different values of temperature (Lo Frano et al., 2014). The strain and stress analysis by finite element method, as well as their effects on the thermal hydraulics were also checked, although the coupled effects were sometimes not found (Zaccari and Aquaro, 2007). Moreover, the effective thermal conductivity was found not to be influenced by the chemical composition of solid material whereas it can be affected by filling gas such as helium clearly (Pupeschi et al., 2017).
In reactor engineering, the pebble path through the reactor core, the resident time in the reactor core, and the optimization of fuel composition (Ho and Obara, 2016) are of crucial importance. Excessive resident time inside the bed may cause severe irradiation and thermal damage to the pebble, as well as possible escape of fission products. In an extreme case, it may lead to the permanent adhesion between pebbles (Zhao et al., 2015). On the other hand, the resident time as well as axial dispersion coefficient can be used to characterize the axial dispersion and the mixing phenomena of the coolant gas flow, which are useful for safe design and efficient operation of these reactors (Abdulmohsin and Al-Dahhan, 2016). For example, the pebble flow dynamics were studied in a scaled down test reactor by using a non-invasive radioactive particle tracking (RPT) technique (Khane et al., 2016a). They used a Cobalt-60 based tracer to mimic pebbles in terms of shape, size, and density. A crosscorrelation based position reconstruction algorithm and RPT calibration data were used to obtain the Lagrangian trajectories, velocity field, and residence time distributions (Fig. 5). Moreover, the pebble size was found to strongly affect axial dispersion and mixing in the packed pebble-bed reactor (Abdulmohsin and Al-Dahhan, 2016). Accurate flow characterization and description of the dynamics of pebbles in pebble beds are important with respect to the basic reactor design calculations, the optimization of fuel cycle, and the burn-up calculations, as well as the monitoring of fuel integrity over its lifetime (Zhao et al., 2015). In addition, the pebble recirculation in pebble beds depends on the maintenance of satisfactory pebble flow through the vessel, its outlet, and the pebble extractor (Zhao et al., 2015). In particular, the pebble flow pattern in the reactor core is expected to be uniform and consistent (a mass flow regime). The flooding of pebbles, channeling, and flow blockage phenomena should be prevented, and the caking or crystallization of pebbles and the formation of stagnant regions need to be avoided (Li et al., 2016a). Therefore, some kinds of wall structures are designed to affect the overall flow field through avoiding or eliminating the crystallization of pebbles in the near-wall region. This can contribute to the elimination of stagnant region and improve reactor's safety capacity by reducing the probability of radiation leakage. Moreover, such structures also enhance the dispersion degree of pebbles in the peripheral region of the core, so that these fuel pebbles are more likely to be exposed to more neutron flux, which would narrow the burnup level in radial direction . Moreover, the pebble motion in the reactor core is expected to follow "first-in-first-out" sequence to ensure each fuel pebble can reach almost the same burnup level when it is drained out (Li et al., 2016a). In addition, Rycroft et al. (2013) studied the scaling properties of granular flow in pebble-bed reactors via DEM simulation in a supercomputing Centreon Rosa, a Cray XE6 system with 47,872 cores, where the pebbles are scaled down by 3:1 and 6:1. The simulation results showed the feature of pebble stresses and dust generation due to pebble wear (see Fig. 6), e.g., the 6:1 simulation exhibits a very different feature with higher levels of wear, particularly pebble-wall wear, and lower stresses due to the Janssen effect.
Besides, the two-region pebble bed reactor core is expected as a promising technique for pebble bed HTGRs (IAEA, 2001). The two-region pebble bed is composed of two distinct regions, i.e., a central column region consisting of graphite pebbles (called graphite region) and an outer annular region consisting of fuel pebbles (called fuel region). Therefore, the graphite pebbles are loaded into the core from the single central hole and fuel pebbles are loaded from the annular periphery of the core. The two-region pebble bed is an advantageous type as it flattens the neutron flux and consequently allows a significantly higher power output without reducing safety margin. The decay heat also transfers a shorter distance from the core to outside during accidents (Yang et al., 2014b). In addition, thermal analysis of the pebble beds is also of vital importance for the reliable blanket design, e.g., the helium cooled solid tritium breeder blanket in Chinese fusion engineering test reactor (Zhou et al., 2015), which uses the pebbles of lithium ceramics (Li 4 SiO 4 ) and beryllium as tritium breeder and neutron multiplier, respectively.
The numerical algorithms for pebble flow study can be divided into two categories. One is based on the Lagrangian framework to describe the motion of every discrete pebble, and the other one is based on the Eulerian framework to simulate the flow features in analogies to continuum theory. One of the most prevalent Lagrangian approach is the discrete element method and its various extensions Yang et al., 2014a;Li et al., 2015Li et al., , 2016aLi et al., , 2016bLi et al., , 2016cJia et al., 2017a;Al Falahi et al., 2018). For example, numerical simulation of large scale pebble flow by DEM code have been performed (Rycroft et al., 2006(Rycroft et al., , 2013Li et al., 2009;Sun et al., 2017). Particular attentions have been paid to the streamlines, velocity distributions and pebble diffusion characteristics.
The representative approach in the Eulerian framework is the kinematic model (Nedderman and Tüzün, 1979). As it was difficult to apply the kinematic model to analyze pebble velocity for reactors having complex geometries, e.g., the PBMR with an annular core and three defuel chutes, Kim et al. (2013) modified it to improve the reconstruction ability of the pebble velocity profile by using cylindrical core experiments to determine the specific coefficients used in this modified kinematic model and verifying the model by the pebble flow experiment of PBMR core. But, these coefficients are still not widely applicable for other reactor design.

Coolant material and flow
The coolant of pebble beds includes helium, nitrogen, supercritical carbon dioxide (S-CO2) (Latifi and du Toit, 2019), etc. The high temperature gas-cooled reactor (HTGR) is based on a Brayton cycle (e.g., three-shaft recuperative Brayton cycle for the Pebble Bed Micro Model (Greyvenstein et al., 2003)) with helium as gas coolant (Nicholls, 2000). Compared to helium, the supercritical carbon dioxide can be used in both direct and indirect cycles, and it is regarded as a suitable coolant for pebble bed reactor due to its large mass density and heat transfer characteristics (Latifi and du Toit, 2019).
Up to date, most of the experimental studies on the hydrodynamics of coolant were restricted to the effect of operating conditions on the global parameters, such as pressure drop, helium mass flow rates, pebble-bed channel lengths, ball diameters , and overall bed voidage (Abdulmohsin and Al-Dahhan, 2016). For example, it showed that the helium flow pressure drop in pebble bed obeyed the flow characteristics in porous medium, i.e., the helium flow pressure drop was proportional to the channel length . The pressure drop of coolant flow is found to be decreased as the ball diameter is increased in the pebble bed  or vice versa (Abou-Sena et al., 2014), and it can use the Ergun's equation to predict with the pressure drop as follows well (Abou-Sena et al., 2014).
where P D is the pressure drop of pebble bed, ε is the pebble bed porosity, μ is the gas dynamic viscosity,  is the pebble sphericity, ρ is the gas density, U is the superficial velocity, L is the bed length, and p D is the pebble diameter. The fluid flow passing through the gaps between the randomly packed pebbles has complex natures, whereas few studies on the flow field of pebble beds have been carried out. Therefore, the hydrodynamic phenomena are not yet well understood. Amongst these studies, particle image velocimetry (PIV) (Hassan and Dominguez-Ontiveros, 2008;Lee and Lee, 2009) and the matched index of refraction (MIR) technique (Hassan and Dominguez-Ontiveros, 2008) are always utilized to measure the liquid-phase flow field of the interior region of a small packed bed. Similarly, fluorescent particle image velocimetry technique in conjunction with matched index of refraction (MIR) was applied by Northrup et al. (1991Northrup et al. ( , 1993 to study the intrapore mixing within porous medium. Both LES and the particle tracking velocimetry (PTV) and MIR technique experimental results showed very complex flow structures and vortical structures within the gaps between the fuel elements (Hassan and Dominguez-Ontiveros, 2008). Experimental measurement of the axial dispersion coefficients of the gas phase and their residence time distributions (Abdulmohsin and Al-Dahhan, 2016) indicates that the flow pattern of the gas phase does not deviate much from the idealized plug-flow condition at high flow rate, which depends on the gas flow rate and the bed structure of the pebble bed. With regard to the bypass flow in near-wall gaps, Amini and Hassan (2014) utilized a hot wire anemometry (HWA) system coupled with a hot film X-probe to measure the radial and axial velocity components of the coolant flow passing through gaps close to the outer wall of a test section modeling an annular pebble bed reactor.
The coolant flow characteristics can be studied by computational fluid dynamical (CFD) simulation (e.g., k ε model in ANSYS fluent Wu et al., 2018d), CFD-DEM coupled method (Li and Ji, 2013;Chen et al., 2016Chen et al., , 2017a, direct numerical simulation (DNS) (Shams et al., 2013b) or large eddy simulation (LES) (Ebara et al., 2010)) to study the interstitial flow between pebbles, or by thermalfluid network analysis code (e.g., Flownex (Venter and Lamprecht, 2012). In general, ANSYS fluent or most CFD codes, in combination with suitable RANS (Reynolds Averaged Navier-Stokes) and LES models, are suitable to simulate the local details of flow, turbulence and heat transfer, and limited scales whole thermal-hydraulic system, supposing the huge numbers of grids are allowed by computer-clusters or super-computers. Detailed comparison of LES, Hybrid (RANS/LES) and RANS models with the reference q-DNS were performed for a well-defined single face cubic centered pebble configuration (Shams et al., 2012(Shams et al., , 2013a(Shams et al., , 2013b(Shams et al., , 2013c(Shams et al., , 2015. To optimize a pebble bed configuration, DNS is a good method to provide a reference for the validation of different turbulence modeling approaches, such as results for mean, RMS and covariance of velocity field and results related to temperature field of two different cross-sections through the cubic centered pebble structure and 43 profiles at different locations in the computational domain (see Fig. 7) (Shams et al., 2013b(Shams et al., , 2013c. But the highly complex physics of three-dimensional flow behavior between packed pebbles still make it a challenging work for turbulence model validation. Besides, body-centered cubic (BCC) and face-centered cubic (FCC) structures are also simulated to predict higher heat transfer capability and lower pebble temperature. The predicted average Nussel (Nu) number was found to decrease from the first layer and reach an asymptotic value when the gas passes through the sixth layer of pebbles (Fig. 8) (Ferng and Lin, 2013). In a similar manner, Shams et al. (2012) applied also an FCC geometry which contains full spheres at each corners and face centre of the cube. In this computation, the single cubic FCC domain is composed of half spherical pebble at each face and 1/8 at each corner (see Fig. 9). Similar work on CFD simulation of coolant flow within a BCC geometrical model for a pebble bed water cooled reactor was performed by using the standard k ε turbulence model (see Fig. 10  ). In more details, Lee et al. (2007a) studied two kinds of inter-pebble gap and two kinds of direct contact between pebbles, and found that the direct contact between pebbles shows numerous differences in the results of the flow regime around the pebbles as well as in the wake, compared to the cases of the inter-pebble gap. However, no large differences were found between two cases of direct contact. For pebble bed thermal hydraulics, the thermal field heat transfer characteristics, and effective thermal conductivity (k eff ) for the uniform-size (Bauer and Schluender, 1978) or multiple-size (Ades and Peddicord, 1982;Chen et al., 2015) pebble beds are important issues (Chen and Lee, 2017). In particular, it is critical to maintain the reactor core within the thermal limit in the coolant loss accident (Lohnert, 1990). The effect thermal conductivity determines the passive safety characteristic of removal of decay heat via heat conduction from the reactor core to the environment after the reactor is shut down. It is extremely important to identify the local hot spot locations which relies on the understanding of thermodynamics and field distribution inside thereactor core (Chen and Lee, 2017). Therefore, heat transfer characteristics of regular structures, e.g., simple cubic (SC), body-centered cubic (BCC), and face-centered cubic (FCC) structures as well as their combinations, have been studied via both experimental measurement (Chen and Lee, 2017) and CFD simulations (Yesilyurt and Hassan, 2003;Lee et al., 2007aLee et al., , 2007bLaguerre et al., 2008;Kim et al., 2009;Sobes et al., 2011;Li et al., 2012;Shams et al., 2012Shams et al., , 2013cShams et al., , 2015Ferng and Lin, 2013). A test section was designed to measure the    the destruction of a pebble's integrity. The three-dimensional and two-dimensional views of it can be seen in Fig. 11 (Chen and Lee, 2017). It used 2 full spheres, 4 hemispheres, and 8 quarter spheres packed inside in an FCC structure to ensure that they are in mutual contact with each other. It was found that areas with φ = 36 °C and 117 °C from the z-axis to the hole are the strongest heat transfer zones, while areas with φ = 0 °C, 90 °C, and 180 °C are the weakest heat transfer zones. The heat transfer coefficients are correlated with the Nusselt number and the Reynolds number in the form of

= .
, respectively. But these data points still seem to be not enough.
The local convective heat transfer coefficients as well as radial profiles were measured experimentally using a sophisticated noninvasive heat transfer probe (Abdulmohsin and Al-Dahhan, 2015). The effect of gas velocity on the heat transfer coefficient was investigated over a wide range of Reynolds numbers. The local heat transfer coefficient increases from the bed center to the wall due to the change in the bed structure and the gas flow pattern. The effective thermal conductivity is an overall heat transfer parameter of a packed bed of spheres which is crucial in the analysis and design of pebble bed gas-cooled reactors. During depressurized loss of coolant, the dominant heat transfer mechanisms for the passive removal of decay heat are radiation and conduction. Predicting the value of the effective thermal conductivity is complex since it depends on the temperature level and temperature gradient through the bed, as well as the pebble packing structure (de Beer et al., 2017). De Beer et al. (2017) showed an accompanying methodology to combine physical measurements with computational fluid dynamics to separate the contributions of radiation and conduction heat transfer.
For simulation of coolant flow over the whole bed, the drag force and pressure drop of the packed pebbles must be the key issues, where existing various correlations obtained from measuring the flows past packed or fluidized beds of powders and granules were always used analogously. Then, experimental measurement on pressure drop is of crucial importance (Abou-Sena et al., 2013). On the industrial scale, the VSOP code can approximate the pebble bed flow under limited computer capability, and the VSOP treatment of the pebble bed flow model is reasonable. Statistics analysis shows that several millions of Monte Carlo simulation can provide reasonable result and new code based on Monte Carlo method was also developed to analyze the characteristic of the pebble bed flow (Chen and Fu, 2014). In addition, the flownex is regarded as a general simulation tool that can solve steady-state and transient flows, and pressure and temperature distribution in large-scale arbitrary-structured thermal-fluid networks). It is capable to handle a wide variety of network components, such as pipes, pumps, orifices, heat exchangers, compressors, turbines, controllers, and valves (Flownex, 2010).
Besides, the pneumatic transportation of absorber pebbles in the pipes outside the core of pebble-bed reactor was also simulated by CFD methods to analyze the force and motion of pebbles . The bypass flow through the gaps among graphite reflectors, which may be interacting with the main flow through the pebble bed, was simulated by the computational fluid dynamic (CFD) to clarify the resistance coefficients (Sun et al., 2018). In the full power operation condition, the bypass flow rate ratio was estimated between 0.33% and 1.94% of the total helium flow rate when the gap was consistent with 1.6 mm in width. As the bypass flow takes place in a very thin vertical gaps (see Fig. 12 (Sun et al., 2018)), it is also an important fundamental issue for the area of fluid mechanics. The friction may be a dominating issue where the turbulence development might be suppressed significantly. Therefore, some kinds of modeling methods as well as direct simulation techniques should be developed.

Operational and geometric issues
The optimization of operation of pebble bed always aims at reducing the maximum fuel temperature and to increase fuel performance. To this regard, power profile of pebbles can be influenced by the pebble bed design and loading/ discharge strategy. For pebble bed design, two-region pebble bed , or two-graphite-reflector pebble bed (Koster et al., 2003) can modify the peaks of radial power profiles and the temperature of pebbles at outlet. For example, the multi-zone design can result in a reduction of the maximum fuel temperature of 80 °C and 300 °C for normal operation and depressurized loss of coolant accident conditions (Boer et al., 2009).
On the other hand, the loading and discharge strategy can influence the power profile. In the pebble bed, the pebbles are loaded at the top of the core and move downward by gravity. The pebbles can be discharged after one cycle, called the once-through-then-out (OTTO) scheme (Hansen et al., 1972), and also be reloaded at the top depending on the burnup level, called the cyclic scheme, such as the HTR-10 reactor. For the multiple recycling scheme, the state of equilibrium core composition is different from the initial core composition of pebbles with different fuel content. The distribution of inlet position of pebbles can influence the radial distribution of nuclide in the core, such as the breeder pebbles loaded at the center position in the AVR reactor (Bäumer, 1990). It is also possible to modify the radial power profile by placing fresh pebbles in the outer region of the core and omitting the use of pebbles with different enrichments or burnable poison (Kloosterman, 2003). In addition, using the multiple recycling scheme can also modify the axial power profile (Boer et al., 2009). DEM simulations of pebble flow in reactor core under earthquake condition were analyzed for safety analysis (Keppler, 2013;Gui et al., 2016e).
As the static friction between the wall and the particles increases, there is significantly higher head independence of pressure and flow rate of granular flows in beds (McCabe et al., 1985;Luo et al., 2010). It was found that a bed-heightto-bed-diameter ratio (H/D ≥ 2) will yield long duration of radioactive particle tracking (RPT) experiments. At H/D = 1, the bottom cone angle was found to have significant influence on the flow of pebbles and presence of dead zones (Khane et al., 2016b). During evolution of the exit control mechanism, it was found that the bottom cone angle also affects the jamming of pebbles in the bottom section. A half-cone angle of 60 degrees was found to be less prone to jamming problems (Khane et al., 2016b;Jia et al., 2017b).

Recent researches on pebble flow by our group
As a whole, our group has performed phenomenological measurement, analysis, numerical simulation, and modeling of pebble flow and pebble radiation. On one hand, we have tried to improve the DEM/DEM-CFD based methods to well predict the pebble flow and heat transfer in our test experimental pebble flow bed, as well as the HTR-10 reactor, and the HTR-PM nuclear power plant. On the other hand, we have used the experimental data to validate the proposed various models.
Our interests mainly include:  Phenomenological observation and PTV measurement. A test bed (a pseudo-2D test pebble bed with about 10d p in depth) was built to observe the pebble flow within a 2 m and 1 m high pebble bed respectively (see Fig. 13) Yang et al., 2012). The stagnant time, stagnant time ration, and streamlines of pebbles, as well as the diffusion of pebble around the streamlines are analyzed pheomenologically, which serve as a basis for the numerical simulation and further mechanical analysis (see Fig. 14). To quantify the pebble flow, the velocity profiles were all computed through the PTV technique (see Fig. 15) (Jia et al., 2016).


The distribution of void fraction (voidage) in the 3D pebble bed (Yang et al., 2014a) and near the wall (Yin et al., 2018). The distribution of void fraction is either averaged over the horizontal plane to show the distribution vertically, or averaged circumferentially to show the radial distribution. The main cylindrical body and the conical base of the bed are separately analyzed but combined together to get the statistical empirical correlations.


The geometrical optimization.
(1) Optimize the geometrical base . To make the pebble flow more uniform, the various arc shapes and the Brachistochrone curve were used to form bed base. After comparison, it was shown by both phenomenological comparison and quantitative evaluation (through the recyling and stagnant rates) that the Brachistochrone shape is the best. Otherwise the R 2 arc is the closest shape to the Brachistochrone, followed by the direct line, i.e. (Fig. 16): (2) Optimize the base angle of the bed. The beds of 30°, 45°, and 60° are compared to see which one is the best for reducing or eliminating the stagnant region. The total mass level follows the sequence of 2 3 1 60 45 30 (See in total mass flow level for all bed configurations in Table 1 and Fig. 17) (Jia et al., 2017b). (3) Optimize the structure effects on the cylindrical main body to avoid the formation of stagnant (dead) zone of pebble flow in the bed corner . Some various geometric triangular, arc, and sawtooth shapes were used to see which shape and what parameters are the best for getting a uniform pebble flow. It was shown that, in the two-dimensional case, the right angled triangular shape of d p in short side is the best option among the investigated configurations since it reaches the lowest rate of stagnation. (4) Optimize the pebble flow pattern inside and decrease the stagnant zone by a flowcorrective insert. It is indicated that the key of insert design is to balance the flows inside and outside the insert. A proper insert can restrict central flow and promote side flows to create a relatively uniform flow field. The outlet diameter of insert is the critical parameter that determines the balancing mechanism. The insert should be located close to the discharge hole and its inlet should not exceed the hopperbin junction surface too much to avoid adverse effect on upper flow field (see Fig. 18).    Parameter optimization of pebble flows. Can particleparticle and particle-wall frictional coefficients affect the flow uniformity of pebbles? Yes, we have demonstrated the cross-sectional velocity distributions under different particle-wall and particle-particle frictions (Figs. 19(a) and 19(b)). In particular, which part is more important to get uniform pebble flow? The cylindrical body or the conical base? It is shown that the former is less important than the latter, i.e., the pebble-to-wall friction coefficient in the conical base is of crucial importance to affect the flow regime and geometric shape of pebble stripes (Fig. 19(c)) Li et al., 2016b).  (1) When the pebbles move from the some start point on the bed top to the bottom silo, their trajectories are a bit diffusive, where we name all of them a pebble spindle ( Fig. 20(a)). Then pebble diffusion in the horizontal direction can be evaluated by the probability density function. As it is a extremely slow flow, the diffusion probability can be well predicted by the Pascal trangle ( Fig. 20(b)). Moreover, the high value zone (> 95%) of the Pascal trangle is a long narrow spingle-like region, which explains the diffusion mechanism well (Fig. 20(c)). (2) The flow intermittency and correlation of pebble flow: the pebble flow seems intermittent in time with sometimes high large values of mean force in combination with large scale bulk movement of pebbles inside (Fig. 21)  . The arc connection of pebbles inside the bed was analyzed. Arch size spans from two to six in the detection window with 8d p wide and 6d p high. The arch size distribution obeys the second-order polynomial distribution in the semi-logarithmic scale for all cases. The probability distribution function of the lifetime of arches was also computed. The arching particles show a higher correlation with the mean velocity of total particles too. Not only the blocking arches but also bulk arches should be responsible for the velocity fluctuations. Additional parameters, the correlation time τ c and the intermittency index C 2 were proposed based on the multifractal analysis with the WTMM method (Jia et al., 2016(Jia et al., , 2017c.  Flow regime of pebble flows. We have tried to categorize the pebble flow regime through three view-points.
(1) Correlation of the velocity (kinematic viewpoint) and force (dynamic viewpoint). The correlation analysis indicates that gravity-driven granular flow can be characterized into two classes: the kinematic flow regime (fast dense flow) and the kinetic flow (slow granular flow). In the former, the flow is fluid-like, dynamically stationary, with flow features dominated by kinematic variables, e.g., velocities. Meanwhile, there is significant correlation between mean force and velocity. In the latter, it is a slow quasi-static considerably intermittent flow, with internal sudden "bulk" motion and sudden change in structures. This implies that the transition from slow to fast regime can be characterized by the transition of characteristic frequency of kinetic variable (e.g., contact force) to kinematic variable (e.g., velocity) (see Fig. 22). (2) Using the standard deviation σ and the span of energy magnitude SOE (see Table 2 and Fig. 23), it is also feasible to categorize the pebble flow in the Fig. 19 Distribution of cross-sectional vertical velocity under different particle-particle (a) and particle-wall (b) friction coefficients, and (c) the pebble stripes in the beds of different particle-cylindrical-body-wall friction and particle-conical-base-wall coefficients. three regimes: the extremely slow flow (σ < 0.1, SOE > 1), the consistent flow (σ < 0.1, SOE < 1), and the intermittent flow (σ > 0.1). (3) Intermittency index. It is defined based on the multiplicative cascade method (MCM) (Jia et al., 2017a). Together with the relative fluctuation of kinetic energy, the flow regime can be categorized into intermittent flow, consistent flow and transition fluctuation flow (as seen in Fig. 24).  Fundamentals of particle flow modeling. With regard to the particle phase, our group has proposed a series of models, including the generalized hard particle models for arbitrary shape (GHPM) (Gui et al., , 2016d, the soft-sphere imbedded pseudo-hard-particle model (SIPHPM) , and coupled approach of soft-and hard models (EHPM-DEM) . These models are capable for simulation of a number of generalized gas-particle flows, including the pebble bed reactors.       Fundamentals of radiation modeling in the pebble beds. Our group has also proposed a series of numerical radiative models, including the short range model, long range model, microscopic model, the sub-cell model, semi-empirical model, and smoothed void fraction method, coupled with CFD-DEM simulation framework to study the effective thermal conductivity coefficient of the packed pebble bed (Wu et al., 2016(Wu et al., , 2017(Wu et al., , 2018a(Wu et al., , 2018b(Wu et al., , 2018c.

Brief summary
In this work, the fundamental aspects relevant to the pebble bed high temperature gas-cooled reactor were reviewed, focusing on the pebble bed flow. In our opinion, particular attention should be paid to the following points on pebble bed.


Fast simulation methods and models. The most numerical simulations of pebble beds are based on the computational fluid dynamics simulation and discrete element methods. These two methods all cost a huge amount of computational time. It is better to develop fast numerical models to fit for the industrial demand and requirements for fast design and safety analysis.  Accurate prediction. Both the pebble bed pressure drop and effective conduction coefficient need to be predicted more accurately. As these two parameters are all dependent on time and various operational conditions. It is hard to give better temporary predictions under practical conditions.