This study attempted to conduct a topological data analysis of groups of particles in packed-bed structures to provide quantitative evaluations of void shapes. This study examined the spatial correlation between the packed bed structure and the holdup sites through geometric data obtained from the coordinates of the various particles composing the packed bed to isolate characteristic structural data for liquid holdup sites in a packed bed. When the study defined “bottlenecks” as narrow areas of a scale below capillary length, it was discovered that, in packed beds consisting of particles of a single diameter, the contribution to holdup was related to the number density of bottlenecks. Regarding the dependence on the void fraction of the holdup sites, as trends were demonstrated that differed from the continuous change that accompanies changes in the modified Capillary number, the difference from the dimensionless correlation occurs. When particles of differing diameters intermingle, the bottleneck number density increases near particles with small diameter, and the percentage of blockages from droplets increases. As position and density differ depending on the particle packing method, it is suggested that holdup sites decrease in number when particles with small diameter are appropriately dispersed.

In the lower part of a blast furnace, the molten slag’s occupation of the coke bed causes the gas to drift and the pressure to increase, presenting a potential operational issue.^1,2,3,4) In chemical engineering, the liquid volume fraction (the “holdup,” or hs) is used as an index of the degree of accumulation of the molten slag in the packed bed. Generally, for holdup, methods are widely used where estimates are made after a dimensionless correlation formula reflecting the characteristics of the molten slag and the packed bed structure is obtained via experiments using room-temperature models.^{5,6,7,8,9,10,11,12,13,14,15)} However, it has been indicated that, if a dimensionless formula obtained via room-temperature model experiments is applied to high-temperature systems, it will yield results that show different trends, for reasons that are as yet unclear.^10,12,14) The surface tension of a molten slag differs greatly from that of a room-temperature liquid; hence, it is believed that its downward flow behavior within a packed bed would differ from that of the room-temperature model. Furthermore, in recent years, it has become apparent that not only the melt’s wettability with a solid but also its repellence from solid surfaces affects the detachment at the solid–liquid interface.^16,17,18) The authors used a numerical simulation to demonstrate analytically that the spatial distribution of the solid particles affects holdup.^19,20,21,22) Therefore, the contribution to holdup from the localized void information in a packed bed, including averaged values such as solid particle diameter and sphericity, must be ascertained.

In recent years, a method has been proposed to discuss the impact of spatial heterogeneity on macroscopic physical quantities by describing short-range and long-range correlations between large, complex, and disparate groups of data. This study focused on the geometric shapes formed between the data points, and data analysis on the correlation data was performed to evaluate quantitatively the differences using topology.^23,24,25) Hiraoka et al. used a mathematical approach called “persistent homology” to describe the geometric characteristics of the differences between the liquid state of SiO₂ and the internal structure of glass.²⁶⁾ This method is capable of describing and characterizing the shapes of objects when analyzing the structures and networks created by groups of particles. When applied to the internal structure of a packed bed, this type of geometric analysis is expected to be capable of isolating characteristic structural information at the holdup sites, including void fraction and average particle size.

This study converted information from discrete point clouds representing the packed structure into geometric data representing the void shapes, and then examined the validity of the analysis and the spatial correlation with the holdup droplets. By showing the geometric data shapes obtained from the coordinates of the various particles that form the packed bed, this study attempted to clarify the spatial correlation between the packed structure and the holdup sites.

In persistent homology, the inputted discrete point cloud data are replaced with each control volume at multiscale resolutions. Based on the overlap between the control volume spheres around the discrete points, geometric information is obtained for points, edges, and triangles; voids in the data are detected, and based on their size, shape, and generation and vanishing behavior, the shape of the input data is characterized. This study follows this approach, using the coordinates of the center of each particle in the packed bed as representative discrete points, then continuously changing the effective radius of the control volume, r_e. The discrete element method (DEM) was used in generating the digital data for the packed bed. Details are provided in the section below. If the control volumes overlapped between two points, then an edge was marked; if the control volumes overlapped between three points, then a triangle was marked. Figure 1 depicts the calculation process and an example calculation. The edges created between particles were marked as black lines; the triangles were depicted as gray planes. When the effective radii r_e of the spherical control volumes that used the coordinates of the packed-bed particle centers as discrete points were set and expanded, the edges between the discrete points began to take shape at the point r_e = d_p/2. Moreover, when r_e was expanded, triangular planes were formed in the voids between the discrete points. These triangular planes represent voids where three particles are in very close contact.

