Geometrical structure and thermal conductivity of dust aggregates formed via ballistic cluster-cluster aggregation

We herein report a theoretical study of the geometrical structure of porous dust aggregates formed via ballistic cluster-cluster aggregation (BCCA). We calculated the gyration radius $R_{\rm gyr}$ and the graph-based geodesic radius $R_{\rm geo}$ as a function of the number of constituent particles $N$. We found that $R_{\rm gyr} / r_{0} \sim N^{0.531 \pm 0.011}$ and $R_{\rm geo} / r_{0} \sim N^{0.710 \pm 0.013}$, where $r_{0}$ is the radius of constituent particles. Furthermore, we defined two constants that characterize the geometrical structure of fractal aggregates: $D_{\rm f}$ and $\alpha$. The definition of $D_{\rm f}$ and $\alpha$ are $N \sim {( R_{\rm gyr} / r_{0} )}^{D_{\rm f}}$ and ${R_{\rm geo}} / {r_{0}} \sim {\left( {R_{\rm gyr}} / {r_{0}} \right)}^{\alpha}$, respectively. Our study revealed that $D_{\rm f} \simeq 1.88$ and $\alpha \simeq 1.34$ for the clusters of the BCCA. In addition, we also studied the filling factor dependence of thermal conductivity of statically compressed fractal aggregates. From this study, we reveal that the thermal conductivity of statically compressed aggregates $k$ is given by $k \sim 2 k_{\rm mat} {( r_{\rm c} / r_{0} )} \phi^{(1 + \alpha) / (3 - D_{\rm f})}$, where $k_{\rm mat}$ is the material thermal conductivity, $r_{\rm c}$ is the contact radius of constituent particles, and $\phi$ is the filling factor of dust aggregates.


