Geometric modeling of aerosol surface roughness and numerical simulation of its impact on aerosol optical properties

TL;DR

This paper introduces a stochastic geometric modeling method for aerosol surface roughness and evaluates its impact on optical scattering properties using advanced computational techniques, improving understanding of aerosol-light interactions.

Contribution

It presents a novel stochastic modeling approach for aerosol surface roughness and applies it to compute optical properties with enhanced accuracy over smooth models.

Findings

Rough surface models outperform smooth particles in scattering simulations.

Entropy measures effectively distinguish scattering similarities.

The method covers particles from Rayleigh to geometric optics regimes.

Abstract

A stochastic method is introduced for geometric modeling aerosol surface roughness with random field and discrete differential geometry theory. Optical scattering properties are computed for randomly oriented spheroidal particles with uniformly random surface roughness. Invariant imbedding T-Matrix and geometric optics method are applied to compute light scattering properties of aerosol particles covered from Rayleigh to geometric optics region. We simulatied scattering phase matrix of model particles, better performace is obtained than smooth particles. We aslo introduce introduce entropy and relative entropy as measures to distinguish the optical scattering similarities between different particles.

Equations

$$h(\mathbf{x})=\sum _{\mathbf{k}}A(|\mathbf{k}|)\cos (\mathbf{k}\cdot \mathbf{x}+{\varphi }_{\mathbf{k}})$$

$$\mathcal{H}=\frac{1}{2}\nu (\nabla h{)}^2+V(h)$$

$$\frac{\partial }{\partial t}h(\mathbf{x},t)=\nu {\nabla }^{2}h-ah+\eta (\mathbf{x},t)$$

$$\mathbb{E}[\eta (\mathbf{x},t)]$$

$$\mathbb{E}[\eta (\mathbf{x},t)\eta ({\mathbf{x}}^{\prime },t^{\prime })]$$

$$h(\mathbf{k},t)=h(\mathbf{k},0)e^{-(\nu k^2+a)t}+{\int }_0^{t}d\tau e^{-(\nu k^2+a)(t-\tau )}\eta (\mathbf{k},\tau )$$

$$\mathbb{E}[h(\mathbf{k},t{)}^{\ast }h({\mathbf{k}}^{\prime },t)]=(2\pi {)}^2{\delta }^2(\mathbf{k}-{\mathbf{k}}^{\prime })(1-e^{-2(\nu k^2+a)t})\frac{T}{\nu k^2+a}$$

$$\mathbb{E}[h(0,t)h(\mathbf{x},t)]=2\pi T\int \frac{J_0(kr)(1-e^{-2(\nu k^2+a)t})}{\nu k^2+a}kdk$$

$$\mathbb{E}[h(\mathbf{x},t)h(\mathbf{x},t)]$$

$$=$$

$$h(\mathbf{k},\omega )=\chi (\mathbf{k},\omega )\eta (\mathbf{k},\omega )$$

$$C(\mathbf{k},\omega )=\frac{2T}{\omega }\Im \chi (\mathbf{k},\omega )$$

$$\tilde{h}=\sum _j{h}_j{\varphi }_j(\mathbf{x})$$

$$⟨\nabla {\varphi }_i,\nabla \tilde{h}⟩=\sum _j{h}_j⟨\nabla {\varphi }_i,\nabla {\varphi }_j⟩$$

$$⟨\nabla {\varphi }_i,\nabla {\varphi }_j⟩={\int }_{\Omega }dx\nabla {\varphi }_i(x)\cdot \nabla {\varphi }_j(x)$$

$$⟨\nabla {\varphi }_k,\nabla {\varphi }_k⟩=\frac{1}{2}(\cot \alpha +\cot \beta )$$

$$⟨\nabla {\varphi }_i,\nabla {\varphi }_j⟩=-\frac{1}{2}\cot \gamma$$

$$({\nabla }^{2}h{)}_i=\frac{1}{2}\sum _j(\cot {\alpha }_j+\cot {\beta }_j)(h_j-h_i)$$

$$h_i^{n+1}=h_i^n+\frac{1}{2}\nu \delta t\sum _j(\cot {\alpha }_j+\cot {\beta }_j)(h_j^n-h_i^n)+{\eta }_i\delta t$$

$$\frac{x^2+y^2}{b^2}+\frac{z^2}{a^2}=1$$

$$a>b:\text{prolate spheroid}$$

