Topological analysis of pattern formation in cooling granular gases confined by elastic wall

In this paper, we investigate the topological characteristics of the pattern formation in the cooling granular gases confined by the elastic wall. The persistent homology and Voronoi’s analysis and its derivative analyses are applied to accomplish our aim. The growth of the pattern formation can be identified by the switch between the logarithmic concave and logarithmic convex in the life-span-distribution obtained using the persistence diagram. Furthermore, three phases are identified by the zeroth or first order Betti number, when a form of the wall is the square. Finally, the characteristics of the coordination of granular particles condensing around the elastic wall are investigated by the Voronoi’s analysis, bond-angle analysis, and polyhedral template matching. We confirm that some clusters of the granular particles condensing around the elastic spherical-wall certainly attribute to their crystallization categorized as the typical coordination.


Introduction
The accurate classification of the phase of matter [1] is significant to characterize the physical property of matter such as solid (crystal structure), liquid, or gases. Meanwhile, such a classification of the phase of matter sometimes involves difficulties as observed in the classification of the phase of glass [2]. For example, we are still unable to answer to the question whether glass corresponds to the liquid, solid or glassy state. The topological characterization of the phase of matter on the basis of locations of atoms or particles has been used in order to overcome such difficulties involved with the specification of the phase of matter. For example, the topological characterization of the crystal has been widely used, as represented by the Voronoi's analysis (VA) using coordination of atoms [3]. On the other hand, the specification of the coordination of atoms (particles) from their locations becomes difficult, when thermal fluctuations of atoms (particles) become significant. Therefore, the specification of the coordination of thermal fluctuating particles involves difficulties even with the VA. Meanwhile, the recent application of the persistent homology (PH) [4] to numerical data (i.e., locations of atoms) of glass by Hiraoka et al [5] succeeded the discrimination of the glassy state from the liquid state, where the coordination of atoms in glass seems to be random, because thermal fluctuations are significant at glance. Then, the PH enables us to find the ordered structure even under the marked thermal fluctuations. The analogy between glass and granular matter is sometimes indicated by the fact that sheared granular media forms the disordered glassy state via the jamming transition [6]. In previous studies of the dense granular packing [7], the characteristics of the force-chain-network was analyzed by the PH, whereas the crystal characteristics of the densely packed granular particles was investigated by the VA.
As a peculiar characteristics of the dilute granular gases, the cluster formation from the initially homogeneous cooling state (HCS) [8,9] has been well studied by Brey [10] and his coworkers, whereas the characteristics of the HCS of the granular gases has been studied in detail by Santos [11], Brilliantov [8] and their coworkers or Yano [12]. The mode-analysis [8,10] of the hydrodynamics equation of the granular gases can demonstrate such an instability of the HCS, which attributes to the pattern formation via aggregations (clusters) of the granular particles. Of course, the granular gases correspond to the status, in which the volume fraction of Betti number on the basis of locations of granular particles is calculated using the monotone PNG format with the CHomP [21]. Figure 2 shows the zeroth and first order Betti numbers, which are calculated using monotone PNG. We can readily understand that both B 0 (connection) and B 1 (hole) are equal to '4', as shown in figure 2.

Persistence diagram (PD)
For a nonnegative integer p, the p-th Betti number counts 'p-dimensional holes' of a simplicial complex, as we explain in section 2.1. A persistence diagram (PD) is a diagram, which records birth and death of 'p-dimensional holes', i.e. we consider 'time-evolution' of a simplicial complex and are interested in the change of its topology.
The time-evolution, formally speaking, corresponds to filtration of a simplicial complex. A filtration of a simplicial complex X is given by a sequence A filtration can be understood as a movie film which records how a simplicial complex grows as time evolves. We give an example in figure 3; the filtration Î X t t { } is an empty complex for sufficiently small t, say t=t −1 ; the filtration Î X t t { } is constant except for t=t 0 , t 1 , t 2 , where it grows up, as shown in figure 3. The p-th PD of a filtration Î X t t { } is induced by the data of when the p-dimensional holes are born and dead as time evolves. Informally speaking, for a p-dimensional hole h of some X t , we mark a point b d , , resp.) is the time when the hole h is born (dead, resp.). The p-th PD is obtained by plotting such b d , h h ( ). The PD is given by multiple set in general since there is a possibility that some holes have the same birth time and death time. We give an example of the PD in figure 4, which is obtained from the filtration in figure 3. Note that the pointt t , 1 0 ( ) in the PD has multiplicity of '2'. There is an obvious problem with respect to the definition of such a hole h in the time-evolution. The notion of persistent homology and some decomposition theorems are necessary to define PDs in a formal way. Their brief overview is given in appendix C.