Introduction
The study of the aggregation of small dust particles into larger aggregates is crucial for understanding the fundamental processes in astrophysics and geophysics. For example, the growth of aerosol or haze particles in the atmosphere causes the scattering and absorption of the sunlight [e.g., 1]. In addition to this, the aggregation of dust particles also occurs in the mineral clouds of exoplanets and hence understanding the fundamental processes involved in the formation process of dust aggregates in exoplanets is imperative to interpret the transmission spectra [e.g., 2,3]. In addition, the aggregation of dust particles in the solar nebula is the first step towards the formation of the planets [e.g., [4][5][6][7], and the density evolution of dust aggregates is considered to be the key to understanding the evolution from nm-or µm-sized dust grains to km-sized small bodies [e.g., [8][9][10][11]. The resulting aggregates frequently have a complex and fractal structure with an extremely low filling factor [e.g., [12][13][14]. Therefore, understanding the geometrical structure of these fractal aggregates and its influence on the physical properties such as the thermal conductivity and the compressive strength is of immense current interest.
c The Author(s) 2012. Published by Oxford University Press on behalf of the Physical Society of Japan. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by-nc/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
For porous dust aggregates composed of micron-sized SiO 2 glass grains, the thermal conductivity is obtained by several experimental studies [e.g., 15,16], and it is empirically known that the thermal conductivity is approximately proportional to the square of the filling factor [17][18][19]. However, the theoretical explanation of the dependence of thermal conductivity on the filling factor is still lacking.
In contrast to the thermal conductivity, the filling factor dependence of the tensile strength of porous dust aggregates is well understood [20]. The tensile strength of porous dust aggregates is evaluated from the fractal structure of dust aggregates and the maximum force required to separate two sticking particles. The compression strength would also be evaluated from the fractal structure of dust aggregates and the rolling energy needed to rotate a constituent particle around its connecting points [21].
In this paper, we describe the geometrical structure of porous dust aggregates formed in astrophysical environments. We present the calculation of the gyration radius (which is defined in Section 2.1) and the graph-based geodesic radius (which is defined in Section 2.2) of porous dust aggregates. Subsequently, we present the interpretations of the filling factor dependence of thermal conductivity from the geometrical structure, and we also confirm the validity of our theoretical understanding by comparing it with the result of direct numerical calculations in Section 3. Finally, we also present the modified interpretation of the filling factor dependence of compression strength and the average coordination number of dust aggregates in Section 4.

Ballistic cluster-cluster aggregation
In the early stage of dust growth in astrophysical environments such as protoplanetary disks and circumplanetary disks, the collision velocity is sufficiently low to avoid collisional compaction and collisions between similar-sized dust aggregates are dominant [e.g., 22,23]. Therefore, the shape of dust aggregates in astrophysical environments is expected to resemble the clusters of ballistic cluster-cluster aggregation (BCCA). The BCCA clusters are formed by the sticking of two equal sized BCCA clusters with no restructuring (see Figure 3(a) of Okuzumi et al. [24]). We prepare BCCA clusters in the following procedure: (1) prepare an aggregate composed of 2 i particles (initially i = 0), (2) copy this aggregate and change the orientation of the copy randomly, (3) by ballistic sticking of these two aggregates with a randomly chosen offset, make a new aggregate composed of 2 i+1 particles, (4) continue the procedure (1)- (3). As shown in Figure 1, a BCCA cluster has a highly porous structure. In this study, we assumed that all constituent particles (hereinafter referred to as "monomers") are spherical and have the same radius r 0 .

Gyration radius
For the quantitative analysis of the structure of porous aggregates, we need to define a typical cluster radius. Here, we use the gyration radius R gyr , which is customary in aggregation studies [e.g., [23][24][25] defined by  where N is the number of constituent monomers, x i is the coordinate of the i-th monomer, and x o is the coordinate of the center of mass. We carried out 20 growth sequences of N -body simulations of BCCA as previous studies [24,25]. Here we show the geometric mean of the gyration radius R gyr as the function of the number of monomers N in Figure 2(a). We found that the gyration radius R gyr is given by and given uncertainties are the standard errors. The structure of BCCA clusters is therefore described in terms of the fractal dimension D f , which is defined as Our numerical data shows that D f is or 1.85 < D f < 1.92 when we take uncertainty into consideration. This result is consistent with previous studies [24][25][26]. The effective volume of BCCA clusters V is evaluated as V ∼ 4πR gyr 3 /3 and the volume of monomers is V 0 = 4πr 0 3 /3. The filling factor of BCCA clusters φ is given by Therefore, we can calculate the filling factor of BCCA clusters from the number of monomers.

Graph-based geodesic radius
Granular materials and dust aggregates transmit compressive stresses via a network of force chains [e.g., 27]. Further, the chains of monomers also conduct heat [e.g., 28]. Therefore, understanding the chain structure within dust aggregates is essential, especially for the study of the mechanical and heat transfer properties. Here, we introduce the graph geodesic and the 3/14 Schematic description of the distance between i-th and j-th particles, |x i − x j |, and the graph geodesic between i-th and j-th particles, d i,j . It is clear that the graph geodesic d i,j is larger than the distance |x i − x j |.
graph-based geodesic radius. Figure 3 schematically illustrates a BCCA cluster. The distance between i-th and j-th particles is |x i − x j |, and we define the graph geodesic between i-th and j-th particles as d i,j in Figure 3. Considering the graph structure of BCCA cluster, which is a tree (i.e., a connected acyclic graph), we can uniquely determine the graph geodesic d i,j . The typical length of the chain of monomers is obtained using the same method as the definition of R gyr . We define the graph-based geodesic radius R geo as It is clear that d i,j 2 ≥ (x i − x j ) 2 and then R geo ≥ R gyr by definition. Here we show the geometric mean of the graph-based geodesic radius R geo as a function of the number of monomers N over 20 runs in Figure 2(b). We found that the graph-based geodesic radius 4/14 R geo is given by log 10 R geo r 0 = (0.710 ± 0.013) log 10 N + (0.034 ± 0.007), (7) and given uncertainties are the standard errors.
Here, we consider the ratio of R geo and R gyr . The ratio of R geo and R gyr is given by where α is the dimensionless constant and the constant α must depend on the aggregation process of clusters. For BCCA clusters, we found that or 1.29 < α < 1.39 when we take uncertainty into consideration.

Bifractality of statically compressed BCCA clusters
In the early stage of dust growth, the fractal dimension of dust aggregates is D f 1.9. When the dust aggregates grow into cm-sized cluster, BCCA clusters are dynamically compressed by dust-dust collisions [e.g., 29] and/or statically compressed by ram pressure of the disk gas [e.g., 30]. Although it depends on the physical properties of the disk and the monomers, the compression mechanism for icy aggregates composed of submicron-size monomers in the minimum mass solar nebula [31,32] is the static compression by ram pressure [30]. In this study, we focus on the geometrical structure of statically compressed BCCA clusters. The geometrical structure of statically compressed BCCA clusters is characterized by bifractality [21]. Kataoka et al. [21] calculated the average number of particles in spheres of radii r in , N in . For statically compressed BCCA clusters, N in is approximately given by and the transition radius r tr is evaluated as r tr ∼ φ 1/(Df −3) r 0 . In other words, the fractal dimension becomes three on a large scale, while it remains 1.9, i.e., D f of BCCA, on a small scale. This structure evolution suggests that the static compression reconstructs the fractal aggregate first on a large scale, because of the weak compressive strength on a large scale [21]. Therefore, we can imagine that it is possible to understand the physical properties of statically compressed BCCA clusters from the geometrical structure of small BCCA clusters preserved in compressed aggregates (hereinafter referred to as "BCCA cells"). Figure 4 shows the schematic description of a BCCA cell in the compressed BCCA cluster.
It is important to note that the geometrical structure of dynamically compressed BCCA clusters is also characterized by bifractality [e.g., 6]. The resulting fractal dimension is approximately 2.5 on a large scale and it remains D f of BCCA on a small scale. Therefore, bifractality is a common characteristic of compressed BCCA clusters. We also hypothesize that this bifractality is a common structure for compressed fractal aggregates although the initial cluster is not originated from BCCA but other aggregation processes, for example, diffusion-limited cluster aggregation [e.g., 33] or reaction-limited cluster aggregation [e.g., 34]. We will, however, need to confirm this hypothesis in the future.

Thermal conductivity
In this section, we calculate the thermal conductivity of compressed BCCA clusters and demonstrate the manner in which the geometrical structure affects the thermal conductivity.

Methods
We calculate the thermal conductivity of compressed BCCA clusters composed of 16384 (= 2 14 ) monomers. The snapshots used in this study and used in our previous study (Arakawa et al. [19]) are the same and were prepared by Tatsuuma et al. [20]. The methods of the thermal conductivity calculation are described in our previous studies [18,19], which we briefly summarize it here.
Dust aggregates are statically compressed in a cubic periodic boundary (see Figure 1 of Arakawa et al. [18]). We consider one-directional heat flow from the lower to the upper boundary plane. The thermal conductivity of a dust aggregate in a cubic periodic boundary k is given by where f is a dimensionless function of φ, k mat is the material thermal conductivity, and r c is the contact radius of monomer grains. The normalized thermal conductivity f is given by where L is the length of the side of the cube, S = L 2 is the area of the upper and lower boundaries, T i is the temperature of i-th monomer, and ∆T is the temperature difference between the upper and lower boundaries. We took the sum of contacts between the i-th grain on the upper boundary and j-th grain inside the boundaries (see Arakawa et al. [18] for details).
In this study, we also consider the series connection of dust aggregates in a cubic periodic boundary (see Figure 5). It is expected that the series connection of dust aggregates would reduce the artificial effects of the boundary condition on the thermal conductivity calculations.   The aggregate used in this calculation is the same that used in Figure 6. We calculated f ∞ and σ ∞ from three directions (x, y, and z).

Filling factor dependence
It is predicted that the normalized thermal conductivity of dust aggregates in a series connection of n cubes, f n , is given by where f ∞ is defined as We can rewrite Equation (14) as and we found that (f 1 /f n ) − 1 is approximately proportional to 1 − (1/n). In Figure 7(a), we confirmed that the relation between (f 1 /f n ) − 1 and 1 − (1/n) works well. Therefore, we can evaluate f ∞ using f 4 and f 8 as follows: Figure 8(a) shows the normalized thermal conductivity for the limiting case of n → ∞, f ∞ , as a function of the filling factor φ. We used 10 snapshot data for each φ obtained from different compression simulations [20] and calculated f ∞ from three directions. The circles represent the geometric mean of 30 calculation results of the temperature structure, with vertical error bars of twice the standard error. We found that f ∞ is given by log 10 f ∞ = (2.068 ± 0.034) log 10 φ + (−0.022 ± 0.007), (18) and given uncertainties are the standard errors. This result is consistent with those of previous studies [18,19].

Surface density of heat paths
The thermal conductivity of dust aggregates must be affected by the typical length of the chain of monomers R geo and the surface density of heat paths σ. Here, we introduce the number of heat paths at the temperature T , N path (T ). We define N path (T ) as the number of contacts between two monomers whose temperatures are T i and T j with T i < T < T j . Thereafter, the average number of heat paths N path is given by where the temperature at the upper and lower boundaries are −∆T /2 and +∆T /2, respectively. The surface density of heat paths σ is given by The average number of heat paths N path depends on the number of connected cubes, n. In Figure 5, N path = 3/2 for the case of n = 1 and N path = 409/281 for the case of n = 2. As well as f ∞ , we evaluated the surface density of heat paths for the limiting case of n → ∞, σ ∞ . In Figure 7(b), we confirmed that σ ∞ can be predicted as follows: where σ n is the surface density of heat paths within a dust aggregate in a series connection of n cubes. Figure 8(b) shows the surface density of heat paths for the limiting case of n → ∞, σ ∞ , as a function of the filling factor φ. We found that log 10 σ ∞ r 0 −2 = (1.775 ± 0.025) log 10 φ + (0.076 ± 0.005), (22) and given uncertainties are the standard errors.

9/14
For a BCCA cell, we can imagine that the average number of heat paths is Thereafter, the surface density of heat paths within a BCCA cell is given by The relation between the gyration radius R gyr and the filling factor φ is and we obtain the relation between σ and φ: We found that 1.85 < D f < 1.92 in Section 2.1, therefore, we obtain 1. We note that the tensile strength of compressed BCCA clusters, P t , is approximately given by P t ∼ σF c , where F c = 3πγr 0 /2 is the maximum force required to separate two sticking monomers and γ is the surface energy [5]. Therefore, the tensile strength is given by which is consistent with the numerical result of Tatsuuma et al. [20], P t 0.6γr 0 −1 φ 1.8 . The coincidence of the filling factor dependence may indicate not only the number of heat paths but the number of force chains is also on the order of unity for BCCA cells.

Understanding the filling factor dependence of the thermal conductivity
Here, we demonstrate the manner in which the filling factor dependence of the thermal conductivity is derived from the geometrical structure. For compressed BCCA clusters, the fractal dimension is three on a large scale, then the thermal conductivity of compressed BCCA clusters should be the same as the thermal conductivity of BCCA cells. The spatial scale of BCCA cells is L ∼ R gyr and the area of the BCCA cells is S ∼ R gyr 2 , where R gyr ∼ N 1/Df r 0 is the gyration radius of a BCCA cell and N is the number of monomers in a BCCA cell. The surface density of heat paths is approximately given by σ ∼ R gyr −2 . The typical temperature difference between two contacting monomers, δT , is also given by where ∆T is the temperature difference between the upper and lower region of a BCCA cell.
The heat conductance at the contact of two monomers, H, is [e.g., 35,36] and the heat flow at the contact of two monomers, I, is I ∼ HδT . Therefore, the heat flow density within the BCCA cell is 10/14 and the thermal conductivity k is rewritten as follows: The normalized thermal conductivity f of the BCCA cell (and the compressed BCCA cluster) is therefore given by The relation between N and φ is N ∼ φ Df /(Df −3) , then we obtain the following equation: The derived relation shows excellent coincidence with our numerical result, f ∼ φ 2.068±0.034 .

Reinterpretation of the filling factor dependence of the compressive strength
We can also derive the filling factor dependence of the compressive strength of compressed BCCA clusters from the geometrical structure. In this section, we evaluate the compressive strength P c as Kataoka et al. [21] did. The compressive force on the surface area of the BCCA cell F c and the compressive strength P c are given by The length of the force chain within the BCCA cell is R geo . Since the compression is accompanied by the rolling of pairs of monomers in the force chain, the work required for compression can be given by where E roll = 6π 2 γr 0 ξ cr is the energy needed to rotate a monomer around its connection point by π/2 rad called the rolling energy, and ξ cr is the critical rolling displacement [5]. Subsequently, we found that the compressive strength P c is given by and (2 + α)/(3 − D f ) 2.99 for compressed BCCA clusters. The derived relation shows an excellent agreement with the numerical results of Kataoka et al. [21], i.e., P c ∼ (E roll /r 0 3 )φ 3 .
We note that the original explanation by Kataoka et al. [21] might not be accurate. Kataoka et al. [21] evaluated the work required for compression as and the filling factor dependence of the compressive strength was obtained as This estimate was based on the assumption that the compression is accompanied by the rolling of single pair of monomers in a BCCA cell. In this derivation, 3/(3 − D f ) 2.69 and it might not reproduce their numerical results. Although our findings suggest that the α parameter associated with the chain length plays a significant role on the compression of dust aggregates, further studies on the force distribution within compressed fractal aggregates are required. 11/14

Revisiting the average coordination number of compressed aggregates
The average coordination number (i.e., the average number of contacts per monomer) Z increases as an aggregate is compressed. Arakawa et al. [19] found that, for compressed BCCA clusters, the filling factor dependence of Z is given by Z = 2 + 9.38φ 1.62 . Here, we derive this equation from the geometrical structure.
Considering the graph structure of BCCA cluster, which is a tree (i.e., a connected acyclic graph), the average coordination number of a compressed BCCA clusters is where the constant C is the number of the inter-cell contacts per BCCA cell. The number of the faces, edges, and corners within a cube is 6, 12, and 8, respectively. Subsequently, we assume that the number of the inter-cell contacts per BCCA cell is C ∼ 9.
We define the deviation of the coordination number from two, ζ ≡ Z − 2. The deviation ζ is given by and 1.60 < D f /(3 − D f ) < 1.79 when we take the uncertainty of D f into consideration. Therefore, we obtain the following equation: which is consistent with the numerical result of Arakawa et al. [19], although the uncertainty of D f /(3 − D f ) is large and future studies on both the fractal dimension analysis and the average coordination number are essential. We note that the fractal dimension D f depends on the formation process of dust aggregates. Thereofre, the average coordination number Z also depends on the formation process of dust aggregates, as reported in Seizinger and Kley [37]. The compressive strength P c is also affected by the average coordination number Z. If the average coordination number is Z 2, nearly all the monomers can roll when they are compressed. Therefore, the interparticle force is close to the rolling friction force and the compressive strength is given by Equation (36). On the other hand, in the high-density region (φ 0.1 and Z 2), most of the particles cannot roll freely and the compressive strength is larger than the evaluated value for the case of Z 2 [e.g., [38][39][40][41]. Then, we expect that the compressive strength would be given by the sliding friction force in the high-density limit [40], although future studies are required to understand this in detail.

Summary
In this study, we conducted the numerical simulations of the BCCA of small dust particles and calculated the geometrical structure of the fractal dust aggregates. Additionally, we derived the filling factor dependence of the physical properties of porous dust aggregates.
Our key findings are summarized as follows.
(1) We calculated the gyration radius R gyr and the graph-based geodesic radius R geo as the functions of the number of constituent particles N . We found that R gyr /r 0 ∼ N 0.531±0.011 and R geo /r 0 ∼ N 0.710±0.013 , where r 0 is the radius of constituent particles. Thereafter, we defined two constants which characterize the geometrical structure 12/14 of fractal aggregates: D f and α. The definition of D f and α are N ∼ (R gyr /r 0 ) Df and R geo /r 0 ∼ (R gyr /r 0 ) α , respectively. We revealed that D f 1.88 and α 1.34 for BCCA clusters. (2) Kataoka et al. [21] found that the geometrical structure of statically compressed BCCA clusters is characterized by bifractality. This structure evolution suggests that the static compression reconstructs the fractal aggregate first on a large scale because of the weak compressive strength on a large scale. Therefore, we can imagine that it is possible to understand the physical properties of statically compressed BCCA clusters from the geometrical structure of small BCCA clusters preserved in compressed aggregates ("BCCA cells"). (3) We investigated the filling factor dependence of thermal conductivity of statically compressed aggregates. We found that the filling factor dependence can be interpreted from the geometrical structure of dust aggregates. The thermal conductivity of statically compressed aggregates k is given by k ∼ 2k mat (r c /r 0 )φ (1+α)/(3−Df ) , where k mat is the material thermal conductivity, r c is the contact radius of constituent particles, and φ is the filling factor of dust aggregates. (4) The compressive strength P c is also derived from the geometrical structure as P c ∼ (E roll /r 0 3 )φ (2+α)/  , where E roll is the energy needed to rotate a monomer around its connection point by π/2 rad. Our finding suggests that the α parameter associated with the chain length plays a significant role in the compression of dust aggregates. In addition, the average coordination number Z is derived from the geometrical structure as Z = 2 + Cφ Df /(3−Df ) , where C ∼ 9 is the number of the inter-cell contacts per BCCA cell.