$$a<b:\text{ oblate spheroid}$$

$$A_r=\mathbb{E}[\int (1+|\nabla h{|}^2{)}^{1/2}d\mathbf{x}]$$

$$P(\mu )=\left[\begin{array}{cccc}{P}_{11}& P_{12}& 0& 0\\ P_{12}& P_{22}& 0& 0\\ 0& 0& P_{33}& P_{34}\\ 0& 0& -P_{34}& P_{44}\end{array}\right]$$

$$\frac{1}{2}\int P_{11}(\mu )d\mu =1$$

$$S=-\int d\mu P_{11}(\mu )\ln [P_{11}(\mu )]$$

$$KL(P_{11}||P_{11}^{\prime })=\int P_{11}(\mu )\ln [P_{11}(\mu )/P_{11}^{\prime }(\mu )]d\mu$$

$$P(\mu )=\frac{{\sum }_{s,X}{n}_s{n}_X{\sigma }_{sc}(s,X)P(\mu ,s,X)}{{\sum }_{s,X}{n}_s{n}_X{\sigma }_{sc}(s,X)}$$

Full text

1 Introduction

Light scattering by small particles is of great importance both in scientific research[1], and in industrial technology[2][3], examples include dusts, cloud particles, ocean particles, nano metal particles et. al.. In order to characterize these particles, especially dusts, it is necessary to know how they interact with light.

To obtain dust optical properties, light measurements and numerical simulations can be used. In terms of numerical scattering studies, perhaps the simplest geometry is isotropic homogeneous sphere, in which Lorenz-Mie theory[1] is a standard tool for analysis and still very active. Not only has qualitative insights about scattering properties been acquired with this method, the corresponding analytical solutions serve as the benchmark for testing the accuracy and convergence behaviors of various computation software packages of light scattering. Based on Lorenz-Mie theory, light scattering properties of a large class of natural or artificial small particles[1] has been modeled assuming these scatterers to be homogeneous spheres suspending in some medium.

Despite the tremendous success of sphere in modeling of dust optical properties, any extension to non-spherical geometry has been proven to be quite useful[4]. For example, spheroids are frequently used for modeling dust, among all non-spherical model particles. Researchers[5, 6, 7, 8, 9, 10] are fascinated by spheroids due to their ability in recovering light scattering properties of dust aerosols with such dissimilar appearance. Counterintuitively, other irregular model particles[11, 12] generated in random manners, like Gaussian random spheres[13], do not perform as good as agreement with laboratory measurement as spheroids in optical modeling of small natural particles. On the other hand, all real particles’s surface are roughened to some degree. Compared with other factors like shapes, chemical composition, less studies[14, 15, 16] focus on surface roughness effects on optical properties of particles. One difficulty is generation of reasonable surface roughness. Most current light scattering researchers construct particle’s surface roughness texture using Gaussian random field(GRF)[17, 18] on spheres or prisms, which are based on random Fourier series expansion on spheres or finite planes. Another noticeable effort in modeling surface roughness is given by the random tilt model[19], in which the local surface normal at the interception of incident ray with the particle surface is randomly tilted with the Gaussian distribution. There are also some surface roughness models that fail to take the correlation of surface roughness into consideration.

In [20], light scattering by spheres with surface roughness has been investigated using Gaussian random sphere from Rayleigh to geometric optics region. Considering the common interest in spheroids and surface roughness effects, light scattering by large spheroids with surface roughness modeled using random tilt model has been studied in [21] with a improved geometric optics method(IGOM)[22]. Wave optics simulation[15] for non-absorptive roughened spheroids scattering was also reported.

Although considered as an important factor, current findings on surface roughness effects were controversial and even conflicting, partly due to lack of rigous roughness model that can characterize roughness effectively and filter out effects from host particles. This paper follows our previous work[15], in which we model rough surface geometry based on the theory of Gaussian random field on Riemann manifolds. The idea of this approach is to link the physical mechanism of surface growth to geometric modeling of surface roughness. To study the physical mechanism of surface roughness, several experimental and numerical researches has been conducted for ice particles[23]. In this paper, we simulate light scattering by randomly rough spheroids in our new Gaussian random field model using both T-Matrix and ray tracing method. We propose the entropy as optical similarity measure between particles. We apply this new rough model particles to simulate ensemble average scattering properties and compare the results with measurement data of dust particles.