Fig. 1.

Schematic diagram of topological data analysis, and a calculation example for a packed bed composed of a single particle diameter.

As stated in the previous section, difficulties have been reported in generating estimates regarding non-wettened holdup using dimensionless correlation. Here, this refers to the holdup characteristics on water-repellent surfaces. In systems where the wettability of the coke and the molten slag is poor, the static holdup mechanism is described via the capillary length, the ratio between gravity and the surface tension, as expressed by κ⁻¹ = σ ρg . If the scale of the void spaces is less than that of the capillary length, then droplets larger than the capillary space would be trapped in these voids, as the impact of the surface tension would exceed that of gravity. In these circumstances, the effective channel is severed, as the droplets accumulate in the funnel-shaped channel—the bottleneck. In contrast, if the scale of the void spaces in the packed bed is larger than the capillary space, then gravity will dominate, and the droplets will fall. The only droplets that will be trapped at this time will be those directly on the surface of the coke, in the direction of gravity. The objective of this study is to quantify the differences between these characteristics. First, by taking the coordinates of the centers of the particles composing the packed bed to form a discrete point cloud, a three-dimensional Delaunay triangulation scheme is used to partition and reconstruct the space using tetrahedrons, and each surface is characterized. Figure 2 depicts a schematic of a droplet trapped by packed particles and the form of the void. The triangular planes forming the tetrahedrons correspond to the void structures formed between three touching particles. More precisely, the shapes are delineated by three fan shapes; there is at minimum a correlation with the formed triangular planes, and the vectors normal to the surfaces match. The size of the void, A_v, created by the space between three particles of known diameters and equal effective radii is calculated geometrically; thus, the equivalent circular area diameter of the void, which is given by the formula d_v = 4 A v π , can be calculated. Considering that the particle resolution is limited in numerical simulation, if d_v < 1.05 κ⁻¹, then the given triangular plane is defined as a “bottleneck.” Furthermore, if the distance between the triangular plane of the bottleneck and the holdup droplet is less than the spherical particle radius, the bottleneck is defined as a “holdup site.”

Fig. 2.

Schematic diagram of the dynamic droplet holdup mechanism in the bottleneck.

The analysis targeted a rectanguloid space measuring W 200 × D 200 × H 400 mm. A simulation was conducted via DEM where spheres of diameter 25.0 mm or 12.5 mm were pseudo-randomly distributed, dropped by gravity, and packed into the space. The coefficient of friction was set as a parameter varying from 0.01 to 5.0; the packed structure was also varied. A volume 200 × 200 × 10 mm of molten slag was positioned above the packed particle bed, and the smoothed-particle hydrodynamics (SPH) approach was used to create a liquid-flow simulation. Assuming a general 45SiO₂-45CaO-10Al₂O₃ composition (at 1873 K), the base conditions were set at ρ = 2600 [kg/m³], μ = 0.1 [Pa·s], and σ = 0.820 [N/m]. Furthermore, the coke was assumed to be a non-wettable surface that would not react with the molten slag. For the specific calculation scheme, please refer to the information provided in previous papers.^19,20,21,22)

This section will focus on packed beds composed of particles of a uniform diameter. First, to identify changes in the holdup characteristics caused by the packed structure, the dimensionless correlation obtained from the simulation and the analysis results were compared. Figure 3 depicts packed structures featuring single-size spherical particles generated via DEM with various coefficients of friction. Here, the diameter of the packed-bed particles was 25.0 mm, and the number of particles was 600. As the coefficient of friction increases, the packed-bed particle number density tended to decrease; in square control volumes established within the packed bed, the porosities were (a) ε = 0.3663, (b) 0.3675, (c) 0.4204, (d) 0.4253, and (e) 0.4580. Figure 4 depicts the relationship between the results of the liquid-phase drop simulation and the bottlenecks. The liquid phase does not pass through all the voids but instead travels through specific passages; the holdup droplets are distributed spatially in the dispersed phase. The blue bottlenecks (triangular planes) and the yellow holdup droplets decrease in number as porosity increases; both numbers are assumed to have a fixed correlation. Figure 5 depicts the relationship between the volume of trapped droplets, the number of droplets, and the void fraction of the number of triangles formed. As porosity (the void fraction) increased, the volume of trapped droplets, the number of droplets, and the number of bottlenecks all tended to decrease. In particular, the droplet volume greatly decreased in a parabolic fashion in opposition to the void fraction. Therefore, when considering the contribution of bottlenecks versus holdup, it is necessary to clarify not only simple numbers but also the spatial relationships. Figure 6 depicts a comparison with the dimensionless correlation regarding the holdups and the spatial correlation regarding the holdup droplets and the positions of the triangular planes. For the dimensionless correlation, the static holdup relationship hs = 9.96 (Cp_m^*)^−1.38 for the molten slag in the coke bed from Ohgusu et al. is used.¹⁰⁾ As the number of capillaries in the immersional wetting mode, Cp_m^* = ρg (ϕε D p ) 2 | σcos θ ¯ | (1-ε) 2 , was used, hs is not dependent on viscosity. Based on Fig. 6(a), as Cpm* increases (as ε increases), the correlation formula demonstrates a trend toward overestimating the simulation results. To explain this, Fig. 6(b) shows the relationship between the angle θ, created by the normal vector to the bottlenecks in the packed bed and the z axis, and the holdup droplets. Here, the bottleneck-holdup ratio, C, yields the bottleneck contribution ratio and is defined as follows:   