From point cloud data to PD
In this subsection, an overview about how the PD is obtained from point cloud data is given. Since the PD is obtained from a filtration of a simplicial complex, it suffices to construct a filtration of a simplicial complex There are several ways to construct a cell complex from the open covering R U : the Čech complex, the Alpha complex, the Vietoris-Rips complex, etc [20,22,23].
If we denote one of such complexes by C R , then, gives a filtration of a cell complex where we consider R as time, i.e.
We note that the PD is determined by the point cloud data S in  N . Figure 5 shows growths of R U of eight point cloud by enlarging R. The upper-left frame of figure 5 shows the birth of the hole, when = 1 R R. The upper-middle frame of figure 5 shows births of eight holes, when = 2 R R . The upper-right frame of figure 5 shows deaths of holes, when = 3 R R . Consequently, the PD of these birth-death-sets of holes , b d R R ( ) are shown in lower-frame of figure 5, in which b R and d R correspond to the radius of B R (i.e., R), which yields the birth and death, respectively. Of course, eight-fold-points are plotted on ) ( )owing to the birth-death of eight holes.

VA, BAA and PTM
The Voronoi's analysis (VA), bond-angle-analysis (BAA) and polyhedral template matching (PTM) have been used to analyze the structure of the crystal on the basis of coordination of atoms. In particular, the analytical tools are freely provided by the codes of Voro++ [24] and Vorotop [25] in Ovito [26]. Then, we can utilize such analyzers in order to analyze the crystal structures of granular particles. Of course, the crystal structure usually postulates the densely packed status of granular particles, where thermal fluctuations of atoms are negligible. We, however, focus on the crystal characteristics, before and after granular particles are densely condensed around the wall owing to their markedly low kinetic energy. Here, schematics of the VA, BAA and PTM are demonstrated, briefly. In order to discuss the VA, BAA and PTM, we start our discussion by mentioning to the Voronoi's diagram (VD).
 is space such as Í   n ( Î  n ). Let's think elements ≔ . We define the distance between the point Î  P and g i as 9), when = 1 R R (upper-left frame), = 2 R R (upper-middle frame) and = 3 R R (upper-right frame) in Čech filtration. PD obtained by birth-death of holes corresponding to above three frames (lower frame).