2 Modeling rough dust surfaces

Natural particles vary in shapes, of which most carry obvious surface roughness. In this section, we will discuss geometric modeling of roughness on arbitrary 2D surfaces based on random field and discrete differential geometry theory. A random field $h(\mathbf{x})$ in $\mathbb{R}^{n}$ is a function whose value are random variables for $\mathbf{x}$. Specifically, a GRF $h(\mathbf{x})$ (height function of the surface) can arise from the random superposition of harmonic waves:

[TABLE]

where $A(|\mathbf{k}|)$ is an amplitude spectrum, sampled from some probability distribution. The random phases $\phi_{\mathbf{k}}$ are uniformly independently distributed within the range $[0,2\pi)$. The power spectrum $C(|\mathbf{k}|)=A^{2}(|\mathbf{k}|)$, containing the two-point correlation information of random fields, entirely code statistical properties of $h(\mathbf{x})$.

For generating GRF, the convolution method above is limited to simple geometry like plane or sphere, since there are well defined harmonic waves(plane, cylindrical or spherical waves). In order to generate GRF on an arbitrary 2D surface, we construct random fields through stochastic differential equation(SDE) approach[24].

In the following, we will briefly review stochastic surface growth theory[25, 26] on 2D plane first and derive a correlation length which is a reasonable approximation to general 2D surfaces with small scale roughness. The energy density functional of a surface can be given by

[TABLE]

where $\nu$ represents ”surface diffusion” and $V(h)$ potential energy. The corresponding linear Langevin equation is

[TABLE]

where $\eta$ is a Gaussian random variable, $V(h)=ah^{2}$ and

[TABLE]

where $\mathbb{E}[\cdot]$ represents expectation value, $\sigma$ is the standard deviation of the surface height function. After making a Fourier transform, the solution at time $t$ can be obtained as

[TABLE]

If $h(\mathbf{k},0)=0$, the average of surface height at time $t$ will be 0, and the covariance function of the surface height $h$ grows as

[TABLE]

After transforming to the coordinate space, we obtain

[TABLE]

Let $\Lambda=1/l_{0}$ be the cutoff. When $a\to 0$, we get

[TABLE]

where we have a correlation length $l_{c}=\sqrt{2\nu t}$[25].

[TABLE]

and $\chi=\frac{1}{-i\omega+\nu k^{2}+a}$ is the suseptibility. According to the fluctuation-dissipation theorem, we obtain correlation function:

[TABLE]

Now, we extend the SDE approach to generate GRFs on general 2D closed surfaces for the purpose of modeling particle surfaces. For simplicity, we adopt Edwards-Wilkinson model(a=0 in eq. 3). Following the approach given in[27], we give a review on how to approximate the differential equation on a general 2D surface using discrete differential geometry theory. Explicitly, suppose $\tilde{h}$ is the approximation of the true solution $h$, then:

[TABLE]

where $\{\phi_{j}\}$ is a collection of expansion functions. If we apply Green’s identity to the triangle mesh, the integral is breaking up into a sum over individual triangles:

[TABLE]

where each component is given by

[TABLE]

Within a given triangle, each expansion function associated with the corresponding vertex

[TABLE]

where $\alpha$ and $\beta$ are the remaining two angles of that given triangle (see Fig. 1). Similarly, for the expansion functions $\phi_{i}$ and $\phi_{j}$ associated with two connected vertices along the same edge, we have

[TABLE]

Now, we can obtain discretized version of Laplacian of arbitrary 2D surfaces, which can be approximated as

[TABLE]

where the sum is taken over the connected neighbors of vertex i.

[TABLE]

Next, we describe how to generate uniform roughness on spheroid. In our simulation, dust aerosol particles are first modeled with the smooth spheroid, of which with $z$-axis as the symmetry axis the equation is given by

[TABLE]

where $a$ is the distance from centre to either pole along the symmetry axis, and he equatorial radius of the spheroid is the semi-major axis $b$. The aspect ratio is defined as $\epsilon=a/b$. And we have

[TABLE]

Then we discretize the smooth spheroid surface with triangular mesh using DistMesh[28]. The normal direction at each vertex on the mesh is defined as the average of adjacent triangular facets. The surface height function $h$ is defined as the normal distance deviation from smooth surface. The simulation is performed according the above equation. The deviation and correlation of height function can be controlled by total simulation time.