Fig. 3.

Simulation results of packed single-size particles by DEM: (a) f_r = 0.01 (ε = 0.3663), (b) f_r = 0.1 (ε = 0.3675), (c) f_r = 0.5 (ε = 0.4204), (d) f_r = 1.0 (ε = 0.4253), (e) f_r = 5.0 (ε = 0.4580).

Fig. 4.

Calculation results of (a) molten slag trickle flow simulation and bottleneck (blue triangle)-droplet holdup site (yellow) distributions in each void fraction: (b) ε = 0.3663, (c) ε = 0.3675, (d) ε = 0.4204, (e) ε = 0.4253, (f) ε = 0.4580. (Online version in color.)

Fig. 5.

Relationship between holdup droplets, the number of bottlenecks, and void fraction.

Fig. 6.

Spatial characteristics of static holdup droplets in each void fraction. (a) Comparison between non-dimensional correlation and simulation results; (b) bottleneck-holdup ratio as a function of the angle between the normal vector to the triangle and the z-axis. (Increment interval for θ is 10°; for example, bottlenecks in the range of 30° < θ ≤ 40° are counted as 40°.) (Online version in color.)

Fig. 7.

Correlation rate of holdup site as a function of void fraction.

This section will maintain a fixed packed-bed structure and consider changes in holdup generated by changing the surface tension σ. Using the packed-bed structure depicted in Fig. 3(c) and varying the σ of the molten slag from 0.50 to 1.50 yields the results shown in Fig. 8. Here, κ⁻¹, calculated from σ, is normalized via the minimum void area-equivalent diameter d_v, created by three particles. As κ⁻¹/d_v increases, the trapped droplet volume, the number of droplets, and the number of bottlenecks all increase. When compared with the curves in Fig. 5, the impact of σ is parabolic for droplet volume and the number of bottlenecks but logarithmic for the number of droplets. To explain this, Fig. 9 compares the dimensionless correlation with the effects of the simulation and depicts the spatial correlation between the holdup droplets and the positions of the triangular planes. Figure 9(a) illustrates that, as Cp_m^* increases (and κ⁻¹/d_v decreases), the correlation equation tends to overestimate the simulation results. At first glance, this appears to be a trend similar to that shown in Fig. 6(a); however, as shown in Fig. 9(b), the change in the holdup mechanism owing to the change in surface tension is somewhat different from cases when the void fraction is changed. In other words, in the relationship between angle θ, formed by the normal vector to the bottleneck’s triangular plane and the z-axis, and the holdup droplets, while no significant difference can be observed if κ⁻¹/d_v > 1 for the value C on the low-angle side, if κ⁻¹/d_v < 1, C suddenly approaches zero. This corresponds to the loss of the contribution of the bottleneck if κ⁻¹/d_v < 1. Figure 10 shows the relationship between κ⁻¹/d_v and λ. In the range of κ⁻¹/d_v > 1, as capillary length increases, it approaches λ = 1, and it is evident that almost all bottlenecks have droplets trapped within. When κ⁻¹/d_v < 1, λ suddenly decreases, and it is illustrated that holdup is negligible owing to the bottleneck mechanism. It has been numerically demonstrated that, in systems dominated by surface tension, droplet volume in specific holdup sites greatly increases and that the number of bottlenecks contributing to holdup increases.