. We call these divided structures of  by  g R , i ( ) as the VD.
The VD is described by the Delaunay tetrahedralization in the case of Í   3 . The significant problem in the VA is classifications of the VD using the typical crystal-coordination such as the face-centered-cubic (FCC), body-centered-cubic (BCC), hexagonal close-packed (HCP), icosahedral (ICO) and hybrid of FCC and HCP (FCC-HCP). The VA by Voro++ applies Weinberg's algorithm [27] to find the most appropriate crystalcoordination among the FCC, BCC, ICO, FCC-HCP and Other from the VD, in which Other indicates the status of coordination of granular particles, from which any crystal structure is not specified with the FCC, BCC, ICO and FCC-HCP.
Of course, there are other methods in order to categorize the crystal structure such as the Common Neighbor Analysis (CNA) [28], PTM [15], BAA [14] other than the VA. The PTM was proposed by Larsen-Schmidt-Schiøtz [15]. In the PTM, the similarity between two diagrams is evaluated by the Root-Mean-Square Deviation (RMSD), which is defined by [15] where Q is the right handed orthogonal matrix , = À å = À å  are barycenter of v and w. s is the optimal scaling of v. The methods of findings of s and Q are proposed by Horn [29] and Theobald [30]. ℵis the number of the neighboring atoms, and À = 6 for the simple cubic (SC), 12 for the FCC, 12 for the HCP, 12 for the icosahedral (ICO) and 14 for the BCC are used. The detail of algorithm of the PTM is demonstrated by Larsen-Schmidt-Schiøtz [15]. Finally, we find the smallest RMSD among them which are calculated using the SC, FCC, HCP, ICO, and BCC. Additionally, Other corresponds to the state that any structure is not identified even with the SC, FCC, HCP, ICO, and BCC, when the minimum RMSD is larger than its critical value. Actually, v is the vector which indicates the vertex-set of the convex hull formed by N-neighboring atoms and w is the vector, which indicates the vertex-set of the convex hull of the reference templates, namely, the SC, FCC, HCP, ICO, and BCC (see figure 6 for convex-fulls of the SC, BCC, FCC, HCP and ICO). Hence, we can calculate the distribution of the RMSD ( , : time), when we specify the crystal structure of the atoms using two categories Other and XX (XX ≔ SC, FCC, HCP, ICO and BCC). Finally, we mention to the BAA [14], briefly. The BAA was proposed by Ackland and Jones [14]. The BAA identifies the crystal structure of atom-i with the FCC, HCP, ICO, and BCC from angles (θ jik ) between two bonds (r ij and r ik ), which connect two sets of two neighboring atoms (i-j and i-k). The neighboring atoms are searched by the mean square distance of six neighboring atoms around the atom-i. Once the calculation of the angle between two bonds for several sets of two bonds is finished, functional of χ, which are calculated by angles θ jik , identifies the crystal structure of atom-i with the FCC, HCP, ICO, BCC or Other.

Event-Driven method to simulate granular particles
The event-driven (ED) method has been frequently used in order to simulate a large number (N) of granular particles, whose simulation is difficult using the discrete element method (DEM), when the calculation of the -N N 1 2 ( ) paired force-chain is beyond the computational resource. In particular, the ED method is suitable to simulate the dilute granular gases, in which the effects of deformations of contacting granular particles on the collective motion of all the granular particles are negligible.
The algorithm of the ED method is so simple. Provided that the k-th collision occurs at t=t k ( Î +  t k ), we search for the occurrence-time of k+1-th collision, namely, t k+1 , which is calculated by Figure 6. Schematic of convex hulls of SC, BCC, FCC, HCP and ICO used in PTM [15].
is the location of the granular particle indexed by i at Î  v t, i d is the velocity of the i-th granular particle, r i is the radius of the i-th granular particle.
Provided the boundary (wall) is considered, the occurrence-time of the collision between the granular particle and wall is calculated in a similar way to equation (3). Afterward, we compare the next collisional time calculated by two granular particles with that calculated by the granular particle and wall and select smaller t k+1 as next collisional time.
Once t k+1 is determined, we revise the locations of all the granular particles using ). Finally, we revise the velocities of colliding granular particles or velocity of the granular particle colliding with the wall. Provided that the i-th and j-th granular particles collide, the velocities of colliding granular particles change at t=t k+1 as follows: ] is the restitution coefficient, g ij is the relative velocity and Î W n ij 2 is the relative location-unitvector.
The calculation of the velocity of the granular particle post-collision with the elastic wall is calculated in a similar way to equation (4) [17].

Persistent homology and their characteristics
In this section, we apply the PH to the cooling process of two dimensional granular gases confined by the elastic wall. Firstly, we investigate the PH of the granular discs in the case of the square (SQ) boundary. Next, we investigate the PH in the case of the circular (CI) boundary to consider the effects of form of the elastic wall on the PH. Finally, we investigate the PH of the granular spheres in the case of the spherical (SP) boundary.