Fig. 2 illustrates a typical sample of rough spheroids. And if we keep the correlation length fixed, the standard deviation of surface height $w=\sqrt{\mathbb{V}ar(h)}$ can be used to characterize the degree of surface roughness. The average surface area of a particle is written as

[TABLE]

So we introduce another characteristic of roughness, which is defined as the ratio of rough and smooth surface area. In particular, it has been noted that the above equation becomes $A_{r}\simeq A_{s}+\frac{1}{2}\mathbb{E}[\int|\nabla h|^{2}d\mathbf{x}]$ for small roughness $|\nabla h|\ll 1$ ($A_{s}$ is the area of smooth base particle.), which unveils the connection between the average surface area of rough surface and its average local slope. This roughness characteristic also connects to physical properties of natural particles. For example, the surface area of a large particle is proportional to its radiation power according to the Stefan-Boltzmann law.

3 Relevant light scattering theory and information entropy

Before illustrating numerical scattering simulation with our new rough particles, we would like to discuss the relevant scattering theory first. We will suggest to introduce entropy concept to measure the similarity between particles’ scattering properties. The Stokes parameters give a complete description of light fields in the sense of principle of optical equivalence. The presence of obstacles in the path of light fields produces light scattering. At the large distance, the scattered field moving in a given direction can be considered as a local plane wave field. The Stokes parameters of the incident field and scattered field are related by scattering phase matrix. Under random orientation and mirror reflection symmetry assumptions[29], the scattering matrix is block diagonal (only 6 elements are independent)

[TABLE]

where $\mu=\cos\theta_{sc}$, $\theta_{sc}$ is the scattering angle between incident direction and scattered direction of light. Phase matrix $P$ are defined with the normalization condition:

[TABLE]

Given the frequency of light, the phase function $P_{11}$ of a scatterer can be interpreted the angular distribution of light field intensity scattered by the object.

As we had mentioned above, random irregular particles usually do not fit measured scattering phase matrix of natural particles like dust aerosols as well as spheroids. There are long and severe debates on why spheroids model optical properties of real particles so well. In this paper, we propose a potential tool to solve this debate by introducing the concept of scattering information entropy. In addition, we note that the entropy of 2D image was proposed for characterizing particle roughness in light scattering experiments[30]. Generally, Shannon’s information entropy for any probability distribution of random variable is defined as

[TABLE]

Since scattering phase function is indeed a probability density distribution of scattered directions. For each phase function, we can calculate its entropy, the degree of disorder. Overall, for larger particles, the higher peak in forward direction will correspond to lower entropy. To quantify optical similarity between different particles, we introduce the distance between the two distributions $P_{11}(\mu)$ and $P^{\prime}_{11}(\mu)$ is given by the relative entropy or more formally, Kullback-Leibler divergence[31]

[TABLE]

As an example, the scattering phase matrix of small spheroidal particles with refractive index $n=1.31$ is computed with the invariant imbedding T-matrix method (II-TM) [32, 33]. Fig. 3 (top) shows the Shannon entropy as a function of size parameter $X=ka$(k is the wavenumber, and a is maximum dimension of the corresponding particle.). The pattern shows the entropy will decrease in trend but fluctuate around certain sizes. The amplitude of fluctuations decay gradually as particles’ size parameter increases. It seems reasonable that entropy curves will become much smoother when particles’ size go to the geometrical optics region. For $X>20$, the curve of Shannon entropy begin to split. The surface roughness increase the Shannon entropy for each size.

In Fig. 3 (bottom), it shows that as standard deviation $w$ (roughness metric) increases, so does the KL-divergence between rough spheroids and the corresponding smooth ones. It also indicates that as size parameter increases, the KL-divergence grows too. Therefore, KL-divergences may help us understanding the mysterious success of spheroids, optical distance(KL-divergence) among spheroid ensemble are generally larger than random irregular particle ensemble, detailed study about this will be presented in our future work. In other words, functional spaces of phase matrix elements for spheroids span larger than that of random irregular particles.

4 Numerical simulation