Fig. 8.

Calculation results of bottleneck (blue triangle)-droplet holdup site (yellow) distributions in each surface tension condition. (a) σ = 0.50; (b) σ = 0.82; (c) σ = 1.00; (d) σ = 1.50. Note that the packed structures are similar to Fig. 3(c) (ε = 0.4204). Relationships between holdup droplets, the number of bottlenecks, and normalized capillary length shown. (Online version in color.)

Fig. 9.

Spatial characteristics of static holdup droplet in each κ⁻¹/d_v value. (a) Comparison between non-dimensional correlation and simulation results, and (b) bottleneck-holdup ratio as a function of the angle between the normal vector of triangle and z-axis. (Online version in color.)

Fig. 10.

Correlation rate of holdup site as a function of normalized capillary length.

In the blast furnace, a large distribution of particle diameters may be created through the use of nut coke and volume breakage of pellet coke. In the dimensionless correlation, when holdup and spatial information are considered, there are concerns that divergence from the actual phenomenon may occur, as the process uses information—such as particle diameter and void fraction—that ends up spatially averaged. This section uses geometric analysis to examine in detail the characteristics of holdup in systems with a mix of particle sizes.

To compare mixed-particle packed beds with single-size spherical-particle packed beds, particles having a large diameter of 25.0 mm (500 particles) and particles having a small diameter of 12.5 mm (800 particles; N_p2) were used for packing. Furthermore, for this section, the surface tension was adjusted so that bottlenecks involving only particles of small diameter were formed, with σ = 0.50 [N/m]. Figure 11 depicts volumetric data for packed beds with differing packing systems, as well as bottleneck and holdup particle distributions obtained through geometric analysis. In (a)’s large-diameter particle packed-bed structure (#1), bottlenecks do not form, whereas in (b) through (d) (#2, 3, and 4), bottlenecks would form near the interfaces of the particles with large and small diameters. The positions of the bottlenecks greatly differ depending on the particle packing system, and the holdup positions also adjust roughly in response to this phenomenon. In the case of #2, it is believed that the penetration of the particles of small diameter among the particles of large diameter via the action of gravity can be confirmed; it is believed that local structural changes with #3 can be quantified via geometric data analysis. In contrast, in the case of #4 and its random structure, a trend was demonstrated where bottlenecks increased at parts where particles of small and large diameters were in contact. These different packing structures affect the average void fraction (#2: 0.460, #3: 0.491, #4: 0.474). Figure 12 depicts a comparison of the dimensionless correlation with the simulation results, as well as the spatial correlation regarding holdup droplets and the positions of bottlenecks. Here, the average particle size was used to yield Cp_m*. Though the prediction accuracy was high in #4 (hs = 0.290), the changes in position in #2 (hs = 0.454) and 3 (hs = 0.498) indicated that hs nearly doubled. The packed bed structure in which the small particles are concentrated forms a number of void sizes smaller than κ⁻¹. At this time, hs increases. Generally, the void fraction is decreased by randomly mixed particles with different diameter ratio, but this time, the hold-up amount is reduced. This depends on the ratio of the number of particles to be mixed. From the distribution in Fig. 12(b), in #2 and 3, the value of C is large throughout the domain of θ, whereas in #4, as in the previous sections, on the low-angle side, C is drastically larger. The droplets trapped in the bottlenecks formed between the small grains occupy all the voids. Therefore, if there is a high spatial density of bottlenecks, then regardless of the angle θ, the bottlenecks and the droplets will come in contact. When the particles are distributed randomly, as in #4, then the lower the bottleneck density, the more dominant droplet blockages are through the bottleneck mechanisms shown on the low-angle side. Regarding the correlation of the holdup sites, while #2 and 3 exhibited high conformity, with λ = 0.90, #4 showed fairly low conformity, with λ = 0.65. With #2 and 3, where particles of small diameter are intensively packed, as with cases where surface tension predominates in the single-size grain cases in the previous section, droplets are trapped in almost all the bottlenecks. With #4, however, it is believed that a different trapping mechanism is at work.

Fig. 11.

Calculation results of packed beds from a two-particle-diameter system (f = 0.5). (a): Packed bed consisting of only large particles (#1), (b) upper small particles-lower large particles structure (#2); (c) upper large particles-lower small particles structure (#3); (d) random mixed structure (#4). (e)–(h) show the bottleneck (blue triangle)-droplet holdup sites (yellow) distributions in each packed structure, respectively. (Online version in color.)

Fig. 12.

Spatial characteristics of static holdup droplets in various packed-bed structures formed from a two-particle-diameter system. (a) Comparison between non-dimensional correlation and simulation results; (b) bottleneck-holdup ratio as a function of the angle between the normal vector of the triangles and the z-axis. (Online version in color.)

Through topological analysis to evaluate the structure of the voids from the geometric shapes created between the discrete points, the characteristics of the droplet holdup sites that accumulate in a packed bed with poor wettability composed of spherical particles were analyzed. The following conclusions were drawn.

Packed beds composed of single-size particles demonstrated a trend where the number of bottlenecks—narrow areas formed between three particles—decreased in a linear manner as the void fraction increased. In contrast, according to the numerical fluid simulation, holdup drastically decreased as the void fraction increased. This led to the new discovery that holdup contribution is related to bottleneck density. As the influence of holdup site void fraction demonstrated a trend that differed from the continuous change that accompanies variability in revised capillary numbers, divergence from the dimensionless correlation also proved to be a factor. A similar phenomenon is evident in cases where surface tension declines and the scale of the capillary space is smaller than the equivalent circular area diameter of the narrow voids; the holdup sites dramatically decrease. In contrast, when the scale of the capillary length increases, the probability of droplet obstruction in the bottleneck approaches 1.

In cases featuring particles of differing diameters, numerous bottlenecks were formed by particles of small diameters, and a trend was demonstrated where the blockage rate would increase when the bottleneck density was high. Positions and densities differ depending on the particle packing method, suggesting that holdup decreases in arrangements where particles of small diameters were appropriately dispersed.

A part of this research was supported by the Steel Foundation for Environmental Protection Technology (SEPT) of Japan, and ISIJ Research Promotion Grant.

This study focused on packed beds composed of spherical particles. By using the coordinates of the centers of spheres as representative discrete points, the relationship between void shape and holdup characteristics was clarified. In an actual blast furnace, packed-bed structures are composed of aspherical coke pellets with a range of grain sizes. The method used here presents the results of experimental analysis of the effect of aspherical elements on bottleneck formation behavior. To ensure that the packed particles are of approximate equivalent volume, data were generated for two types of packed beds: a bed composed of spherical particles of uniform size (particle diameter d_p = 25 mm; number of particles = 600; ε = 0.4204) and a bed composed of aspherical coke pellets (particle diameter = 15 to 35 mm; number of coke pellets = 600; ε = 0.4380). For the aspherical coke pellets, isovolumetric-equivalent diameters were calculated from the actual volume, and effective radii about the centers of gravity were derived from the product of the isovolumetric-equivalent diameters and the effective radius scaling factor. Figure A1 shows the changes in the number of triangular planes in each packed bed when the effective radius of the control volume r_e = sd_p⁄2 is changed. Here, the amount of increase in the number of triangular planes when the effective radius of the control volume is changed in 0.1-mm increments is displayed. For spherical particles, for an effective radius of r_e = 14.5 mm (s = 1.16), triangular planes suddenly formed in the voids in the packed bed, appearing as a peak. In contrast, the coke-pellet packed bed did not yield a sharp peak like the single-size spherical particles did; instead, the formation of a broad curve was evident. This is due to the presence of complex void structures from the form characteristics of the coke pellets such as particle size and sphericity, which did not produce simple void structures as did the single-size spherical-particle packed bed. In contrast, the peak cresting in the neighborhood of s =1.4 and the distribution in the domain of s > 1.4 are approximately the same for both the spherical particles and the coke pellets. Therefore, packed beds composed of aspherical particles such as coke pellets have points of commonality with packed beds consisting of single-size spherical particles with an averaged matching particle diameter in terms of void structure and long-distance structure. Differences in local void structures were identified in the s < 1.4 domain. For coke-pellet packed beds, the presence of narrow parts even in the very short range in the neighborhood of s = 1.0, which cannot be composed of groups of spherical particles, is vital. This area demonstrates extremely tight packing, with characteristic bottlenecks.

Fig. A1.

Characteristics of calculated packed beds consisting of 15 types of coke.

• 1)   M.  Baniasadi and  B.  Peters: ISIJ Int., 58 (2018), 25.
• 2)   S.  Ghosh,  N. N.  Viswanathan and  N. B.  Ballal: Steel Res. Int., 88 (2017), 1600440.
• 3)   G. S.  Gupta,  J. D.  Litster,  V. R.  Rudolph,  E. T.  White and  A.  Domanti: ISIJ Int., 36 (1996), 32.
• 4)   P. R.  Austin,  H.  Nogami and  J.  Yagi: ISIJ Int., 37 (1997), 458.
• 5)   T.  Fukutake and  K.  Okabe: Tetsu-to-Hagané, 66 (1980), 1947.
• 6)   T.  Sugiyama,  T.  Nakagawa,  H.  Sibaike and  Y.  Oda: Tetsu-to-Hagané, 73 (1987), 2044.
• 7)   H.  Kawabata,  Z.  Liu,  F.  Fujita and  T.  Usui: ISIJ Int., 45 (2005), 1466.
• 8)   H.  Kawabata,  K.  Shinmyou,  T.  Harada and  T.  Usui: ISIJ Int., 45 (2005), 1474.
• 9)   S. J.  Chew,  P.  Zulli and  A.  Yu: ISIJ Int., 41 (2001), 1112.
• 10)   H.  Ohgusu,  Y.  Sassa,  Y.  Tomita,  K.  Tanaka and  M.  Hasegawa: Tetsu-to-Hagané, 78 (1992), 1164.
• 11)   K.  Sunahara,  K.  Nakano,  M.  Hoshi,  T.  Inada,  S.  Komatsu and  T.  Yamamoto: Tetsu-to-Hagané, 92 (2006), 875.
• 12)   D.  Jang,  M.  Shin,  J. S.  Oh,  H. S.  Kim,  S. H.  Yi and  J.  Lee: ISIJ Int., 54 (2014), 1251.
• 13)   T.  Sugiyama,  T.  Nakagawa,  H.  Sibaike and  Y.  Oda: Tetsu-to-Hagané, 73 (1987), 2044.
• 14)   W. M.  Husslage,  T.  Bakker,  A. G.  Steeghs,  M. A.  Reuter and  R. H.  Heerema: Metall. Mater. Trans. B, 36 (2005), 765.
• 15)   H. L.  George,  R. J.  Longbottom,  S. J.  Chew and  B. J.  Monaghan: ISIJ Int., 54 (2014), 820.
• 16)   S.  Ueda,  T.  Kon,  S.  Kikuchi,  S.  Natsui,  S.  Sukenaga and  H.  Nogami: ISIJ Int., 55 (2015), 1277.
• 17)   T.  Kon,  S.  Sukenaga and  S.  Ueda: ISIJ Int., 57 (2017), 1166.
• 18)   D.  Nakajima,  T.  Kikuchi,  S.  Natsui and  R. O.  Suzuki: RSC Adv., 8 (2018), 37315.
• 19)   S.  Natsui,  A.  Sawada,  K.  Terui,  Y.  Kashihara,  T.  Kikuchi and  R. O.  Suzuki: Tetsu-to-Hagané, 104 (2018), 347.
• 20)   S.  Natsui,  A.  Sawada,  T.  Kikuchi and  R. O.  Suzuki: ISIJ Int., 58 (2018), 1742.
• 21)   S.  Natsui,  K.  Ohno,  S.  Sukenaga,  T.  Kikuchi and  R. O.  Suzuki: ISIJ Int., 58 (2018), 282.
• 22)   S.  Natsui,  A.  Sawada,  K.  Terui,  Y.  Kashihara,  T.  Kikuchi and  R. O.  Suzuki: Chem. Eng. Sci., 175 (2018), 25.
• 23)   M.  Kimura,  I.  Obayashi,  Y.  Takeichi,  R.  Murao and  Y.  Hiraoka: Sci. Rep., 8 (2018), 3553.
• 24)   M.  Gameiro,  Y.  Hiraoka,  S.  Izumi,  M.  Kramar,  K.  Mischaikow and  V.  Nanda: Jpn. J. Ind. Appl. Math., 32 (2015), 1.
• 25)   J. M.  Chan,  G.  Carlsson and  R.  Rabadan: Proc. Natl. Acad. Sci. U.S.A., 110 (2013), 18566.
• 26)   Y.  Hiraoka,  T.  Nakamura,  A.  Hirata,  E. G.  Escolar,  K.  Matsue and  Y.  Nishiura: Proc. Natl. Acad. Sci. U.S.A., 113 (2016), 7035.