Results in case of SQ boundary
The restitution coefficient of the granular discs is set as ò=0.85, and the volume fraction of the granular discs is fixed as f=7.85×10 −2 . The length of one side of the square box (SQ boundary) is set as L=300. Three types of the diameter of granular discs (r d ), namely, r d =0.15, 0.21 and 0.47 are considered. As a result, total number of granular discs (N) is set as N=10 4 in the case of r d =0.47, N=5×10 4 in the case of d=0.21 and N=10 5 in the case of r d =0.15 owing to the constant volume fraction (f=7.85×10 −2 ). The initial velocities of the granular discs are randomly distributed in the range of Î - ] and the initial positions of the granular discs are randomly distributed inside the elastic wall, namely, the range of Î - The time-evolution of granular discs is calculated using the event-driven (ED) method [17], because the conventional molecular dynamics requires the vast calculation time, when forces between -N N 2 2 ( ) paired granular discs are calculated. Figures 7-9 show time-evolutions of the granular discs inside the SQ boundary at N=10 4 , 5×10 4 and 10 5 , respectively. Figure 7 indicates that some clusters are formed in the vicinity of the wall and the cluster with the maximum size rotates to clockwise direction along the wall. Such a rotation of granular discs along the wall is similar to the emergence of the ordered collective-motion in biological swarm inside the wall [31]. On the other hand, the granular particles tend to move away from the wall, when the heating wall with the constant temperature is used, as reported by Esipov and Pöschel [13]. The numerical results obtained using N=5×10 4 and 10 5 do not indicate such a rotation of the cluster with the maximum size, as shown in figures 8 and 9. Indeed, it is commonly observed in numerical results of N=10 4 , 5×10 4 and 10 5 that the several small-clusters aggregate toward the larger clusters via their connections during the time-evolution. Of course, the elastic wall is a mathematical toy model to avoid the freeze of the calculation via the excessive accumulation of granular discs with the small kinetic energy in the vicinity of the wall and long range correlation of the granular discs owing to the use of the periodic boundary condition, which is unfavorable to topological analyses, because the distance between two granular discs in the PH must be modified by considering the periodicity at the boundary. The calculations of the PDs are performed using locations of the center of granular discs, namely, x y , ). Therefore, we must remind that the connection in d=0-persistent homology (d: dimension of PH) is not equivalent to the physical connections due to contacting granular discs. Reminding that the domain, which is occupied by neighborhood of i-th point-cloud, whose center is set as x y , is called as the life-span of the hole, and defined by: The image of the life-span (ℓ) is shown in the lower-frame of figure 5.
The open source program homcloud [32] is used to calculate the PD in our study. The PD is plotted using not , R R ( )owing to the specification of the homcloud. Figure 10 shows PDs (left-half frames) and birth of holes with max ℓ (maximum value of ℓ) via connections of some of at N c =0, 2.15M, 3.85M and 11.4M (right-half frames), when N=10 5 (N c : collision number, M: million), respectively. The color of the contour expresses the number of holes, which has , b d R R ( )and length of one side of the square is set as D Î +  . Readers remind that The expeditious method to read the characteristics of the PD is to calculate the life-span-distribution (lsd) by each time step. Then, we approximate the lsd, namely, - where , A B, and Í +  n , , , D C E and m Í. The reason why we use equation (6) is that it can express the hybrid of the logarithmic convex (log-convex: and x<y)) and logarithmic concave (log-concave: and x<y)) using appropriate set of n m , ( ) in equation (6). Similarly, equation (6) is also able to express completely concave or convex for all the range of ℓ (i.e., Î +  ℓ ). Figures 11-16 show f ℓ ( ) versus ℓ together with f ap ℓ ( ) for d=0 (connection) and 1 (hole) in cases of N=10 4 , 5×10 4 and 10 5 , whose n , , , , , , A B C D E F and m in equation (6) are defined in table 1. The authors will consider that one standard for the evaluation of the topological phase-transition of the granular discs is determined by the switch between the log-convex and log-concave of f ℓ ( ) (or f ap ℓ ( )), as discussed later. Figure 11 indicates that f ℓ ( ) (or f ap ℓ ( )) follows the log-concave at N c =0 and 0.25M, whereas f ℓ ( ) (or f ap ℓ ( )) follows the log-convex at 0.5M N c . Consequently, we conclude that the topological phasetransition occurs in the range of 0.25M< N c < 0.5M, when d=0 and N=10 4 . Similarly, figure 12 indicates that f ℓ ( ) (or f ap ℓ ( )) follows the log-convex at  9 ℓ and log-concave at < 9 ℓ in the case of N c =0, whereas  The lsd obtained using d=1-PD, which follows the log-convex, indicates that the number of holes with the large diameter d R and small b R (i.e., long life-span) is larger than that obtained using the d=1-PD, which follows the log-concave. The holes with the long life-span are generated by the increase of vacant space and  connections of point-clouds by small b R , which are caused by the growths of clusters. Consequently, the switch between the log-concave and log-convex becomes the standard to measure the drastic growth of clusters, by which the topological characteristics of granular gases changes, markedly.
From above discussions, we confirmed that the topological phase-transition can be judged by the switch of the life-span-distribution (lsd) (i.e., f ℓ ( )) between the log-convex and log-concave. Here, we investigate the topological phase-transition using the Betti number, B 0 (the zeroth order) and B 1 (the first order). As discussed in section 2, the Betti number is calculated using the PNG versions of figures 7-9. Therefore, the accuracies of the calculations of B 0 and B 1 depend on the number of pixels in PNG versions of figures 7-9. Figure 17 shows B 0 /N and B 1 /N in cases of N=10 4 (upper-left frame), 5×10 4 (upper-right frame) and N=10 5 (lower-left frame). B 0 /N ; 0.96 at N c =0 obtained using N=10 5 indicates that the numerical error exists owing to the insufficient number of pixels, whereas B 0 /N=1 at N c =0 obtained using N=10 4 and 5×10 4 indicate that the number of pixels is adequate to calculate B 0 and B 1 with the good accuracy.
Firstly, we confirm that there exist typical three phases (i.e., phases-I, II and III) in cases of N=5×10 4 and 10 5  We are, however, unable to identify phases-I and II, clearly, when N=10 4 , whereas the phase-III is conformed in N c ä [0.25M, ¥), when N=10 4 . The increase and decrease in B 0 and B 1 in the phase-III are caused by temporal changes in form of rotating cluster along the wall obtained using N=10 4 , so that such changes in phase-III do not attribute to intrinsic changes in the topological phase, which characterize the pattern formation from the HCS. It is the significant result that the topological phase-transition from the phase-I to the phase-II is identified by not the lsd but B 0 and B 1 , whereas the phase-transition due to the topological phase-transition from the phase-II to the phase-III is identified by both the lsd and B 0 and B 1 .

