- Email: [email protected]
Grouping and adding method for calculating light scattering by large fluffy aggregates
ARTICLE IN PRESS
Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80 www.elsevier.com/locate/jqsrt
Grouping and adding method for calculating light scattering by large fluffy aggregates Y. Okadaa,, T. Mukaia, I. Manna, H. Nomurab, T. Takeuchia, I. Sanoc, S. Mukaic a
Graduate School of Science and Technology, Kobe University, 1-1 Rokkodai, Nada, Kobe 657-8501, Japan b School of Mathematics and Physics, Queen’s University Belfast, Belfast BT 1NN, Northern Ireland, UK c Faculty of Science and Technology, Kinki University, 3-4-1 Kowakae, Higashi-Osaka 577-8502, Japan Received 1 August 2006; received in revised form 25 February 2007; accepted 5 April 2007
Abstract We present a method to derive the light scattering properties of very porous fractal aggregates composed of a large number of monomers where the size parameter of monomer is larger than unity. Our new method is based on the grouping of the aggregate: The aggregate is divided into groups, where each group is located along a line of the incident light, and the scattering properties of the group are calculated taking into account multiple scattering with monomers located inside the group, as well as those in a buffer region around the group. The scattering and absorption efficiencies are obtained by adding the resultant scattering properties for all the groups. We have shown that the method effectively works when the monomer scatters the incident lights predominantly in the forward direction, which is the case if the monomer size is large compared to the wavelength of the incident light. The errors in resulting scattering and absorption efficiencies for porous aggregates are investigated for various refractive indices and sizes of monomers. We found that the errors are larger for low absorbing materials and they can be reduced by expanding the buffer region. In the case of the buffer region for each group consisting of 1/8 of the total number of monomers, the results show errors less than 15% and 10% for absorption and scattering, respectively. It is also shown that the errors have a small standard deviation (i.e., 2%) for different directions of the incident light. r 2007 Elsevier Ltd. All rights reserved. Keywords: BCCA; Cluster T-matrix method; Dominant forward scattering; Scattering and absorption efficiencies; Non-Rayleigh region
1. Introduction Porous irregular particles are ubiquitous in space and in the atmospheres of the Earth and other planets. Since the 1970s it has been obvious that data on observations of dust in the solar system cannot be reproduced by applying the light scattering properties of compact spherical particles with the same volume [1]. For estimating the scattering properties of the irregular particles laboratory measurements of irregular particles were taken into account [2–5]. Also, the numerical studies of the formation of aggregates by sticking of particles are studied [6]. Corresponding author. Tel.: +81 78 803 6566; fax: +81 78 803 6483.
E-mail address: [email protected] (Y. Okada). 0022-4073/$ - see front matter r 2007 Elsevier Ltd. All rights reserved. doi:10.1016/j.jqsrt.2007.04.006
ARTICLE IN PRESS 66
Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
The concept of fractal aggregates has been applied to the calculations of the radiation pressure [7] and the temperature [8] of the interplanetary dust particles. For atmospheric models, biomass burning aerosol particles produced from the forest fires are often modeled as fractal aggregates [9,10]. From laboratory experiments of dust aggregation in microgravity environment it is shown that aggregates form from the sticking of particles [11,12]. In many cases, however, the concept of fractal aggregates has been applied to particles that are not produced by such a growth process. The fractal aggregates provide a concept for a reproducible description of irregular particles, which are actually produced by different processes [13]. Since light scattering properties are fundamental for various applications to interpret data observed by such as Earth-observing satellites, space-borne sensors and ground-based instruments observing dust from comets, asteroids and protoplanetary disks, it is important to have a through understanding of the light scattering of fractal aggregates. Many different methods have been proposed for the calculation of light scattering properties of nonspherical particles, such as DDA (discrete dipole approximation) [14,15], a1 term method [16,17], T-matrix for cluster of spheres [18,19], FDTD (finite difference time domain) method [20–22], EMT-Mie (effective medium theories combined with Mie theory) method [7,23], GO (geometrical optics) [24–27]. Calculations are especially difficult for aggregates composed of a large number of spherical monomers when the monomer size is larger than the wavelength of the incident light. The DDA method is based on the calculation of interactions of dipoles which describe the shape of nonspherical particles. To calculate for fractal aggregates composed of monomers, the shape of the aggregate should be described with a number of dipoles. In order to keep the accuracy, the product of dipole separation d, wavenumber k and refractive index jmj (i.e., jmjkd) should be less than 1.0 [14,15]. This requires the subdivision of each monomer into number of dipoles when the monomer size is larger than the wavelength of the incident light as conducted in [28,29] and the memory limitation leads to the calculation of aggregates composed of no more than a few hundred monomers. The FDTD method solves the Maxwell equations by calculating the time evolution of electromagnetic waves (EMs), and the shape of the particles is described using cubic cells. The resultant scattering properties in the time domain are then converted into the frequency domain. The FDTD method can be used to calculate the scattering properties of particles larger than that are computable with DDA [20–22,30]. However for the considered aggregates, the geometric limitations of the FDTD method mean that the shape should be finely discretized into a larger number of cubic cells in order to reproduce the spherical surface of each monomer [30]. This also requires a large amount of computer memory, and hence it is difficult to handle aggregates composed of a large number of spheres with the FDTD method. The EMT-Mie method is easy to implement since it first derives the effective refractive index based on the volume fraction of the material relative to empty space, and then calculates the light scattering properties based on Mie theory [7,23,31]. However, compared with other methods such as DDA, the accuracy of the EMT-Mie method deteriorates very rapidly as the monomers become larger than the wavelength of the incident light (i.e., non-Rayleigh region), especially for highly porous particles [7,31–37]. In contrast to most methods, which have upper limitations for the monomer size, the GO method has a lower limitation for the monomer size. The GO method is applicable only when the circumference of the monomer is much larger than the wavelength of the incident light (i.e., around 100 times) and it is based on treating the incoming light as a bundle of rays. Fresnel’s formulae are used to describe the reflection and refraction of rays [24–26]. However, there is still a large gap in the range of allowable sizes of monomers for porous particles between other methods (around a few times the wavelength) and the GO method (around 100 times) [38], indicating the need to modify the GO method to include the influence of wave optics or the need to expand the size range into larger one for other methods some how. The cluster T-matrix method developed by Mackowski and Mishchenko [18] is based on the superposition of the T-matrix for each sphere. The calculation is accurate as long as higher order of expansions of scattering function are calculated. However, this method is limited by the extensive computing resources required for a larger number of monomers and for larger monomers. The upper limitation with our computing resources (shown in Section 3) is around 600 monomers for the circumference of the monomer of 4.5 times the wavelength.
ARTICLE IN PRESS Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
67
As shown above, it is required to make a convenient method to derive the light scattering properties for the porous fluffy particles composed of thousands (or a larger number) of monomers where the circumference of the monomer is larger than the wavelength of the incident light (i.e., non-Rayleigh region) but smaller than the region of validity for the GO method. In this paper we try to expand the number of monomers in around 10 times with our proposed method applied to the currently available cluster T-matrix method of fixed direction of incident light [18]. In Section 2, scattering characteristics of single compact spheres are presented, then our proposed method is introduced based on the scattering characteristics of each monomer. In Section 3, errors are investigated for different sizes of buffer region where the higher order of multiple scattering is considered. Then, error analysis is conducted as a function of real and imaginary part of refractive indices and sizes of monomers. The error maps of absorption are compared with the maps of scattering properties of a single monomer. Also the errors for different directions of incident lights are investigated. In Section 4, examples of the light scattering properties of very porous fractal aggregates composed of up to 8192 monomers are presented. The result shows the power law relationship between the number of monomers and the efficiencies. In the conclusion, we summarize this study and also we briefly mention some possible future improvements combining with methods other than the cluster T-matrix method used in this study. 2. Method 2.1. Forward scattering of a single compact sphere For a model of the shape and structure of the scattering objects, we consider the BCCA (ballistic cluster–cluster aggregates) as shown in Fig. 1, where the constituent monomers have same radii of rm . The fractal dimension Df determined from the radius of gyration and the number of composing monomers is around 2.0 [7], which shows the high fluffiness of the aggregates. BCCAs have very large porosities (i.e., 99% of vacuum inside the characteristic radius) [7,31]. It is shown that EMT-Mie method yields a few tens of percentages of errors for the light scattering properties of very porous particles such as BCCA, when the monomer size becomes larger than the wavelength of the incident light (i.e., non-Rayleigh region), and the errors tend to be larger for larger particles (e.g., see Fig. 6 in [32], Fig. 3 in [33] and Fig. 9 in [35]). Table 1 summarizes the parameters of the light scattering simulations considered in this study. The size information is described with size parameter x ¼ 2pr/l, where r is the radius of the particle and l is the
Fig. 1. Three-dimensional image of a BCCA composed of 512 monomers viewed from a zenith angle of 151. Note that all the monomers are of the same size and chemical composition.
ARTICLE IN PRESS 68
Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
Table 1 Parameters of the light scattering simulations xm l N mono m ¼ n þ ik Qabs Qsca
Monomer size parameter: Wavelength of the incident light: Number of monomers in the aggregate: Refractive index of all the monomers: Absorption efficiency of the aggregate: Scattering efficiency of the aggregate: The Qabs and Qsca are defined for the volume-equivalent sphere radius.
p11 (Θ) / p11 (0°)
100 xm=0.001 1.0 5.0 10.0 50.0 100.0 500.0
10-2
10-4
10-6 0
30
60 90 120 Scattering angle Θ (degree)
150
180
Fig. 2. Scattering function p11 ðYÞ of a single compact sphere as a function of size parameter of the sphere, xm , where m ¼ 1:5+0:1i. Each scattering function is normalized to its value at Y ¼ 0 .
wavelength of the incident light. The radius of the aggregate is often taken as reff , corresponding to the radius of a sphere with the volume equal to that of the aggregate. Here we follow this custom when to the pffiffiffiffireferring ffi radius of the aggregate. We mainly consider the size parameter of the monomer xm ¼ xag / 3 N mono where xag represents the size parameter of the aggregate calculated from reff and l. N mono is the number of monomers. The chemical composition of the monomers composing the aggregate is considered to be homogeneous, and we study the light scattering as a function of refractive index m ¼ n þ ik. Of the various light scattering properties of the aggregates, we mainly consider the efficiency factors for absorption Qabs and for scattering Qsca . The angular distribution of the scattered light is described with a scattering function p11 ðYÞ, which denotes the intensities of scattered light as a function of scattering angle Y, where Y ¼ 0 represents the forward and Y ¼180 the backward direction relative to the direction of the incident light. Fig. 2 shows the scattering functions of single compact spheres normalized to the intensity at Y ¼ 0 for various particle sizes calculated using Mie theory. This figure clearly shows that the forward scattering component (e.g., Yo30 ) becomes more dominant relative to scattering into side direction (e.g., Y460 ) as the particle becomes larger. This characteristic forms the basis of the method proposed in Section 2.2.
2.2. Grouping and adding method The maximum computable number of monomers for the given size of monomers is limited in the cluster Tmatrix method (and in most other methods). As the number of monomers increases, the memory requirement increases very rapidly. At the same time, a larger number of equations have to be solved for the interactions of electric fields and the calculation can become unstable.
ARTICLE IN PRESS Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
69
Our method proposed in this section, called the grouping and adding method (GAM), has the potential to calculate the scattering properties of large fluffy aggregates beyond the computable limitations for the number of monomers, and involves the following steps: (1) Divide the aggregate into N DIV groups. (2) Calculate efficiencies (i.e., Qabs and Qsca ) of each group. (3) Sum the efficiencies of all the groups. The complete aggregate is first divided into groups, which consist of monomers close to one another on the two-dimensional plane (standard plane) perpendicular to the direction of the incident light as shown in Fig. 3a and a0 , where all constituent monomers are projected on this standard plane. In this study, the grouping is conducted using the K-means method [40] as follows: (1) Randomly choose N DIV monomers on the standard plane and each of which is used as a seed for each group. (2) Select other members in each group based on the shortest distance from the seed monomers for the group on the standard plane.
Top view 9 group 4 surrouding monomers
8 7
1
Y
Y
10
4 6 2 5
3
EM wave
EM wave X
X
Z
Z
Side view
EM wave X
EM wave X
Fig. 3. Schematics of grouping of the aggregate (a and a’) and of the buffer region (b and b’) for EMs propagating in the Z direction. The aggregates shown in Fig. 1 are grouped with N DIV ¼ 10. In panels b and b’, gray- and black-shaded circles are the N in gr monomers in the group and N out gr monomers in the buffer region, respectively, and the large solid circle and dashed one denote the maximum radius of group 4 and buffer region for the group, respectively. The direction of the incident light is shown in the right bottom of each panel.
ARTICLE IN PRESS Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
70
(3) Calculate the center of each group using the selected members together with the seed monomer for each group. (4) Regroup all the monomers based on the shortest distance from the centers of groups. (5) Repeat steps 3 and 4 until the changes in the locations of the centers of groups for all the groups become smaller than a pre-configured threshold. The result of the grouping method is not reproducible because the K-means method uses a random number generator to locate the first seed monomers. However, we have found that the different groupings generated for the same N DIV do not largely influence on the results. In the calculation of the light scattering properties of each group, a contribution from the monomers belonging to other groups cannot be neglected. To examine such a contribution, we use the N in gr monomers of the group and the N out monomers in the buffer region (see Fig. 3b and b’, and Table 2). The buffer region is set gr as that surrounding the selected group based on x and y coordinates on the standard plane, and not using information of z coordinates (i.e., the direction of the incident light). From the inner part of the buffer region, N out gr monomers are selected until the chosen number (e.g., N mono =8) are selected. When grouping is not used, the light scattering properties of the aggregate are calculated in the presence of all the monomers. In contrast, in our proposed method, only N out gr nearby surrounding monomers in the buffer region are used together with the N in gr monomers to calculate the scattering contributions of these monomers. in The N out gr monomers are used in order to partly account for the light that is scattered sideways by the N gr monomers. The grouping reduces the memory requirement for the calculation, from being proportional to out N mono to being proportional to (N in gr + N gr ). The light scattering properties of each group are calculated by considering a contribution from the monomers included in the group and its buffer region, as follows: out (1) The monomers classified into the group (N in gr ) and in the buffer region (N gr ) are used for the light scattering calculations for a particular group. (2) Only the light scattering properties of the N in gr monomers are stored as the result for the selected group. (3) Steps 1 and 2 are repeated for all the groups. (4) The scattering properties from all the groups are simply added to deduce the total light scattering properties of the aggregate.
In the Superposition T-matrix method for clusters of spheres, cross sections (i.e., C ext and C abs ) are formulated as the superposition of those of composing spheres as the following: C ext ¼
Ns X
C iext ,
(1)
C iabs ,
(2)
i¼1
C abs ¼
Ns X i¼1
respectively from Eqs. (24) and (43) of [39]. Where C ext , C abs are the cross sections of extinction and absorption of total cluster while C iext and C iabs are those of ith composing sphere. Ns is the number of particles. In the above equations, multiple interactions of the EM wave between the particles are included into C iext and C iabs . Table 2 Parameters used in the GAM Number of groups: Number of monomers in each group: Number of monomers in the buffer region:
N DIV N in gr N out gr
ARTICLE IN PRESS Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
71
In the Fortran code developed by Mackowski et al., C ext and C abs are calculated in two ways (see ftp:// ftp.eng.auburn.edu/pub/dmckwski/scatcodes/index.html): (1) summation of the individual sphere results, (2) analytical integration of the cluster scattered field expansion. Then the comparison of calculations (1) and (2) have been checked to agree if the interaction equations are solved to converge. In our paper, we use (1) summation of the individual sphere results. Currently only Qabs and Qsca are calculated based on our method. The phase matrix are not just additive from those of each composing particle. Therefore in this paper, we will not go further to use our method for phase matrix. The cluster T-matrix code outputs the efficiency factors for extinction and absorption of the TM (transverse magnetic) and TE (transverse electric) mode of the EMs for each monomer. These efficiencies are influenced by (a) the incident electric field and (b) the electric fields induced by monomers other than the selected one, due to the multiple scattering by monomers. For each monomer, we average the efficiencies for the TM and TE modes to calculate efficiencies for non-polarized incident light. The grouping increases the computing time N DIV fold. However the grouping also works for the reduction out of time by the reduction of monomers from N mono to (N in gr + N gr ) for each group. Moreover, once the grouping is conducted, the calculation for each group can be parallelly performed in different computers, resulting in the faster calculations. 3. Error analysis The GAM is based on disregarding the scattering by monomers outside the combined area of the group plus its buffer region. Neglecting the contribution from outer monomers in this way will introduce errors into the scattering properties (i.e., Qsca and Qabs ). To ensure the efficacy of the GAM, we need to investigate the associated errors caused under various conditions (i.e., various ðn; kÞ and xm ). Although the GAM would be applied for a large fluffy aggregate, where the total number of monomers is beyond the limit of available computer system, we need to conduct the investigation of its applicability under the computable number of monomers. In this section, the results of the GAM for small aggregates composed of 512 monomers are compared to those calculated without grouping. The errors in Qabs and Qsca are defined as Errorð%Þ ¼
jQGAM Qorg j 100, Qorg
(3)
where QGAM denotes the value based on the GAM while Qorg is the value without grouping. The 600 monomers for xm ¼ 4:5 is the upper limit for our computing resources (i.e., Pentium 4 2.8 GHz CPU with 2 GB of RAM). Therefore, we investigate the errors under this upper limitation. 3.1. Errors for different sizes of buffer region for silicate, pyroxene and amorphous carbon Fig. 4 shows the errors in Qabs and Qsca for different sizes of buffer region in the cases of five different shapes of aggregates. Only one direction of incident light is considered. We investigate errors for silicate, pyroxene and amorphous carbon, which exhibit low, moderate and high absorbing materials, respectively. Although pyroxene belongs to the silicate group, as a matter of practical convenience, we differentiate them as ‘‘silicate’’ and ‘‘pyroxene’’ in this paper. We use the refractive indices at a wavelength of 0:55 mm for these materials (see Table 3). For the pyroxene, we use the data for Mgx Fe1x SiO3 , where x ¼ 0:7 [42]. For the silicate and pyroxene (see Fig. 4a–d), the errors in Qabs decrease with increasing buffer region (i.e., N out gr to N mono ratio), whereas the errors in Qsca for the two materials do not show clear trends. However, the errors in Qsca are already small (i.e., less than 10%) for most of the results. The absorption occurs when the light interacts with each monomer, and hence the errors are directly influenced by the number of side monomers outside the buffer region that are disregarded.
ARTICLE IN PRESS Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
20
20
10
0
0 1/2
1
0 1/16 1/8
5
1/4
5
Nout gr /Nmono
Nout gr /Nmono 20 shape 1 2 3 4 5
10
0 1
1/2
1/2
0
0 1/16 1/8
5
1/4
5
1/4
Err (Qsca) (%)
15
10
0 1/16 1/8
Nout gr /Nmono
Nout gr /Nmono
20 shape 1 2 3 4 5
15
10
Nout gr /Nmono
1
0 1/2
0 1/4
5
0 1/16 1/8
5
1/4
10
0 1/16 1/8
Err (Qsca) (%)
15
shape 1 2 3 4 5
1
20
1/2
Err (Qabs) (%)
15
shape 1 2 3 4 5
1
20
Err(Qabs) (%)
1/2
Err (Qsca) (%)
15
10
0 1/16 1/8
Err (Qabs) (%)
15
shape 1 2 3 4 5
1
shape 1 2 3 4 5
1/4
72
Nout gr /Nmono
Fig. 4. Errors in Qabs and Qsca calculations as functions of the number of monomers (N out gr ) included in the buffer region relative to the total number of monomers (N mono ). Results are for silicate (a and b), pyroxene (c and d), and amorphous carbon (e and f), where N mono ¼ 512, xm ¼ 3:0 and N DIV ¼ 60.
ARTICLE IN PRESS Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
73
Table 3 Refractive indices of the materials at a wavelength of 0:55 mm Silicate Pyroxene Amorphous carbon
m ¼ 1:48 þ 0:000028i [41] m ¼ 1:635 þ 0:0042i [42] m ¼ 2:03 þ 0:77i [43]
In contrast, Qsca corresponds to values integrated for solid angles, and so the errors depend on the directional distribution of the scattered energy. For example, when the light is scattered more in the forward direction, side monomers contribute less to the scattering efficiencies. Therefore, once sufficient monomers are considered in the buffer region, the errors in Qsca are little affected by the size of buffer region, as shown in Fig. 4b and d. Fig. 4a–d shows that the errors are higher in Qabs than in Qsca when the buffer region is small. Therefore, the N out gr -to-N mono ratio can be set based mainly on the errors in Qabs for materials with low to moderate absorption. The errors in both Qabs and Qsca for highly absorbing amorphous carbon (see Fig. 4e and f) are less than 10% for the sizes of buffer region investigated. The errors in Qabs does not decrease when the size of the buffer region is increased. This is because the lights interacting with amorphous carbon are highly absorbed by each monomer. This results in multiple scattering with higher energies being confined to a smaller region compared to those of silicate and pyroxene. It is considered that the buffer regions investigated include the large region where most multiple scatterings with higher energies for amorphous carbon take place. Fig. 4 indicates that the N out gr -to-N mono ratio needs to be set as at least 1/8 in order to keep the errors in Qabs and Qsca mostly below 10% for xm ¼ 3:0. This allows us to treat 8 times 600 monomers with errors of less than 10% for these three materials. 3.2. Error maps for ðn; kÞ and xm Fig. 4 shows the results of the errors of different sizes of buffer regions for only three materials. In order to systematically investigate the errors caused by the grouping, we have investigated the errors for various sets of ðn; kÞ and xm . Fig. 5 shows error maps for n vs k for different xm . The N out gr -to-N mono ratio for the buffer region is set as 1/8 for this figure. Only one direction of the incident light is considered. The variations in errors for different directions of the incident light are discussed in Section 3.4. For the results of xm ¼ 1:0, errors are around 5% for all sets of n and k. The errors of xm ¼ 1:0 are smaller than those of xm ¼ 3:0 and xm ¼ 4:5. The reason of this smaller errors for xm ¼ 1:0 will be shown in Section 3.3. For xm ¼ 3:0 and xm ¼ 4:5, the errors in Qabs show certain trends for different ðn; kÞ and xm . For smaller n (i.e., np1:30), the errors in Qabs is around 5% for all sets of xm and k. The real part of refractive index (i.e., n) denotes how much the light is scattered. The smaller n results in less scattering of the light and that corresponds to the smaller errors in Qabs . On the other hand, for nX1:4, errors in Qabs become larger for smaller k (i.e., ko0:01). This is because a smaller k indicates less absorption by each monomer, with the resulting higher order of multiple scattering into the side direction outside the buffer region causing larger errors. The error maps of Fig. 5 show that the grouping with an N out gr -to-N mono ratio of 1/8 causes errors in Qabs of less than 15% for ko0:01 and nX1:4, and less than 10% for the remaining parameter combinations. The errors in Qsca for all sets of n and k in Fig. 5 are less than 10%. This indicates that the use of the GAM allows us to expand the number of monomers eight-fold within these errors for all sets of n and k for xm X1:0. 3.3. Errors related to the scattering properties of each monomer In Fig. 5, the errors in Qabs show certain trends (i.e., higher errors for larger n and for smaller k). These trends in errors are considered to have relations with how each monomer scatterers the light.
ARTICLE IN PRESS 74
Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
Errors in Qabs (%)
Errors in Qsca (%)
xm=1.0
0.00001 0.00001 0.00005 0.00005 0.0001 0.0001 0.0005 0.0005 0.001 0.001 1.1 1.1 0.005 0.005 1.2 1.2 0.01 0.01 k k 1.3 1.3 1.4 1.4 0.05 0.05 1.5 1.5 0.1 0.1 1.6 1.6 1.7 1.7 0.5 0.5 n n 1.8 1.8 1.9 1.0 1.9 1.0 2.0 2.0
xm=3.0
0.00001 0.00001 0.00005 0.00005 0.0001 0.0001 0.0005 0.0005 0.001 0.001 1.1 1.1 0.005 0.005 1.2 1.2 0.01 k 1.3 0.01 k 1.3 1.4 1.4 0.05 0.05 1.5 1.5 0.1 1.6 0.1 1.6 n 1.7 n 1.7 0.5 0.5 1.8 1.8 1.9 1.0 1.9 1.0 2.0 2.0
xm=4.5
0.00001 0.00001 0.00005 0.00005 0.0001 0.0001 0.0005 0.0005 0.001 0.001 1.1 1.1 0.005 0.005 1.2 1.2 0.01 k 1.3 1.3 0.01 k 1.4 1.4 0.05 0.05 1.5 1.5 0.1 1.6 1.6 0.1 n 1.7 1.7 n 0.5 0.5 1.8 1.8 1.9 1.0 1.9 1.0 2.0 2.0 Errors (%)
0
5
10
15
20
Fig. 5. Maps of errors in Qabs and Qsca for n vs k for different xm . The buffer region has an N out gr -to-N mono ratio of 1/8, and N mono ¼ 512 and N DIV ¼ 60.
ARTICLE IN PRESS Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
xm=1.0
1.1
1.2
1.3
1.4
1.5
n
1.6
1.7
1.8
1.9
2.0
75
xm=3.0
0.00001 0.00005 0.0001 0.0005 0.001 0.005 0.01 k 0.05 0.1 0.5 1.0
1.1
1.2
1.3
1.4
1.5
n
1.6
1.7
1.8
1.9
2.0
0.00001 0.00005 0.0001 0.0005 0.001 0.005 0.01 k 0.05 0.1 0.5 1.0
xm=4.5
1.1
1.2
1.3
1.4
1.5
n
1.6
1.7
1.8
1.9
2.0
0.00001 0.00005 0.0001 0.0005 0.001 0.005 0.01 k 0.05 0.1 0.5 1.0
mono Qmono sca /g
0
2
4
6
8
10
Fig. 6. Scattering properties of a monomer for n vs k for different xm . Qmono and gmono , which are the scattering efficiency and asymmetry sca parameters of the monomer, respectively, are calculated using Mie theory.
Fig. 6 presents the scattering efficiencies of the single monomer (Qmono sca ) divided by the asymmetry parameter of the monomer (gmono ), based on Mie theory. The asymmetry parameter g is defined as [44] Z p11 ðyÞ cos y dO, (4) g ¼ o cos y4 ¼ 4p
where p11 ðyÞ is the scattering function for the scattering angle y. O is the solid angle for the integration. The value of g approaches 1 as forward scattering becomes dominant, and hence the 1=g represents how much light is scattered sideways. Qmono represents how effectively the incident light is scattered, sca mono and so the product Qmono indicates how widely the light is distributed after scattering by each sca =g monomer. Therefore, the trends seen in Fig. 6 coincide with those in Fig. 5; namely, larger errors in Qabs by the GAM are caused by a lower degree of forward scattering by the monomer. This indicator can therefore be used to quantify which combination of n, k, and xm results in larger errors. For xm ¼ 1:0, mono the Qmono shows smaller values without trends for n vs k, corresponding to the smaller errors shown in sca =g Fig. 5.
ARTICLE IN PRESS 76
Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
20
20
15
mean : 6.71(%), SD : 2.44(%) Err (Qsca) (%)
Err (Qabs) (%)
mean : 9.09(%), SD : 2.09(%)
10 5 0
15 10 5 0
0
100 200 ID of direction
300
0
100 200 ID of direction
300
Fig. 7. Errors in Qabs and Qsca of a BCCA for 324 directions of incident light, where xm ¼ 3:0, N mono ¼ 512, and N DIV ¼ 60 for silicate (i.e., m ¼ 1:48 þ 0:000028i). The N out gr -to-N mono ratio is set at 1/8.
3.4. Errors for the different directions of incident light The error analysis presented above considered incident light from only one direction for each shape of aggregate. The errors differ for different directions of incident light, because this causes differences in the alignment of monomers. Therefore, we also investigate variations in errors for different directions of incident light. The BCCA consisting of 512 monomers is rotated in azimuth and in polar angles at regular intervals to give a total of 324 directions. Fig. 7 shows the errors in Qabs and Qsca for 324 directions of incident light for xm ¼ 3:0. The average errors in Qabs and Qsca are less than 10% with a standard deviation of 2%. In general, the values of Qabs and Qsca for aggregates are deduced from averaging over those values obtained at different directions [8]. 3.5. Selection of N DIV The division of the aggregate into groups requires the light scattering of N DIV groups to be computed. We recommend that the N out gr -to-N mono ratio is set as 1/8, especially for materials with a low absorption. Then in N DIV can be taken so that (computable number of monomers—N out gr ) is assigned to N gr . In this paper, we have investigated errors for N DIV ¼ 60. The errors for smaller N DIV are either larger or smaller due to the following two reasons: out (1) The maximum radius of the N in gr monomers of each group increases as N DIV decreases. In this case the N gr monomers are assigned from monomers located over a larger radius with the thinner width, leading to multiple scattering between outer of N in gr monomers and those of buffer region being treated less. This results in the increase in the errors. (2) The increase of N in gr monomers for smaller N DIV causes higher-order multiple scattering in the group monomers being considered. This results in the decrease in errors.
Fig. 8 shows the errors for different numbers of groups (i.e., N DIV ¼ 72, 30, and 16) where the N out gr -toN mono ratio is fixed at 1/8. However, the errors in Qabs do not show clear trends for the increase of N DIV because of the complicated combination of the item 1 and 2 of the above. Then the N DIV can be selected based on the computing time. 4. Application to large aggregates composed of a number of monomers beyond the limit of computing resources The results of error analysis above indicates that the GAM can be utilized to increase the number of monomers beyond the computable number for single aggregate within certain errors. In this section, we
ARTICLE IN PRESS Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
Error in Qsca (%)
Error in Qabs (%) 15 10 5
15 NDIV=72 30 16
10
n=1.1
n=1.6
n=1.2
n=1.7
0 15 n=1.2
10
5
5
n=1.7
0 15
0 15 10
n=1.1
5
n=1.6
0 15 10
n=1.3
n=1.3
10
5
5
n=1.8
n=1.8 0 15
0 15 n=1.4
n=1.4
10
10 5
5
n=1.9
n=1.9
0 15
0 15
10
10
n=1.5
5 0 10-5
77
n=1.5 10-3
5
n=2.0 10-1 1 10-5 k
10-3
n=2.0
10-1 1
0
10-5
10-3
10-1 1 10-5 k
10-3
10-1 1
Fig. 8. Errors in Qabs and Qsca for different sets of n and k, where N mono ¼ 512 and xm ¼ 3:0. Different number of grouping (i.e., N DIV ¼ 72, 30 and 16) are investigated. The N in gr -to-N mono ratio is set as 1/8.
present the light scattering properties of large aggregates composed of up to thousands of monomers calculated with the GAM. Fig. 9 shows the Qabs and Qsca of very porous BCCAs composed of up to N mono ¼ 4096 for silicate and pyroxene, and those up to 8192 for amorphous carbon where xm is 4.5. The results of N mono p512 are calculated without grouping while N mono X1024 based on the GAM. The results of up to N mono ¼ 4096 are calculated with the N out gr -to-N mono ratio of 1/8. The results for the amorphous carbon for N mono ¼ 8192 are in out calculated with an N out gr -to-N mono ratio of 1/16 due to the upper limitation of (N gr + N gr ) being 600 for our computing resources. In this case of N mono ¼ 8192, the errors in Qabs and Qsca are expected to be around 10% for amorphous carbon (see Fig. 4e and f). In Fig. 9, only one direction of incident light is considered for each of the aggregates with different N mono . Therefore, the deviation of Qsca from others seen for N mono ¼ 128 and 256 are attributed to this rather than being considered representative of random orientation results. However, all of the results show the characteristic of increasing efficiencies for increasing N mono . We note that in contrast with the case of compact sphere, Qabs plus Qsca in Fig. 9 does not converge to 2 when the aggregate size increases. It is because the efficient factors of BCCA are defined using the radius of volume-equivalent sphere.
ARTICLE IN PRESS 78
Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
50
101
Silicate Pyroxene Amorphous carbon
30
10-1
Silicate Pyroxene Amorphous carbon
Qsca
Qabs
100
10-2 10-3
10 5 2
64
256 1024 4096 Number of monomers
64
256 1024 4096 Number of monomers
Fig. 9. Log–log plots of Qabs and Qsca of aggregates composed of up to 8192 monomers, where xm ¼ 4:5 for silicate, pyroxene and amorphous carbon. The lines are the regression line for each material. The aggregates with N mono p512 are calculated without grouping while those with N mono X1024 are calculated with the GAM for one direction of the incident light. Table 4 lists the equations of the regression lines.
Table 4 Equations of regression lines in Fig. 9 for silicate (Si), pyroxene (Py) and amorphous carbon (Ac)
Si Py Ac
Qabs
Qsca
y ¼ 0:000811 x0:292 (0.998) y ¼ 0:133 x0:287 (0.997) y ¼ 0:831 x0:298 (0.995)
y ¼ 2:35 x0:282 (0.990) y ¼ 1:82 x0:293 (0.994) y ¼ 0:864 x0:304 (0.996)
The numbers in brackets represent the correlation coefficient of the regression, and the x and y represent N mono and Qabs (or Qsca ), respectively.
Table 4 lists the equations for the regression lines in Fig. 9 and their corresponding correlation coefficients, all of which exceed 0.99. It is also notable that the regression lines have similar exponential slopes of around 0.293, commonly seen not only for all three materials but also both for Qabs and Qsca . The increase in the efficiencies with increasing N mono is related to the increase in multiple scattering by the monomers. Therefore, the structure of the fractal aggregates is considered to be related to the exponential slope of efficiencies for an increasing of number of monomers. However the slopes for aggregates with different fractal dimensions are likely to differ from those of BCCAs. 5. Conclusions We have proposed the GAM to derive the light scattering properties of a large fractal aggregate composed of spherical monomers where two conditions are met: (1) the scattering properties of each monomer can be derived and (2) in the light scattering by the monomer, the forward scattering is dominant. The studies for errors in the results computed by the GAM have been examined for the aggregates consisting of available number of monomers in our computer system, as functions of the number of groups, the sizes of neighbor region influenced on certain group (buffer region), sizes and optical constants of a constituent monomer. It is found that the GAM can keep the errors in Qabs and Qsca less than 15% and 10%, respectively, when the effect of 1/8 monomers of total constituent monomers N mono in the aggregate is taken into account for computing scattering properties of each monomer. This implies that the aggregate consisting of monomers eight times of available number of monomers in our computer system can be treated by GAM within these errors in Qabs and Qsca . In general, the GAM allows us to expand the number of monomers beyond the limit of computer resources. We have also found the power law relationship of Qabs and Qsca for N mono with the similar slopes for three different materials and correlation coefficients larger than 0.99.
ARTICLE IN PRESS Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
79
In this study, we have investigated the errors for BCCA. There is another fractal aggregate named BPCA (ballistic particle-cluster aggregate) often used in the light scattering calculation. BPCA has the fractal dimension of 2.980.02 [7]. Since the composing particles are located closer to each other than BCCA, the side scattering becomes more important. We have also investigated our method for BPCA. The errors are larger as expected (not shown here). We consider that our method works effectively only for porous particle described as BCCA. In this study, we used the cluster T-matrix code developed by Mackowski and Mishchenko [18]. However, the combinations with other methods are possible as long as the above conditions (1) and (2) are met. For example, the combination with the DDA method may allow us to calculate the asymmetry parameter g [45]. The combination with a1 term method may allow us to treat aggregates composed of a larger number than cluster T-matrix method for xm around 1.0 with errors by a1 term method and those by the GAM method [16,17]. The combination with the GMM method may allow us the calculation of aggregates of spheroidal monomers even though the applicability of the GAM method should be investigated for randomly directed spheroidal monomers [46,47]. Then these methods receive advantages in the increased number of monomers as conducted in this paper for cluster T-matrix method. The sample computer code of the GAM for cluster T-matrix method can be found at http:// harbor.scitec.kobe-u.ac.jp/okada/GAM/index.html. Acknowledgments This study used the cluster T-matrix code developed by Dr. Mackowski and Dr. Mishchenko (obtained from http://www.giss.nasa.gov/crmim/). Some of the data analyses were performed on the general common use computer system at the Astronomical Data Analysis Center (ADAC) of the National Astronomical Observatory of Japan. We thank Prof. Yoshitsugu Nakagawa at Kobe University, for comments that helped to improve this study. This research is supported by ‘‘The 21st Century COE Program of Origin and Evolution of Planetary Systems’’ of the Ministry of Education, Culture, Sports and Technology (MEXT).
References [1] Giese RH, Weiss K, Zerull RH, Ono T. Large fluffy particles—a possible explanation of the optical properties of interplanetary dust. A&A 1978;65:265–72. [2] Zerull RH, Giese RH, Weiss K. Scattering functions of nonspherical dielectric and absorbing particles VS Mie theory (E). Appl Opt 1977;16:777–8. [3] Weiss-Wrana K. Optical properties of interplanetary dust—Comparison with light scattering by larger meteoritic and terrestrial grains. A&A 1983;126:240–50. [4] Gustafson BA˚S. Microwave analog to light scattering measurements: a modern implementation of a proven method to achieve precise control. JQSRT 1996;55:663–72. [5] Gustafson BA˚S. Microwave analog to light scattering measurements. In: Mishchenko MI, Hovenier JW, Travis LD., editors. Light scattering by nonspherical particles: theory, measurements, and geophysical applications. New York: Academic Press; 1999. p. 367–90 [Chapter 13], 1999. ISBN 0-12-498660-9. [6] Meakin P. Effects of cluster trajectories on cluster–cluster aggregation: a comparison of linear and Brownian trajectories in two- and three-dimensional simulations. Phys Rev 1984;A29:997–9. [7] Mukai T, Ishimoto H, Kozasa T, Blum J, Greenberg JM. Radiation pressure forces of fluffy porous grains. A&A 1992;262:315–20. [8] Mann I, Okamoto H, Mukai T, Kimura H, Kitada Y. Fractal aggregate analogues for near solar dust properties. A&A 1994;291:1011–8. [9] Manickavasagam S, Menguc MP. Scattering matrix elements of fractallike soot agglomerates. Appl. Opt. 1997;36:1337–51. [10] Mackowski DW. A simplified model to predict the effects of aggregation on the absorption properties of soot particles. JQSRT 2006;100:237–49. [11] Wurm G, Blum J. Experiments on preplanetary dust aggregation. Icarus 1998;132:125–36. [12] Blum J, Wurm G. Experiments on sticking, restructuring, and fragmentation of preplanetary dust aggregates. Icarus 2000;143:138–45. [13] Mann I, Kimura H, Kolokolova L. A comprehensive model to describe light scattering properties of cometary dust. JQSRT 2004;89:291–301. [14] Draine BT. The discrete dipole approximation and its application to interstellar graphite grains. ApJ 1988;333:848–72.
ARTICLE IN PRESS 80
Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
[15] Draine BT, Flatau PJ. User guide to the discrete dipole approximation code DDSCAT 6.1, 2004. hhttp://arxiv.org/abs/astro-ph/ 0409262i. [16] Okamoto H. Light scattering by clusters: the A1-term method. Opt Rev 1995;2:407–12. [17] Okamoto H, Xu Y. Light scattering by irregular interplanetary dust particles. Earth Planets Space 1998;50:577–85. [18] Mackowski DW, Mishchenko MI. Calculation of the T matrix and the scattering matrix for ensembles of spheres. J Opt Soc Am A 1996;13:2266–78. [19] Mishchenko MI, Mackowski DW. Electromagnetic scattering by randomly oriented bispheres: comparison of theory and experiment and benchmark calculations. JQSRT 1996;55:683–94. [20] Sun W, Loeb NG, Fu Q. Finite-difference time-domain solution of light scattering and absorption by particles in an absorbing medium. Appl Opt 2002;41:5728–43. [21] Sun W, Loeb NG, Videen G, Fu Q. Examination of surface roughness on light scattering by long ice columns by use of a twodimensional finite-difference time-domain algorithm. Appl Opt 2004;43:1957–64. [22] Yang P, Liou KN. Finite difference time domain method for light scattering by nonspherical and inhomogeneous particles. In: Mishchenko MI, Hovenier JW, Travis JD, editors. Light scattering by nonspherical particles. San Diego: Academic Press; 2000. p. 173–221. [23] Chylek P, Srivastava V. Dielectric constant of a composite inhomogeneous medium. Phys Rev B 1983;27:5098–106. [24] Macke A. Scattering of light by polyhedral ice crystals. Appl Opt 1993;32:2780–8. [25] Muinonen K, Nousiainen T, Fast P, Lumme K, Peltoniemi JI. Light scattering by Gaussian random particles: ray optics approximation. JQSRT 1996;55:577–601. [26] Okada Y, Nakamura AM, Mukai T. Light scattering by particulate media of irregularly shaped particles: laboratory measurements and numerical simulations. JQSRT 2006;100:295–304. [27] Mukai T, Okada Y. Optical properties of large aggregates. In: Kru¨ger H, Graps A., editors. Dust in planetary systems, Workshop, September 26-30, 2005, Kauai, Hawaii, ESA Publications, SP-643; January 2007. [28] Andersen AC, Mutschke H, Posch Th, Min M, Tamanai A. Infrared extinction by homogeneous particle aggregates of SiC, FeO and SiO2 : Comparison of different theoretical approaches. JQSRT 2006;100:4–15. [29] Koehler M, Kimura H, Mann I. Applicability of the discrete-dipole approximation to light-scattering simulations of large cosmic dust aggregates. A&A 2006;448:395–9. [30] Taflove A, Hagness SC. Computational electrodynamics: The finite-difference time-domain method (Artech House Antennas and Propagation Library). MA: Artech House; 2005. 1038p. [31] Kozasa T, Blum J, Mukai T. Optical properties of dust aggregates, I Wavelength dependence. A&A 1992;263:423–32. [32] Wolff MJ, Clayton GC, Martin PG, Schulte-Ladbeck RE. Modeling composite and fluffy grains: the effects of porosity. ApJ 1994;423:412–25. [33] Wolff MJ, Clayton GC, Gibson SJ. Modeling composite and fluffy grains II Porosity and phase functions. ApJ 1998;503:815–30. [34] Voshchinnikov NV, Mathis JS. Calculating cross sections of composite interstellar grains. ApJ 1999;526:257–64. [35] Voshchinnikov NV, Il’in VB, Henning Th. Modelling the optical properties of composite and porous interstellar grains. A&A 2005;429:371–81. [36] Voshchinnikov NV, Il’in VB, Henning Th, Dubkova DN. Dust extinction and absorption: the challenge of porous grains. A&A 2005;445:167–77. [37] Koehler M, Mann I. Light-scattering models applied to circumstellar dust properties. JQSRT 2004;89:453–60. [38] Liou KN, Hansen JE. Intensity and polarization for single scattering by polydisperse spheres: a comparison of ray optics and Mie theory 1971;28:995–1004. [39] Mackowski DW. Calculation of total cross sections of multiple-sphere clusters. J Opt Soc Am A 1994;11:2851–61. [40] McQueen J. Some methods for classification and analysis of multivariate observations. In: Proceedings of the fifth Berkeley symposium on mathematical statistics and probability, 1967. p. 281–97. [41] Mukai T. Cometary dust and interplanetary particles. In: Bonetti A, Greenberg JM, Aiello SA, editors. Evolution of interstellar dust and related topics. Amsterdam: Elsevier; 1989. p. 397–445. [42] Dorschner J, Begemann B, Henning T, Jaeger C, Mutschke H. Steps toward interstellar silicate mineralogy. II Study of Mg–Fesilicate glasses of variable composition. A&A 1995;300:503–19. [43] Preibisch T, Ossenkopf V, Yorke HW, Henning Th. The influence of ice-coated grains on protostellar spectra. A&A 1993;279:577–88. [44] Bohren CF, Huffman DR. Absorption and scattering of light by small particles. New York: Wiley; 1983 530p. [45] Hoekstra AG, Frijlink M, Waters LBFM, Sloot PMA. Radiation forces in the discrete-dipole approximation. J Opt Soc Am A 2001;18:1944–53. [46] Xu Y-l. Electromagnetic scattering by an aggregate of spheres. Appl Opt 1995;34:4573–88. [47] Xu Y-l. Electromagnetic scattering by an aggregate of spheres: far field. Appl Opt 1997;36:9496–508.
Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80 www.elsevier.com/locate/jqsrt
Grouping and adding method for calculating light scattering by large fluffy aggregates Y. Okadaa,, T. Mukaia, I. Manna, H. Nomurab, T. Takeuchia, I. Sanoc, S. Mukaic a
Graduate School of Science and Technology, Kobe University, 1-1 Rokkodai, Nada, Kobe 657-8501, Japan b School of Mathematics and Physics, Queen’s University Belfast, Belfast BT 1NN, Northern Ireland, UK c Faculty of Science and Technology, Kinki University, 3-4-1 Kowakae, Higashi-Osaka 577-8502, Japan Received 1 August 2006; received in revised form 25 February 2007; accepted 5 April 2007
Abstract We present a method to derive the light scattering properties of very porous fractal aggregates composed of a large number of monomers where the size parameter of monomer is larger than unity. Our new method is based on the grouping of the aggregate: The aggregate is divided into groups, where each group is located along a line of the incident light, and the scattering properties of the group are calculated taking into account multiple scattering with monomers located inside the group, as well as those in a buffer region around the group. The scattering and absorption efficiencies are obtained by adding the resultant scattering properties for all the groups. We have shown that the method effectively works when the monomer scatters the incident lights predominantly in the forward direction, which is the case if the monomer size is large compared to the wavelength of the incident light. The errors in resulting scattering and absorption efficiencies for porous aggregates are investigated for various refractive indices and sizes of monomers. We found that the errors are larger for low absorbing materials and they can be reduced by expanding the buffer region. In the case of the buffer region for each group consisting of 1/8 of the total number of monomers, the results show errors less than 15% and 10% for absorption and scattering, respectively. It is also shown that the errors have a small standard deviation (i.e., 2%) for different directions of the incident light. r 2007 Elsevier Ltd. All rights reserved. Keywords: BCCA; Cluster T-matrix method; Dominant forward scattering; Scattering and absorption efficiencies; Non-Rayleigh region
1. Introduction Porous irregular particles are ubiquitous in space and in the atmospheres of the Earth and other planets. Since the 1970s it has been obvious that data on observations of dust in the solar system cannot be reproduced by applying the light scattering properties of compact spherical particles with the same volume [1]. For estimating the scattering properties of the irregular particles laboratory measurements of irregular particles were taken into account [2–5]. Also, the numerical studies of the formation of aggregates by sticking of particles are studied [6]. Corresponding author. Tel.: +81 78 803 6566; fax: +81 78 803 6483.
E-mail address: [email protected] (Y. Okada). 0022-4073/$ - see front matter r 2007 Elsevier Ltd. All rights reserved. doi:10.1016/j.jqsrt.2007.04.006
ARTICLE IN PRESS 66
Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
The concept of fractal aggregates has been applied to the calculations of the radiation pressure [7] and the temperature [8] of the interplanetary dust particles. For atmospheric models, biomass burning aerosol particles produced from the forest fires are often modeled as fractal aggregates [9,10]. From laboratory experiments of dust aggregation in microgravity environment it is shown that aggregates form from the sticking of particles [11,12]. In many cases, however, the concept of fractal aggregates has been applied to particles that are not produced by such a growth process. The fractal aggregates provide a concept for a reproducible description of irregular particles, which are actually produced by different processes [13]. Since light scattering properties are fundamental for various applications to interpret data observed by such as Earth-observing satellites, space-borne sensors and ground-based instruments observing dust from comets, asteroids and protoplanetary disks, it is important to have a through understanding of the light scattering of fractal aggregates. Many different methods have been proposed for the calculation of light scattering properties of nonspherical particles, such as DDA (discrete dipole approximation) [14,15], a1 term method [16,17], T-matrix for cluster of spheres [18,19], FDTD (finite difference time domain) method [20–22], EMT-Mie (effective medium theories combined with Mie theory) method [7,23], GO (geometrical optics) [24–27]. Calculations are especially difficult for aggregates composed of a large number of spherical monomers when the monomer size is larger than the wavelength of the incident light. The DDA method is based on the calculation of interactions of dipoles which describe the shape of nonspherical particles. To calculate for fractal aggregates composed of monomers, the shape of the aggregate should be described with a number of dipoles. In order to keep the accuracy, the product of dipole separation d, wavenumber k and refractive index jmj (i.e., jmjkd) should be less than 1.0 [14,15]. This requires the subdivision of each monomer into number of dipoles when the monomer size is larger than the wavelength of the incident light as conducted in [28,29] and the memory limitation leads to the calculation of aggregates composed of no more than a few hundred monomers. The FDTD method solves the Maxwell equations by calculating the time evolution of electromagnetic waves (EMs), and the shape of the particles is described using cubic cells. The resultant scattering properties in the time domain are then converted into the frequency domain. The FDTD method can be used to calculate the scattering properties of particles larger than that are computable with DDA [20–22,30]. However for the considered aggregates, the geometric limitations of the FDTD method mean that the shape should be finely discretized into a larger number of cubic cells in order to reproduce the spherical surface of each monomer [30]. This also requires a large amount of computer memory, and hence it is difficult to handle aggregates composed of a large number of spheres with the FDTD method. The EMT-Mie method is easy to implement since it first derives the effective refractive index based on the volume fraction of the material relative to empty space, and then calculates the light scattering properties based on Mie theory [7,23,31]. However, compared with other methods such as DDA, the accuracy of the EMT-Mie method deteriorates very rapidly as the monomers become larger than the wavelength of the incident light (i.e., non-Rayleigh region), especially for highly porous particles [7,31–37]. In contrast to most methods, which have upper limitations for the monomer size, the GO method has a lower limitation for the monomer size. The GO method is applicable only when the circumference of the monomer is much larger than the wavelength of the incident light (i.e., around 100 times) and it is based on treating the incoming light as a bundle of rays. Fresnel’s formulae are used to describe the reflection and refraction of rays [24–26]. However, there is still a large gap in the range of allowable sizes of monomers for porous particles between other methods (around a few times the wavelength) and the GO method (around 100 times) [38], indicating the need to modify the GO method to include the influence of wave optics or the need to expand the size range into larger one for other methods some how. The cluster T-matrix method developed by Mackowski and Mishchenko [18] is based on the superposition of the T-matrix for each sphere. The calculation is accurate as long as higher order of expansions of scattering function are calculated. However, this method is limited by the extensive computing resources required for a larger number of monomers and for larger monomers. The upper limitation with our computing resources (shown in Section 3) is around 600 monomers for the circumference of the monomer of 4.5 times the wavelength.
ARTICLE IN PRESS Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
67
As shown above, it is required to make a convenient method to derive the light scattering properties for the porous fluffy particles composed of thousands (or a larger number) of monomers where the circumference of the monomer is larger than the wavelength of the incident light (i.e., non-Rayleigh region) but smaller than the region of validity for the GO method. In this paper we try to expand the number of monomers in around 10 times with our proposed method applied to the currently available cluster T-matrix method of fixed direction of incident light [18]. In Section 2, scattering characteristics of single compact spheres are presented, then our proposed method is introduced based on the scattering characteristics of each monomer. In Section 3, errors are investigated for different sizes of buffer region where the higher order of multiple scattering is considered. Then, error analysis is conducted as a function of real and imaginary part of refractive indices and sizes of monomers. The error maps of absorption are compared with the maps of scattering properties of a single monomer. Also the errors for different directions of incident lights are investigated. In Section 4, examples of the light scattering properties of very porous fractal aggregates composed of up to 8192 monomers are presented. The result shows the power law relationship between the number of monomers and the efficiencies. In the conclusion, we summarize this study and also we briefly mention some possible future improvements combining with methods other than the cluster T-matrix method used in this study. 2. Method 2.1. Forward scattering of a single compact sphere For a model of the shape and structure of the scattering objects, we consider the BCCA (ballistic cluster–cluster aggregates) as shown in Fig. 1, where the constituent monomers have same radii of rm . The fractal dimension Df determined from the radius of gyration and the number of composing monomers is around 2.0 [7], which shows the high fluffiness of the aggregates. BCCAs have very large porosities (i.e., 99% of vacuum inside the characteristic radius) [7,31]. It is shown that EMT-Mie method yields a few tens of percentages of errors for the light scattering properties of very porous particles such as BCCA, when the monomer size becomes larger than the wavelength of the incident light (i.e., non-Rayleigh region), and the errors tend to be larger for larger particles (e.g., see Fig. 6 in [32], Fig. 3 in [33] and Fig. 9 in [35]). Table 1 summarizes the parameters of the light scattering simulations considered in this study. The size information is described with size parameter x ¼ 2pr/l, where r is the radius of the particle and l is the
Fig. 1. Three-dimensional image of a BCCA composed of 512 monomers viewed from a zenith angle of 151. Note that all the monomers are of the same size and chemical composition.
ARTICLE IN PRESS 68
Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
Table 1 Parameters of the light scattering simulations xm l N mono m ¼ n þ ik Qabs Qsca
Monomer size parameter: Wavelength of the incident light: Number of monomers in the aggregate: Refractive index of all the monomers: Absorption efficiency of the aggregate: Scattering efficiency of the aggregate: The Qabs and Qsca are defined for the volume-equivalent sphere radius.
p11 (Θ) / p11 (0°)
100 xm=0.001 1.0 5.0 10.0 50.0 100.0 500.0
10-2
10-4
10-6 0
30
60 90 120 Scattering angle Θ (degree)
150
180
Fig. 2. Scattering function p11 ðYÞ of a single compact sphere as a function of size parameter of the sphere, xm , where m ¼ 1:5+0:1i. Each scattering function is normalized to its value at Y ¼ 0 .
wavelength of the incident light. The radius of the aggregate is often taken as reff , corresponding to the radius of a sphere with the volume equal to that of the aggregate. Here we follow this custom when to the pffiffiffiffireferring ffi radius of the aggregate. We mainly consider the size parameter of the monomer xm ¼ xag / 3 N mono where xag represents the size parameter of the aggregate calculated from reff and l. N mono is the number of monomers. The chemical composition of the monomers composing the aggregate is considered to be homogeneous, and we study the light scattering as a function of refractive index m ¼ n þ ik. Of the various light scattering properties of the aggregates, we mainly consider the efficiency factors for absorption Qabs and for scattering Qsca . The angular distribution of the scattered light is described with a scattering function p11 ðYÞ, which denotes the intensities of scattered light as a function of scattering angle Y, where Y ¼ 0 represents the forward and Y ¼180 the backward direction relative to the direction of the incident light. Fig. 2 shows the scattering functions of single compact spheres normalized to the intensity at Y ¼ 0 for various particle sizes calculated using Mie theory. This figure clearly shows that the forward scattering component (e.g., Yo30 ) becomes more dominant relative to scattering into side direction (e.g., Y460 ) as the particle becomes larger. This characteristic forms the basis of the method proposed in Section 2.2.
2.2. Grouping and adding method The maximum computable number of monomers for the given size of monomers is limited in the cluster Tmatrix method (and in most other methods). As the number of monomers increases, the memory requirement increases very rapidly. At the same time, a larger number of equations have to be solved for the interactions of electric fields and the calculation can become unstable.
ARTICLE IN PRESS Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
69
Our method proposed in this section, called the grouping and adding method (GAM), has the potential to calculate the scattering properties of large fluffy aggregates beyond the computable limitations for the number of monomers, and involves the following steps: (1) Divide the aggregate into N DIV groups. (2) Calculate efficiencies (i.e., Qabs and Qsca ) of each group. (3) Sum the efficiencies of all the groups. The complete aggregate is first divided into groups, which consist of monomers close to one another on the two-dimensional plane (standard plane) perpendicular to the direction of the incident light as shown in Fig. 3a and a0 , where all constituent monomers are projected on this standard plane. In this study, the grouping is conducted using the K-means method [40] as follows: (1) Randomly choose N DIV monomers on the standard plane and each of which is used as a seed for each group. (2) Select other members in each group based on the shortest distance from the seed monomers for the group on the standard plane.
Top view 9 group 4 surrouding monomers
8 7
1
Y
Y
10
4 6 2 5
3
EM wave
EM wave X
X
Z
Z
Side view
EM wave X
EM wave X
Fig. 3. Schematics of grouping of the aggregate (a and a’) and of the buffer region (b and b’) for EMs propagating in the Z direction. The aggregates shown in Fig. 1 are grouped with N DIV ¼ 10. In panels b and b’, gray- and black-shaded circles are the N in gr monomers in the group and N out gr monomers in the buffer region, respectively, and the large solid circle and dashed one denote the maximum radius of group 4 and buffer region for the group, respectively. The direction of the incident light is shown in the right bottom of each panel.
ARTICLE IN PRESS Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
70
(3) Calculate the center of each group using the selected members together with the seed monomer for each group. (4) Regroup all the monomers based on the shortest distance from the centers of groups. (5) Repeat steps 3 and 4 until the changes in the locations of the centers of groups for all the groups become smaller than a pre-configured threshold. The result of the grouping method is not reproducible because the K-means method uses a random number generator to locate the first seed monomers. However, we have found that the different groupings generated for the same N DIV do not largely influence on the results. In the calculation of the light scattering properties of each group, a contribution from the monomers belonging to other groups cannot be neglected. To examine such a contribution, we use the N in gr monomers of the group and the N out monomers in the buffer region (see Fig. 3b and b’, and Table 2). The buffer region is set gr as that surrounding the selected group based on x and y coordinates on the standard plane, and not using information of z coordinates (i.e., the direction of the incident light). From the inner part of the buffer region, N out gr monomers are selected until the chosen number (e.g., N mono =8) are selected. When grouping is not used, the light scattering properties of the aggregate are calculated in the presence of all the monomers. In contrast, in our proposed method, only N out gr nearby surrounding monomers in the buffer region are used together with the N in gr monomers to calculate the scattering contributions of these monomers. in The N out gr monomers are used in order to partly account for the light that is scattered sideways by the N gr monomers. The grouping reduces the memory requirement for the calculation, from being proportional to out N mono to being proportional to (N in gr + N gr ). The light scattering properties of each group are calculated by considering a contribution from the monomers included in the group and its buffer region, as follows: out (1) The monomers classified into the group (N in gr ) and in the buffer region (N gr ) are used for the light scattering calculations for a particular group. (2) Only the light scattering properties of the N in gr monomers are stored as the result for the selected group. (3) Steps 1 and 2 are repeated for all the groups. (4) The scattering properties from all the groups are simply added to deduce the total light scattering properties of the aggregate.
In the Superposition T-matrix method for clusters of spheres, cross sections (i.e., C ext and C abs ) are formulated as the superposition of those of composing spheres as the following: C ext ¼
Ns X
C iext ,
(1)
C iabs ,
(2)
i¼1
C abs ¼
Ns X i¼1
respectively from Eqs. (24) and (43) of [39]. Where C ext , C abs are the cross sections of extinction and absorption of total cluster while C iext and C iabs are those of ith composing sphere. Ns is the number of particles. In the above equations, multiple interactions of the EM wave between the particles are included into C iext and C iabs . Table 2 Parameters used in the GAM Number of groups: Number of monomers in each group: Number of monomers in the buffer region:
N DIV N in gr N out gr
ARTICLE IN PRESS Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
71
In the Fortran code developed by Mackowski et al., C ext and C abs are calculated in two ways (see ftp:// ftp.eng.auburn.edu/pub/dmckwski/scatcodes/index.html): (1) summation of the individual sphere results, (2) analytical integration of the cluster scattered field expansion. Then the comparison of calculations (1) and (2) have been checked to agree if the interaction equations are solved to converge. In our paper, we use (1) summation of the individual sphere results. Currently only Qabs and Qsca are calculated based on our method. The phase matrix are not just additive from those of each composing particle. Therefore in this paper, we will not go further to use our method for phase matrix. The cluster T-matrix code outputs the efficiency factors for extinction and absorption of the TM (transverse magnetic) and TE (transverse electric) mode of the EMs for each monomer. These efficiencies are influenced by (a) the incident electric field and (b) the electric fields induced by monomers other than the selected one, due to the multiple scattering by monomers. For each monomer, we average the efficiencies for the TM and TE modes to calculate efficiencies for non-polarized incident light. The grouping increases the computing time N DIV fold. However the grouping also works for the reduction out of time by the reduction of monomers from N mono to (N in gr + N gr ) for each group. Moreover, once the grouping is conducted, the calculation for each group can be parallelly performed in different computers, resulting in the faster calculations. 3. Error analysis The GAM is based on disregarding the scattering by monomers outside the combined area of the group plus its buffer region. Neglecting the contribution from outer monomers in this way will introduce errors into the scattering properties (i.e., Qsca and Qabs ). To ensure the efficacy of the GAM, we need to investigate the associated errors caused under various conditions (i.e., various ðn; kÞ and xm ). Although the GAM would be applied for a large fluffy aggregate, where the total number of monomers is beyond the limit of available computer system, we need to conduct the investigation of its applicability under the computable number of monomers. In this section, the results of the GAM for small aggregates composed of 512 monomers are compared to those calculated without grouping. The errors in Qabs and Qsca are defined as Errorð%Þ ¼
jQGAM Qorg j 100, Qorg
(3)
where QGAM denotes the value based on the GAM while Qorg is the value without grouping. The 600 monomers for xm ¼ 4:5 is the upper limit for our computing resources (i.e., Pentium 4 2.8 GHz CPU with 2 GB of RAM). Therefore, we investigate the errors under this upper limitation. 3.1. Errors for different sizes of buffer region for silicate, pyroxene and amorphous carbon Fig. 4 shows the errors in Qabs and Qsca for different sizes of buffer region in the cases of five different shapes of aggregates. Only one direction of incident light is considered. We investigate errors for silicate, pyroxene and amorphous carbon, which exhibit low, moderate and high absorbing materials, respectively. Although pyroxene belongs to the silicate group, as a matter of practical convenience, we differentiate them as ‘‘silicate’’ and ‘‘pyroxene’’ in this paper. We use the refractive indices at a wavelength of 0:55 mm for these materials (see Table 3). For the pyroxene, we use the data for Mgx Fe1x SiO3 , where x ¼ 0:7 [42]. For the silicate and pyroxene (see Fig. 4a–d), the errors in Qabs decrease with increasing buffer region (i.e., N out gr to N mono ratio), whereas the errors in Qsca for the two materials do not show clear trends. However, the errors in Qsca are already small (i.e., less than 10%) for most of the results. The absorption occurs when the light interacts with each monomer, and hence the errors are directly influenced by the number of side monomers outside the buffer region that are disregarded.
ARTICLE IN PRESS Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
20
20
10
0
0 1/2
1
0 1/16 1/8
5
1/4
5
Nout gr /Nmono
Nout gr /Nmono 20 shape 1 2 3 4 5
10
0 1
1/2
1/2
0
0 1/16 1/8
5
1/4
5
1/4
Err (Qsca) (%)
15
10
0 1/16 1/8
Nout gr /Nmono
Nout gr /Nmono
20 shape 1 2 3 4 5
15
10
Nout gr /Nmono
1
0 1/2
0 1/4
5
0 1/16 1/8
5
1/4
10
0 1/16 1/8
Err (Qsca) (%)
15
shape 1 2 3 4 5
1
20
1/2
Err (Qabs) (%)
15
shape 1 2 3 4 5
1
20
Err(Qabs) (%)
1/2
Err (Qsca) (%)
15
10
0 1/16 1/8
Err (Qabs) (%)
15
shape 1 2 3 4 5
1
shape 1 2 3 4 5
1/4
72
Nout gr /Nmono
Fig. 4. Errors in Qabs and Qsca calculations as functions of the number of monomers (N out gr ) included in the buffer region relative to the total number of monomers (N mono ). Results are for silicate (a and b), pyroxene (c and d), and amorphous carbon (e and f), where N mono ¼ 512, xm ¼ 3:0 and N DIV ¼ 60.
ARTICLE IN PRESS Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
73
Table 3 Refractive indices of the materials at a wavelength of 0:55 mm Silicate Pyroxene Amorphous carbon
m ¼ 1:48 þ 0:000028i [41] m ¼ 1:635 þ 0:0042i [42] m ¼ 2:03 þ 0:77i [43]
In contrast, Qsca corresponds to values integrated for solid angles, and so the errors depend on the directional distribution of the scattered energy. For example, when the light is scattered more in the forward direction, side monomers contribute less to the scattering efficiencies. Therefore, once sufficient monomers are considered in the buffer region, the errors in Qsca are little affected by the size of buffer region, as shown in Fig. 4b and d. Fig. 4a–d shows that the errors are higher in Qabs than in Qsca when the buffer region is small. Therefore, the N out gr -to-N mono ratio can be set based mainly on the errors in Qabs for materials with low to moderate absorption. The errors in both Qabs and Qsca for highly absorbing amorphous carbon (see Fig. 4e and f) are less than 10% for the sizes of buffer region investigated. The errors in Qabs does not decrease when the size of the buffer region is increased. This is because the lights interacting with amorphous carbon are highly absorbed by each monomer. This results in multiple scattering with higher energies being confined to a smaller region compared to those of silicate and pyroxene. It is considered that the buffer regions investigated include the large region where most multiple scatterings with higher energies for amorphous carbon take place. Fig. 4 indicates that the N out gr -to-N mono ratio needs to be set as at least 1/8 in order to keep the errors in Qabs and Qsca mostly below 10% for xm ¼ 3:0. This allows us to treat 8 times 600 monomers with errors of less than 10% for these three materials. 3.2. Error maps for ðn; kÞ and xm Fig. 4 shows the results of the errors of different sizes of buffer regions for only three materials. In order to systematically investigate the errors caused by the grouping, we have investigated the errors for various sets of ðn; kÞ and xm . Fig. 5 shows error maps for n vs k for different xm . The N out gr -to-N mono ratio for the buffer region is set as 1/8 for this figure. Only one direction of the incident light is considered. The variations in errors for different directions of the incident light are discussed in Section 3.4. For the results of xm ¼ 1:0, errors are around 5% for all sets of n and k. The errors of xm ¼ 1:0 are smaller than those of xm ¼ 3:0 and xm ¼ 4:5. The reason of this smaller errors for xm ¼ 1:0 will be shown in Section 3.3. For xm ¼ 3:0 and xm ¼ 4:5, the errors in Qabs show certain trends for different ðn; kÞ and xm . For smaller n (i.e., np1:30), the errors in Qabs is around 5% for all sets of xm and k. The real part of refractive index (i.e., n) denotes how much the light is scattered. The smaller n results in less scattering of the light and that corresponds to the smaller errors in Qabs . On the other hand, for nX1:4, errors in Qabs become larger for smaller k (i.e., ko0:01). This is because a smaller k indicates less absorption by each monomer, with the resulting higher order of multiple scattering into the side direction outside the buffer region causing larger errors. The error maps of Fig. 5 show that the grouping with an N out gr -to-N mono ratio of 1/8 causes errors in Qabs of less than 15% for ko0:01 and nX1:4, and less than 10% for the remaining parameter combinations. The errors in Qsca for all sets of n and k in Fig. 5 are less than 10%. This indicates that the use of the GAM allows us to expand the number of monomers eight-fold within these errors for all sets of n and k for xm X1:0. 3.3. Errors related to the scattering properties of each monomer In Fig. 5, the errors in Qabs show certain trends (i.e., higher errors for larger n and for smaller k). These trends in errors are considered to have relations with how each monomer scatterers the light.
ARTICLE IN PRESS 74
Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
Errors in Qabs (%)
Errors in Qsca (%)
xm=1.0
0.00001 0.00001 0.00005 0.00005 0.0001 0.0001 0.0005 0.0005 0.001 0.001 1.1 1.1 0.005 0.005 1.2 1.2 0.01 0.01 k k 1.3 1.3 1.4 1.4 0.05 0.05 1.5 1.5 0.1 0.1 1.6 1.6 1.7 1.7 0.5 0.5 n n 1.8 1.8 1.9 1.0 1.9 1.0 2.0 2.0
xm=3.0
0.00001 0.00001 0.00005 0.00005 0.0001 0.0001 0.0005 0.0005 0.001 0.001 1.1 1.1 0.005 0.005 1.2 1.2 0.01 k 1.3 0.01 k 1.3 1.4 1.4 0.05 0.05 1.5 1.5 0.1 1.6 0.1 1.6 n 1.7 n 1.7 0.5 0.5 1.8 1.8 1.9 1.0 1.9 1.0 2.0 2.0
xm=4.5
0.00001 0.00001 0.00005 0.00005 0.0001 0.0001 0.0005 0.0005 0.001 0.001 1.1 1.1 0.005 0.005 1.2 1.2 0.01 k 1.3 1.3 0.01 k 1.4 1.4 0.05 0.05 1.5 1.5 0.1 1.6 1.6 0.1 n 1.7 1.7 n 0.5 0.5 1.8 1.8 1.9 1.0 1.9 1.0 2.0 2.0 Errors (%)
0
5
10
15
20
Fig. 5. Maps of errors in Qabs and Qsca for n vs k for different xm . The buffer region has an N out gr -to-N mono ratio of 1/8, and N mono ¼ 512 and N DIV ¼ 60.
ARTICLE IN PRESS Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
xm=1.0
1.1
1.2
1.3
1.4
1.5
n
1.6
1.7
1.8
1.9
2.0
75
xm=3.0
0.00001 0.00005 0.0001 0.0005 0.001 0.005 0.01 k 0.05 0.1 0.5 1.0
1.1
1.2
1.3
1.4
1.5
n
1.6
1.7
1.8
1.9
2.0
0.00001 0.00005 0.0001 0.0005 0.001 0.005 0.01 k 0.05 0.1 0.5 1.0
xm=4.5
1.1
1.2
1.3
1.4
1.5
n
1.6
1.7
1.8
1.9
2.0
0.00001 0.00005 0.0001 0.0005 0.001 0.005 0.01 k 0.05 0.1 0.5 1.0
mono Qmono sca /g
0
2
4
6
8
10
Fig. 6. Scattering properties of a monomer for n vs k for different xm . Qmono and gmono , which are the scattering efficiency and asymmetry sca parameters of the monomer, respectively, are calculated using Mie theory.
Fig. 6 presents the scattering efficiencies of the single monomer (Qmono sca ) divided by the asymmetry parameter of the monomer (gmono ), based on Mie theory. The asymmetry parameter g is defined as [44] Z p11 ðyÞ cos y dO, (4) g ¼ o cos y4 ¼ 4p
where p11 ðyÞ is the scattering function for the scattering angle y. O is the solid angle for the integration. The value of g approaches 1 as forward scattering becomes dominant, and hence the 1=g represents how much light is scattered sideways. Qmono represents how effectively the incident light is scattered, sca mono and so the product Qmono indicates how widely the light is distributed after scattering by each sca =g monomer. Therefore, the trends seen in Fig. 6 coincide with those in Fig. 5; namely, larger errors in Qabs by the GAM are caused by a lower degree of forward scattering by the monomer. This indicator can therefore be used to quantify which combination of n, k, and xm results in larger errors. For xm ¼ 1:0, mono the Qmono shows smaller values without trends for n vs k, corresponding to the smaller errors shown in sca =g Fig. 5.
ARTICLE IN PRESS 76
Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
20
20
15
mean : 6.71(%), SD : 2.44(%) Err (Qsca) (%)
Err (Qabs) (%)
mean : 9.09(%), SD : 2.09(%)
10 5 0
15 10 5 0
0
100 200 ID of direction
300
0
100 200 ID of direction
300
Fig. 7. Errors in Qabs and Qsca of a BCCA for 324 directions of incident light, where xm ¼ 3:0, N mono ¼ 512, and N DIV ¼ 60 for silicate (i.e., m ¼ 1:48 þ 0:000028i). The N out gr -to-N mono ratio is set at 1/8.
3.4. Errors for the different directions of incident light The error analysis presented above considered incident light from only one direction for each shape of aggregate. The errors differ for different directions of incident light, because this causes differences in the alignment of monomers. Therefore, we also investigate variations in errors for different directions of incident light. The BCCA consisting of 512 monomers is rotated in azimuth and in polar angles at regular intervals to give a total of 324 directions. Fig. 7 shows the errors in Qabs and Qsca for 324 directions of incident light for xm ¼ 3:0. The average errors in Qabs and Qsca are less than 10% with a standard deviation of 2%. In general, the values of Qabs and Qsca for aggregates are deduced from averaging over those values obtained at different directions [8]. 3.5. Selection of N DIV The division of the aggregate into groups requires the light scattering of N DIV groups to be computed. We recommend that the N out gr -to-N mono ratio is set as 1/8, especially for materials with a low absorption. Then in N DIV can be taken so that (computable number of monomers—N out gr ) is assigned to N gr . In this paper, we have investigated errors for N DIV ¼ 60. The errors for smaller N DIV are either larger or smaller due to the following two reasons: out (1) The maximum radius of the N in gr monomers of each group increases as N DIV decreases. In this case the N gr monomers are assigned from monomers located over a larger radius with the thinner width, leading to multiple scattering between outer of N in gr monomers and those of buffer region being treated less. This results in the increase in the errors. (2) The increase of N in gr monomers for smaller N DIV causes higher-order multiple scattering in the group monomers being considered. This results in the decrease in errors.
Fig. 8 shows the errors for different numbers of groups (i.e., N DIV ¼ 72, 30, and 16) where the N out gr -toN mono ratio is fixed at 1/8. However, the errors in Qabs do not show clear trends for the increase of N DIV because of the complicated combination of the item 1 and 2 of the above. Then the N DIV can be selected based on the computing time. 4. Application to large aggregates composed of a number of monomers beyond the limit of computing resources The results of error analysis above indicates that the GAM can be utilized to increase the number of monomers beyond the computable number for single aggregate within certain errors. In this section, we
ARTICLE IN PRESS Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
Error in Qsca (%)
Error in Qabs (%) 15 10 5
15 NDIV=72 30 16
10
n=1.1
n=1.6
n=1.2
n=1.7
0 15 n=1.2
10
5
5
n=1.7
0 15
0 15 10
n=1.1
5
n=1.6
0 15 10
n=1.3
n=1.3
10
5
5
n=1.8
n=1.8 0 15
0 15 n=1.4
n=1.4
10
10 5
5
n=1.9
n=1.9
0 15
0 15
10
10
n=1.5
5 0 10-5
77
n=1.5 10-3
5
n=2.0 10-1 1 10-5 k
10-3
n=2.0
10-1 1
0
10-5
10-3
10-1 1 10-5 k
10-3
10-1 1
Fig. 8. Errors in Qabs and Qsca for different sets of n and k, where N mono ¼ 512 and xm ¼ 3:0. Different number of grouping (i.e., N DIV ¼ 72, 30 and 16) are investigated. The N in gr -to-N mono ratio is set as 1/8.
present the light scattering properties of large aggregates composed of up to thousands of monomers calculated with the GAM. Fig. 9 shows the Qabs and Qsca of very porous BCCAs composed of up to N mono ¼ 4096 for silicate and pyroxene, and those up to 8192 for amorphous carbon where xm is 4.5. The results of N mono p512 are calculated without grouping while N mono X1024 based on the GAM. The results of up to N mono ¼ 4096 are calculated with the N out gr -to-N mono ratio of 1/8. The results for the amorphous carbon for N mono ¼ 8192 are in out calculated with an N out gr -to-N mono ratio of 1/16 due to the upper limitation of (N gr + N gr ) being 600 for our computing resources. In this case of N mono ¼ 8192, the errors in Qabs and Qsca are expected to be around 10% for amorphous carbon (see Fig. 4e and f). In Fig. 9, only one direction of incident light is considered for each of the aggregates with different N mono . Therefore, the deviation of Qsca from others seen for N mono ¼ 128 and 256 are attributed to this rather than being considered representative of random orientation results. However, all of the results show the characteristic of increasing efficiencies for increasing N mono . We note that in contrast with the case of compact sphere, Qabs plus Qsca in Fig. 9 does not converge to 2 when the aggregate size increases. It is because the efficient factors of BCCA are defined using the radius of volume-equivalent sphere.
ARTICLE IN PRESS 78
Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
50
101
Silicate Pyroxene Amorphous carbon
30
10-1
Silicate Pyroxene Amorphous carbon
Qsca
Qabs
100
10-2 10-3
10 5 2
64
256 1024 4096 Number of monomers
64
256 1024 4096 Number of monomers
Fig. 9. Log–log plots of Qabs and Qsca of aggregates composed of up to 8192 monomers, where xm ¼ 4:5 for silicate, pyroxene and amorphous carbon. The lines are the regression line for each material. The aggregates with N mono p512 are calculated without grouping while those with N mono X1024 are calculated with the GAM for one direction of the incident light. Table 4 lists the equations of the regression lines.
Table 4 Equations of regression lines in Fig. 9 for silicate (Si), pyroxene (Py) and amorphous carbon (Ac)
Si Py Ac
Qabs
Qsca
y ¼ 0:000811 x0:292 (0.998) y ¼ 0:133 x0:287 (0.997) y ¼ 0:831 x0:298 (0.995)
y ¼ 2:35 x0:282 (0.990) y ¼ 1:82 x0:293 (0.994) y ¼ 0:864 x0:304 (0.996)
The numbers in brackets represent the correlation coefficient of the regression, and the x and y represent N mono and Qabs (or Qsca ), respectively.
Table 4 lists the equations for the regression lines in Fig. 9 and their corresponding correlation coefficients, all of which exceed 0.99. It is also notable that the regression lines have similar exponential slopes of around 0.293, commonly seen not only for all three materials but also both for Qabs and Qsca . The increase in the efficiencies with increasing N mono is related to the increase in multiple scattering by the monomers. Therefore, the structure of the fractal aggregates is considered to be related to the exponential slope of efficiencies for an increasing of number of monomers. However the slopes for aggregates with different fractal dimensions are likely to differ from those of BCCAs. 5. Conclusions We have proposed the GAM to derive the light scattering properties of a large fractal aggregate composed of spherical monomers where two conditions are met: (1) the scattering properties of each monomer can be derived and (2) in the light scattering by the monomer, the forward scattering is dominant. The studies for errors in the results computed by the GAM have been examined for the aggregates consisting of available number of monomers in our computer system, as functions of the number of groups, the sizes of neighbor region influenced on certain group (buffer region), sizes and optical constants of a constituent monomer. It is found that the GAM can keep the errors in Qabs and Qsca less than 15% and 10%, respectively, when the effect of 1/8 monomers of total constituent monomers N mono in the aggregate is taken into account for computing scattering properties of each monomer. This implies that the aggregate consisting of monomers eight times of available number of monomers in our computer system can be treated by GAM within these errors in Qabs and Qsca . In general, the GAM allows us to expand the number of monomers beyond the limit of computer resources. We have also found the power law relationship of Qabs and Qsca for N mono with the similar slopes for three different materials and correlation coefficients larger than 0.99.
ARTICLE IN PRESS Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
79
In this study, we have investigated the errors for BCCA. There is another fractal aggregate named BPCA (ballistic particle-cluster aggregate) often used in the light scattering calculation. BPCA has the fractal dimension of 2.980.02 [7]. Since the composing particles are located closer to each other than BCCA, the side scattering becomes more important. We have also investigated our method for BPCA. The errors are larger as expected (not shown here). We consider that our method works effectively only for porous particle described as BCCA. In this study, we used the cluster T-matrix code developed by Mackowski and Mishchenko [18]. However, the combinations with other methods are possible as long as the above conditions (1) and (2) are met. For example, the combination with the DDA method may allow us to calculate the asymmetry parameter g [45]. The combination with a1 term method may allow us to treat aggregates composed of a larger number than cluster T-matrix method for xm around 1.0 with errors by a1 term method and those by the GAM method [16,17]. The combination with the GMM method may allow us the calculation of aggregates of spheroidal monomers even though the applicability of the GAM method should be investigated for randomly directed spheroidal monomers [46,47]. Then these methods receive advantages in the increased number of monomers as conducted in this paper for cluster T-matrix method. The sample computer code of the GAM for cluster T-matrix method can be found at http:// harbor.scitec.kobe-u.ac.jp/okada/GAM/index.html. Acknowledgments This study used the cluster T-matrix code developed by Dr. Mackowski and Dr. Mishchenko (obtained from http://www.giss.nasa.gov/crmim/). Some of the data analyses were performed on the general common use computer system at the Astronomical Data Analysis Center (ADAC) of the National Astronomical Observatory of Japan. We thank Prof. Yoshitsugu Nakagawa at Kobe University, for comments that helped to improve this study. This research is supported by ‘‘The 21st Century COE Program of Origin and Evolution of Planetary Systems’’ of the Ministry of Education, Culture, Sports and Technology (MEXT).
References [1] Giese RH, Weiss K, Zerull RH, Ono T. Large fluffy particles—a possible explanation of the optical properties of interplanetary dust. A&A 1978;65:265–72. [2] Zerull RH, Giese RH, Weiss K. Scattering functions of nonspherical dielectric and absorbing particles VS Mie theory (E). Appl Opt 1977;16:777–8. [3] Weiss-Wrana K. Optical properties of interplanetary dust—Comparison with light scattering by larger meteoritic and terrestrial grains. A&A 1983;126:240–50. [4] Gustafson BA˚S. Microwave analog to light scattering measurements: a modern implementation of a proven method to achieve precise control. JQSRT 1996;55:663–72. [5] Gustafson BA˚S. Microwave analog to light scattering measurements. In: Mishchenko MI, Hovenier JW, Travis LD., editors. Light scattering by nonspherical particles: theory, measurements, and geophysical applications. New York: Academic Press; 1999. p. 367–90 [Chapter 13], 1999. ISBN 0-12-498660-9. [6] Meakin P. Effects of cluster trajectories on cluster–cluster aggregation: a comparison of linear and Brownian trajectories in two- and three-dimensional simulations. Phys Rev 1984;A29:997–9. [7] Mukai T, Ishimoto H, Kozasa T, Blum J, Greenberg JM. Radiation pressure forces of fluffy porous grains. A&A 1992;262:315–20. [8] Mann I, Okamoto H, Mukai T, Kimura H, Kitada Y. Fractal aggregate analogues for near solar dust properties. A&A 1994;291:1011–8. [9] Manickavasagam S, Menguc MP. Scattering matrix elements of fractallike soot agglomerates. Appl. Opt. 1997;36:1337–51. [10] Mackowski DW. A simplified model to predict the effects of aggregation on the absorption properties of soot particles. JQSRT 2006;100:237–49. [11] Wurm G, Blum J. Experiments on preplanetary dust aggregation. Icarus 1998;132:125–36. [12] Blum J, Wurm G. Experiments on sticking, restructuring, and fragmentation of preplanetary dust aggregates. Icarus 2000;143:138–45. [13] Mann I, Kimura H, Kolokolova L. A comprehensive model to describe light scattering properties of cometary dust. JQSRT 2004;89:291–301. [14] Draine BT. The discrete dipole approximation and its application to interstellar graphite grains. ApJ 1988;333:848–72.
ARTICLE IN PRESS 80
Y. Okada et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 108 (2007) 65–80
[15] Draine BT, Flatau PJ. User guide to the discrete dipole approximation code DDSCAT 6.1, 2004. hhttp://arxiv.org/abs/astro-ph/ 0409262i. [16] Okamoto H. Light scattering by clusters: the A1-term method. Opt Rev 1995;2:407–12. [17] Okamoto H, Xu Y. Light scattering by irregular interplanetary dust particles. Earth Planets Space 1998;50:577–85. [18] Mackowski DW, Mishchenko MI. Calculation of the T matrix and the scattering matrix for ensembles of spheres. J Opt Soc Am A 1996;13:2266–78. [19] Mishchenko MI, Mackowski DW. Electromagnetic scattering by randomly oriented bispheres: comparison of theory and experiment and benchmark calculations. JQSRT 1996;55:683–94. [20] Sun W, Loeb NG, Fu Q. Finite-difference time-domain solution of light scattering and absorption by particles in an absorbing medium. Appl Opt 2002;41:5728–43. [21] Sun W, Loeb NG, Videen G, Fu Q. Examination of surface roughness on light scattering by long ice columns by use of a twodimensional finite-difference time-domain algorithm. Appl Opt 2004;43:1957–64. [22] Yang P, Liou KN. Finite difference time domain method for light scattering by nonspherical and inhomogeneous particles. In: Mishchenko MI, Hovenier JW, Travis JD, editors. Light scattering by nonspherical particles. San Diego: Academic Press; 2000. p. 173–221. [23] Chylek P, Srivastava V. Dielectric constant of a composite inhomogeneous medium. Phys Rev B 1983;27:5098–106. [24] Macke A. Scattering of light by polyhedral ice crystals. Appl Opt 1993;32:2780–8. [25] Muinonen K, Nousiainen T, Fast P, Lumme K, Peltoniemi JI. Light scattering by Gaussian random particles: ray optics approximation. JQSRT 1996;55:577–601. [26] Okada Y, Nakamura AM, Mukai T. Light scattering by particulate media of irregularly shaped particles: laboratory measurements and numerical simulations. JQSRT 2006;100:295–304. [27] Mukai T, Okada Y. Optical properties of large aggregates. In: Kru¨ger H, Graps A., editors. Dust in planetary systems, Workshop, September 26-30, 2005, Kauai, Hawaii, ESA Publications, SP-643; January 2007. [28] Andersen AC, Mutschke H, Posch Th, Min M, Tamanai A. Infrared extinction by homogeneous particle aggregates of SiC, FeO and SiO2 : Comparison of different theoretical approaches. JQSRT 2006;100:4–15. [29] Koehler M, Kimura H, Mann I. Applicability of the discrete-dipole approximation to light-scattering simulations of large cosmic dust aggregates. A&A 2006;448:395–9. [30] Taflove A, Hagness SC. Computational electrodynamics: The finite-difference time-domain method (Artech House Antennas and Propagation Library). MA: Artech House; 2005. 1038p. [31] Kozasa T, Blum J, Mukai T. Optical properties of dust aggregates, I Wavelength dependence. A&A 1992;263:423–32. [32] Wolff MJ, Clayton GC, Martin PG, Schulte-Ladbeck RE. Modeling composite and fluffy grains: the effects of porosity. ApJ 1994;423:412–25. [33] Wolff MJ, Clayton GC, Gibson SJ. Modeling composite and fluffy grains II Porosity and phase functions. ApJ 1998;503:815–30. [34] Voshchinnikov NV, Mathis JS. Calculating cross sections of composite interstellar grains. ApJ 1999;526:257–64. [35] Voshchinnikov NV, Il’in VB, Henning Th. Modelling the optical properties of composite and porous interstellar grains. A&A 2005;429:371–81. [36] Voshchinnikov NV, Il’in VB, Henning Th, Dubkova DN. Dust extinction and absorption: the challenge of porous grains. A&A 2005;445:167–77. [37] Koehler M, Mann I. Light-scattering models applied to circumstellar dust properties. JQSRT 2004;89:453–60. [38] Liou KN, Hansen JE. Intensity and polarization for single scattering by polydisperse spheres: a comparison of ray optics and Mie theory 1971;28:995–1004. [39] Mackowski DW. Calculation of total cross sections of multiple-sphere clusters. J Opt Soc Am A 1994;11:2851–61. [40] McQueen J. Some methods for classification and analysis of multivariate observations. In: Proceedings of the fifth Berkeley symposium on mathematical statistics and probability, 1967. p. 281–97. [41] Mukai T. Cometary dust and interplanetary particles. In: Bonetti A, Greenberg JM, Aiello SA, editors. Evolution of interstellar dust and related topics. Amsterdam: Elsevier; 1989. p. 397–445. [42] Dorschner J, Begemann B, Henning T, Jaeger C, Mutschke H. Steps toward interstellar silicate mineralogy. II Study of Mg–Fesilicate glasses of variable composition. A&A 1995;300:503–19. [43] Preibisch T, Ossenkopf V, Yorke HW, Henning Th. The influence of ice-coated grains on protostellar spectra. A&A 1993;279:577–88. [44] Bohren CF, Huffman DR. Absorption and scattering of light by small particles. New York: Wiley; 1983 530p. [45] Hoekstra AG, Frijlink M, Waters LBFM, Sloot PMA. Radiation forces in the discrete-dipole approximation. J Opt Soc Am A 2001;18:1944–53. [46] Xu Y-l. Electromagnetic scattering by an aggregate of spheres. Appl Opt 1995;34:4573–88. [47] Xu Y-l. Electromagnetic scattering by an aggregate of spheres: far field. Appl Opt 1997;36:9496–508.