In the simulation, we adopt rough spheroids described in section 2 as model particles. Light scattering by particles with our new random rough surfaces are simulated using II-TM and IGOM. First, II-TM computes particles’ optical properties based on volume integral equation formalism of Maxwell’s equations[32], in which the invariant imbedding technique is used to solve the volume integration efficiently. Compared with surface integral formalism of T-Matrix method, II-TM can deal with highly non-spherical and inhomogeneous cases. Like other T-matrix methods[4], II-TM can obtain random orientation cross sections and scattering phase matrix analytically. For large size particles, T-Matrix method falls down. Hence, IGOM which hybridizes geometric optics (ray-tracing) and physical optics methods, is utilized to perform the computation. And each ray in IGOM is generated from a uniform random incident direction, random oriented results are obtained after we perform the near field to far field mapping[19].

For particle with size parameter $X=60$, it takes about 25 hours wall time with 100 GB memory usage for a single II-TM simulation. The computing platform used is a IBM NeXtScale nx360 M4 dual socket server based on the Intel Xeon 2.5GHz E5-2670 v2 10-core processor, 2.5 GHz and 256GB RAM, 20 cores on the node. Fig. 4 shows the phase matrix elements for randomly rough spheroids with different roughness degree $w$ ( size parameter $X=60$, aspect ratio $\epsilon=a/b=2$). It shows roughness has little impact on scattering phase function for angles $<5^{\circ}$. While for angles $>40^{\circ}$, roughness will enhance scattering phase function $P_{11}$, especially for the backward direction $180^{\circ}$, which means backscattering is quite sensitive to surface roughness. For other phase matrix elements like $P_{22}/P_{11}$,surface roughness changes the curve dramatically, all local extrema are smoothed out. Generally, increasing surface roughness will make the curve smoother. Increasing surface roughness, variance of the curves of phase matrix elements decreases.

Fig. 5 shows phase matrix elements of the same sample particle as above with IGOM.The figure shows consistent moving trend for these two methods.

We examine the performance of randomly rough spheroids by comparing modeling results with the experimentally measured data for feldspar ensemble[34]. To obtain model scattering phase matrix $P(\mu)$, we compute its ensemble average with

[TABLE]

where $X$ is the particle size parameter, $s$ represents its shape, $n_{X},n_{s}$ are particle size and shape probability mass function and $\sigma(s,X)$ is the scattering cross section for a particle with shape $s$ and size parameter $X$. Model particle size range was set according to the measured data [34]. Furthermore, we assumed that spheroid shape (parameterized with aspect ratio $\epsilon$, $z=\epsilon$ if $\epsilon>1$, and $z=-1/\epsilon$ if $\epsilon<1$) follows a power distribution $p(z)=C|z|^{3}$(C is the normalization constant.), which was also used in [35].

Fig.6 and Fig. 7 indicate that the surface roughness has great impact on phase matrix of spheroidal particles, especially larger size. The following comparison with measurements also shows rough model particles have better performance than smooth ones. The optical scattering of dust particles is size dependent. We simulate light scattering by dust aerosols within size parameter ranging from 1 to 300. In the simulation, II-TM is used for size smaller than 60 and IGOM is used for size larger than 60. The estimated refractive index is set to be $m=1.5+0.001i$, and particle size distribution is statistically estimated by observations given in [34]. Specifically, we use a volume distribution for dust particles in our comparison. We compare our model results with the measured data from the Amsterdam Light Scattering Database[34]. Fig. 8 shows model data versus measured data. The model data show good agreement with measurements. Randomly rough spheroidal particles reproduce most of the phase matrix elements, and give better fit for observed dust scattering data.

5 Conclusion

In conclusion, a method applying discrete differential geometry method to model dust particles’ surface roughness is discussed in detail. This method applies a novel approach for generating rough surfaces based on stochastic differential equation and finite element. When applied to particle surface, correlation of surface random field can be estimated. In particular, the correlation length of an random field on particle surface is obtained with a quasi planar approximation. The mean and variance of the surface random field on each particle can be estimated using classic statistical inference methods.

Furthermore, information entropy and KL-divergence are introduced for quantitatively distinguishing optical differences between two particles, which can be used as a similarity measure. Application to remote sensing of atmospheric particles and optical characterization of small particles are also possible. Light scattering by rough spheroids are simulated with II-TM and IGOM software package. Finally, ensemble modeling results are compared with the laboratory measured scattering phase matrix data of feldspar dust.

6 Acknowledgments

We would like to thank Ping Yang for helpful discussions and reviewers for their useful suggestions for improvements. This study was supported by the National Natural Science Foundation of China (NSFC) (Grant No. 41705010) and “The Fundamental Research Funds for the Central Universities”(DUT16RC(3)121).