Results in case of CI boundary
Now, we investigate the PH in the case of the CI boundary in order to investigate the effects of form of the elastic wall on the PH. The radius of the granular disc (r d ) is set as r d =0. 15 and radius of the CI (R d ) is set as R d =150.
The set of the restitution coefficient and initial velocities of the granular discs are same as those used in calculations in the case of the SQ boundary. The total number of the granular discs (N) is set as = N 10 5 . As a result, the volume fraction of granular discs (f) is calculated as f=0.1, because all the granular discs are randomly distributed inside the CI boundary at t=0. Figure 18 shows snap-shots of the granular discs confined by the CI boundary at N c =0 (upper-left frame), N c =2.5M (upper-center frame), N c =5M (upper-right frame), N c =7.5M (lower-left frame), N c =15M (lower-center frame) and N c =35M (lower-right frame). We can confirm the clear pattern formation at 5MN c . In particular, the clusters seem to grow from the wall at N c =35M. Figure 19 shows the d=1-PDs at N c =0 (top-left frame), N c =2.5M (top-right frame), N c =5M (middleleft frame), N c =7.5M (middle-right frame), N c =15M (bottom-left frame) and N c =35M (bottom-right frame), when N c =10 5 and CI boundary. The tendency of changes in the d=1-PD during the time-evolution is similar to that obtained using the SQ boundary in figure 10. Figures 20 and 21 show f ℓ ( ) versus ℓ together with f ap ℓ ( ) (in equation (6)) versus ℓ when d=0 and d=1, respectively, where n , , , , , , A B C D E F and m in equation (6) are defined in table 2. We can confirm that the topological phase-transition occurs in N c ä [2.5M,5M] in both cases of d=0 and 1, as shown in figures 20 and 21. Figure 22 shows B 0 /N, B 1 /N and . Similarly to B 0 /N obtained using N=10 5 and the SQ boundary, B 0 /N ; 0.96 indicates that the numerical error exists owing to the insufficient number of pixels. We can readily confirm that phases-II and III exist in N c ä [0,2M] and N c ä [2M, ¥). It is not obvious in the present study whether the lack of the phase-I in the case of the CI boundary is caused by the different form of boundary from the SQ boundary or larger volume fraction (f=0.1) in CI boundary than that (f=0.078 5) in the SQ boundary. The phase-transition due to the topological phase-transition from the phase-II to the phase-III is obtained using both the lsd and Betti number in the case of the CI boundary.

Results in case of SP boundary
Next, we investigate the PH in the case of the SP boundary. As discussed in section 2, the 3D physical space postulates the PH with d=2, which corresponds to the spherical surface (void). Here, we consider two types of N. In one case, N=5×10 4 granular spheres with r d =1.5 are calculated. In the other case, N=10 4 granular spheres with r d =1.5 are calculated. In both cases, the radius of the SP boundary is set as R d =150 and restitution coefficient is set as 0.85. Consequently, the volume fraction is 5% in the case of N=5×10 4 and 1% in the case of N=10 4 , because all the granular spheres are randomly (almost homogeneously) distributed inside the SP boundary. The restitution coefficient of the granular spheres is set as 0.85 and the initial velocities of granular spheres are randomly selected in the ranges of Î - Figure 23 shows the snapshots of the granular spheres confined by the SP boundary, when N=10 4 (upper two raw) and N=5×10 4 (lower two raw). The granular spheres cluster in the vicinity of the SP boundary with the band-like form as N c increases, when N=10 4 ,. On the other hand, some band-like distributions of granular spheres emerge along the SP boundary in accordance with the increase in N c when N=5×10 4 . The aggregation of granular spheres in the vicinity of the wall at large N c is similar to those in 2D cases of the SQ and CI boundaries. Figures 24-27 show time-evolutions of d=1-PDs and d=2-PDs obtained using N c =10 4 and 5×10 4 , respectively. The tendency of temporal changes in d=1-PDs is similar to those obtained using the SQ and CI boundaries. d=2-PDs obtained using N=10 4 indicate that the life-spans of most of granular spheres approach to zero at N c =5.5M, when N=10 4 , although the life-spans of some of granular spheres are finite, as observed in d=1-PD at N c =5.5M in figure 24. Meanwhile, the life-spans of granular spheres observed in d=2-PDs are finite at N c =40M, when N=5×10 4 , as shown in figure 27.         topological phase-transition in accordance with the switch of the lsd between the log-concave and log-convex in both cases of N=10 4 and 5×10 5 . The Betti number is not calculated in the case of the SP boundary, because the monotone (0,1) pixels, which are used in the calculation of B 0 and B 1 in the SQ and CI boundaries, are unable to consider the depth-direction in 3D.

3D crystal classification of granular gases by VA BAA and PTM
Based on the numerical results obtained using the SP boundary in section 3, we investigate the 3D crystal classification of the granular gases confined by the SP boundary using the Voronoi's analysis (VA), bond-angle analysis (BAA) and polyhedral template matching (PTM). The numerical schemes of the VA, BAA and PTM have been already demonstrated in section 2, briefly. As discussed in Introduction, there are two goals in our study of the granular gases with the VA, BAA and PTM. One is to confirm whether the VA, BAA or PTM are able to identify the transition of the granular gases from the HCS to the clustering state or not. The other is to investigate which typical coordination (i.e., FCC, BCC, HCP, ICO etc.,) the crystal structures of the granular particles, which are condensed around the wall, densely, are categorized using the VA, BAA and PTM. We proceed our discussions of analytical results in order of the BAA, VA and PTM. Figure 34 shows snapshots of the categorized coordination per a granular particle and its fraction obtained using the BAA at N c =0, 4M, 5M, 6M, 10M, 20M, 30M, 40M and 47.4M. We can confirm that the spatially random locations of granular spheres at N c =0 (t=0) does not obtain any crystal structure. In short, all the coordination of granular particles is categorized as Other. Meanwhile, the fraction of granular spheres with the coordination categorized as the HCP, FCC and BCC (coordination) increases, as N c increases. Actually, 20% of granular spheres are categorized as the structured coordination, namely, FCC, HCP, BCC and ICO at N c =47.7M. Of course, most of granular spheres condensate around the wall at N c =47.7M, so that such a status is not appropriate to be called as the granular gases. Figure 35 shows the fractions of the categories of the coordination versus N c obtained using the BAA (Other: left frame, FCC, HCP, BCC and ICO: right frame). The enlarged figure (log-log plot of the fraction of Other . Then, we can discriminate the topological phase of the granular spheres between N c <4M (HCS) and 4MN c (growth of clusters). Reminding that the lsd for d=0 changes from the log-concave to the log-convex at 3M< N c <4M, as shown in figure 29, the BAA also specifies the topological phase-transition of the granular spheres with a similar accuracy to the lsd obtained using d=0-PD.  Figure 36 shows the snapshots of the categorized coordination per a granular sphere and its fraction at N c =0, 2M 4M, 5M, 6M, 10M, 20M, 30M, 40M and 47.4M obtained using the VA. We can confirm that 3.1% of randomly located granular spheres at N c =0 can be categorized as HCP or FCC (coordination). Finally, 32.8% of granular spheres obtain the crystal structures, namely, FCC, HCP, BCC, FCC-HCP and ICO at N c =47.7M. Figure 37 shows the fractions of the categories of the coordination versus N c obtained using the VA (Other: left frame, FCC, BCC, ICO, FCC-HCP and HCP: right frame). The enlarged figure (log-log plot of the fraction of Other versus N c ) is added to the left frame of figure 34 in order to show that the fraction of Other decreases by the inverse power law function of N c . As shown in figure 34, the fraction of the category Other is almost invariant during N c < 4M and decreases in accordance with

Results of VA
) . Then, the VA-results indicate that the topological phase-transition of the granular spheres occurs between N c <4M (HCS) and 4MN c (growth of clusters) similarly to the PD and BAA. The growths of the fractions of the FCC-HCP and HCP are markedly larger than those of the BCC, ICO and FCC. Then, we can conclude that the HCP type coordination is dominant in both cases of the BAA and VA, when the structured coordination are assigned to granular spheres.

Results of PTM
We consider the categorization of the coordination obtained using the PTM. Figures 38-46 show granular particles inside the SP boundary, which are colored in accordance with the assigned category of the coordination   Figure 38 shows that 48.7% of N is categorized as Other, whereas 45.3% of N is categorized as the SC and 4.9% of N is categorized as the HCP, at N c =0. Then, 50.2% of N is categorized as the structured coordination at N c =0, even when the BAA and VA categorize most of N as Other.

Concluding remarks
We investigated the transition from the homogeneous cooling state (HCS) to the pattern formation of the granular gases confined by the elastic wall using the topological analyses, namely, PH, BAA, VA and PTM. The topological phase-transition is successfully identified by the switch of the life-span-distribution (lsd) between the log-concave and log-convex, when the dimension of homology is zero and unity in the cases of the SQ and CI boundaries or zero in the case of the SP boundary. The zeroth and first order Betti numbers obtained using the SQ boundary indicate that three topological phases exist, whereas the zeroth and first order Betti numbers obtained using the CI boundary indicate that only two topological phases exist. The analyses of coordination of granular spheres on the basis of the BAA, VA and PTM indicate that the time-evolution of fractions of categories also become standard to judge the topological phase-transition as well as the switch of the lsd between the logconcave and log-convex. The HCP coordination is dominant, when the condensed granular spheres around the wall are categorized as the structured coordination after the adequate time passes by. Finally, we confirmed that these topological analyses, which have been applied to crystal or glass, are also useful for the consideration of the phase-transition of the granular gases. { } be an interval module component on an interval I as above. If < a b are endpoints of the interval I, then we say that the p-th dimensional hole is born at t=a and dead at t=b. By using an interval module decomposition, we obtain the PD.