Sand/dust storm is a natural phenomenon that may cause signals of remote sensing and telecommunication being worsened due to the attenuation of electromagnetic waves (EMWs) by sand particles. To understand the attenuation mechanism, much research has been conducted on the scattering of sand storm (Peterson and Strm, 1974; Ansari and Evans, 1982; Haddad et al., 1983; Ghobrial and Jervase, 1997; Islam et al., 2010) with the Lorenz–Mie theory. Notably, particles in the sand storm can be charged due to many reasons, such as friction and collision between particles (Gill, 1948;Gill and Alfrey, 1949; Lowell and Rose-Innes, 1980; Lacks and Mohan Sankaran, 2011; Xie et al., 2013a; Xie et al., 2013b). Existing research has shown that the scattering properties of charged particles are strikingly different from those of neutral ones (Klačka and Kocifaj, 2007; Zhou and Xie, 2011; Kocifaj, 2013; Klačka et al., 2015b; Kocifaj et al., 2015). For example, Klačka et al. found that scattering and absorption will be significantly enhanced by the surface charge for small particles with a size parameter less than 0.01 (Klačka and Kocifaj, 2007). In 1977, Bohren and Hunt put forward a model to solve the scattering problem of a charged particle through revising the boundary condition by introducing an excess charge surface conductivity (Bohren and Hunt, 1977; Bohren and Huffman, 1998). The Drude model that relates the surface conductivity of charged particles to the charge-to-mass ratio of the charge carriers and the driving frequency can be employed to understand the surface currents on the particles (Klačka and Kocifaj, 2007; Kocifaj, 2013; Klačka et al., 2015b).

However, although there is a large number of research conducted on the scattering properties of individual charged particles of different materials and shapes (Klačka and Kocifaj, 2010; Kocifaj et al., 2016; Liao et al., 2019; Yang et al., 2019; Xue et al., 2020), little attention was given on the interactions between charged particles in clusters. It is reported that conventional analytical methods, which assume the granular medium is homogeneous and neglect multiple scattering effects between particles, will lead to significant overestimation of the attenuation through the medium (Huang et al., 2020; Xie et al., 2020). Recently, Xie et al. improved the independent scattering approximation in conventional radiative transfer theory to deal with the multiple scattering effects in a charged granular system (Xie et al., 2020). In 2020, Chai J W et al. analytically studied the optical response of a symmetric circular array of metallic nanoparticles based on the eigen-decomposition method under the dipole approximation (Fung and Chan, 2008; Chai et al., 2020). In this article, we build up a fully analytical theoretical frame based on the T-matrix formulism (Tsang et al., 2001; Mishchenko et al., 2002, 2006; Huang et al., 2020) to study the electromagnetic response of multi-particle systems consisting of electrically charged spherical particles with multiple scattering effect considered. The article is organized as follows: In the Single charged particle section, before the discussion of the multiple scattering effect, we investigated the influence of various factors on the scattering of a single charged particle, such as the charge density, radius, and relative permittivity of the particle. We would like to state that our theoretical frame is under the assumption of uniform charge distribution on the particle surface. Uniform charge distribution is typical for particles of simple geometry such as spherical particles considered in this work, to which separation of variables method is applicable to find scattered electric and magnetic fields (Klačka et al., 2015a). In the Clustered charged particles section, the T-matrix formulism with modified boundary condition to calculate the electromagnetic response of clustered charged spherical particles is described. Then, in the Result section, the discussions on the numerical results of a single charged particle and clustered particles are presented, with comparison between different charge densities and inter-particle distances. Finally, the convergence performance of the proposed theoretical approach is investigated, and conclusions are given in the Conclusion section.

A charged spherical particle with radius a is illuminated by x-polarized EMW propagating along the z-axis. The incident electric field can be expressed as E → i n = e i ( k z − ω t ) x ^ , where ω = 2 π f and f denotes the frequency of the EMW. Here, x ^ is the unit vector of the x-axis. The particle is assumed to be overall and uniformly charged. The material of the particle and ambient medium are set to be isotropic and homogeneous. Let E → 1 and H → 1 denote the field inside the particle, and E → 2 and H → 2 denote the field of the ambient. The electric field E → i (i = 1, 2) and magnetic field H → i (i = 1, 2) satisfy the following wave propagating equations: { ∇ 2 E → 2 + k 2 2 E → 2 = 0 ,     ∇ 2 H → 2 + k 2 2 H → 2 = 0 ∇ 2 E → 1 + k 1 2 E → 1 = 0 ,     ∇ 2 H → 1 + k 1 2 H → 1 = 0 (1) where k ₁ and k ₂ denote the wavenumber inside and outside the particle, which satisfy the dispersion relation: k i = ω ε i μ i (2) where μ i ( i = 1 ,   2 ) and ε i ( i = 1 ,   2 ) are the permeability and permittivity of the particle and the ambient, respectively. We consider the medium outside the particle as a vacuum, i.e., μ 2 = μ 0 and ε 2 = ε 0 . Here, μ 0 = 1.2566 × 10 − 6   H / m and ε 0 = 8.8542 × 10 − 12   F / m are the vacuum permeability and permittivity, respectively. μ 1 = μ r μ 0 and ε 1 = ε r ε 0 are the permeability and permittivity of the particle, respectively, where μ r and ε r are the relative permeability and relative permittivity, respectively. We begin with the Mie scattering theory, by which the problem of plane waves scattered by a spherical particle can be rigorously solved through matching the boundary conditions (Kong, 1986; Bohren and Huffman, 1998; Mishchenko et al., 2002; Lovat et al., 2020). The field outside the particle can be expressed as the superposition of the incident field and the scattered field, { E → 2 = E → i n + E → s H → 2 = H → i n + H → s (3) where E → s and H → s denote the scattered fields of the particle. It is commonly assumed in the Mie theory when studying the scattering of electromagnetic waves by spherical particles that the particle carries a zero charge. The boundary condition can be written as: { ( E → 2 − E → 1 ) × n ^ = 0 ( H → 2 − H → 1 ) × n ^ = 0             a t     r = a (4) where n ^ denotes the outward-directed unit vector normal to the surface of the particle.

However, when it comes to charged particles, the electrons on the surface will be driven by the electric field on the surface. Hence, there exist surface current K → on the surface of the particle due to the movement of excess charge driven by the tangential electric field. The boundary condition is modified as: { ( E → 2 − E → 1 ) × n ^ = 0 ( H → 2 − H → 1 ) × n ^ = K →           a t       r = a (5)

The surface current can be expressed as: K → = σ s E → τ (6) where σ _s denotes the excess charge surface conductivity, and E → τ stands for the tangential component of the electric field on the surface of the particle. σ _s is independent of position but may be complex and dependent on frequency. Under the assumption of free surface charges (i.e., not localized on the surface of the particle), we can adopt the Drude model and the equation of motion for a surface charge with mass m _s and charge q _s as m s   u → ˙ τ = q s E → τ − β s u → τ (7) where u → τ is the component of velocity tangential to the sphere, and β s is the damping coefficient. The oscillatory part of the solution of (Eq. 7) is u → τ = q s / m s − i ω + γ s E → τ (8) where γ s = β s / m s . Since K → = η s u → τ , where η _s is the density of charge, from (Eqs 6, 8), we have (Bohren and Hunt, 1977;Kocifaj, 2013;Klačka et al., 2015b;Holloway et al., 2019;zhang et al., 2020), σ s = η s ( q s / m s ) − i ω + γ s (9) where q _s/m _s is the charge-to-mass ratio of charge carriers. For electrons, we have q s / m s = 1.7587 × 10 11   C / kg . The constant γ _s can be calculated as γ s = k B T / ℏ (Klačka and Kocifaj, 2007), where k _B denotes the Boltzmann constant, and ℏ denotes the Planck’s constant.

The problem of plane waves scattered by a charged particle can then be solved based on a T-matrix-based approach with modified boundary conditions, which is the generalized version of the Lorenz–Mie theory. The scattering of electromagnetic waves by homogeneous, isotropic sphere, which carries excess charge is tackled under the assumption that surface current is linearly related to the electric field parallel to the surface of the particle. The calculated scattering coefficients are, therefore, formally different from those of the standard Mie theory. We first expand the incident wave into spherical harmonics, E → i n c = e i ( k z − ω t ) x ^ = e − i ω t { − ∑ n = 1 ∞ ∑ m =   − 1,1 i n 2 n + 1 n ( n + 1 ) [ x ^ ⋅ C → − m , n ( 0,0 ) γ m n R g M → m n ( k r , θ , ϕ ) − i x ^ ⋅ B → − m , n ( 0,0 ) γ m n R g N → m n ( k r , θ , ϕ ) ] } (10) where i =   − 1 . M → m n and N → m n are unit vectors of spherical waves, which can be written as R g M → m n ( k r , θ , ϕ ) = γ m n j n ( k r ) C → m n ( θ , ϕ ) (11) R g N → m n ( k r , θ , ϕ ) = 1 k ∇ × R g M → m n ( k r , θ , ϕ ) = γ m n { n ( n + 1 ) j n ( k r ) k r P → m n ( θ , ϕ ) + [ k r j n ( k r ) ] ′ k r B → m n ( θ , ϕ ) } (12) where j _n and h _n are the first and third kinds of spherical Bessel function. Rg stands for the regular solution of the spherical harmonics, which means the solution is regular (finite) at the origin. The vector spherical waves without the prefix Rg are the expressions with j _n replaced by h _n. γ m n can be calculated as γ m n =   ( 2 n + 1 ) ( n − m ) ! 4 π n n ( n + 1 ) ( n + m ) ! . B → m n , C → m n , and P → m n can be written as (Tsang et al., 2000) B → m n = [ θ ^ d d θ P n m ( c o s θ ) + ϕ ^ i m s i n θ P n m ( cos ⁡ θ ) ] e i m ϕ   (13) C → m n = [ θ ^ i m s i n θ P n m ( cos ⁡ θ ) − ϕ ^ d d θ P n m ( c o s θ ) ] e i m ϕ (14) P → m n =   r ^ Y n m ( θ , ϕ ) (15) where r ^ , θ ^ , and ϕ ^ are unit vectors under spherical coordinates. P n m ( c o s θ ) is associated Legendre polynomials, and Y n m ( θ , ϕ ) is spherical harmonics, which can be related by   Y n m ( θ , ϕ ) =   P n m ( c o s θ ) e i m ϕ .

Under the classical assumption of independent scattering, the cross section of the multi-particle system can be written as C s c a =   ∑ i C s c a i , where C s c a i is the scattering cross section of the ith particle and can be calculated from the method mentioned in the Single charged particle section. However, in multi-particle systems, particles are electromagnetically coupled by the multiple scattering effect. The incident field of one particle can be expressed as the addition of the original incident field and scattered field from other particles [i.e., E → e x ( l ) ( r → ) =   E → i n c ( r → ) + ∑ j ≠ l E → s ( j ) ( r → ) ].

A general theoretical frame to study the multiple scattering effect in a multi-particle system has been built in previous research (Tsang et al., 2001; Mishchenko et al., 2002, Mishchenko et al. 2006; Altuncu et al., 2019; Huang et al., 2020). Similar to the solid-state analogy, where waves are scattered by crystals lattice, here, lattice sums of spherical harmonic functions are adopted to model the scattering from discrete particles, and the center transformation is implemented with the T-matrix approach combined with the translational addition theorems, given the differences between the Schrödinger and Maxwell equations (Popov, 2012; Tan and Tsang, 2019). Consider N spherical particles with radii a ₁, a ₂, a ₃, …, a _N , respectively. The ith particle has permittivity equal to ε i and permeability equal to μ i . The incident field can be expanded into spherical harmonics with the center at origin as follows: E → i n c = ∑ n ∑ m [ a m n R g M → m n ( k r → ) + b m n R g N → m n ( k r → ) ]   (28)

As mentioned above, for an incident EMW with electric field polarized in x-axis and propagating in z-axis, it has the form E → i n = e i ( k z − ω t ) x ^ . We can then obtain the coefficients for the spherical harmonics for the incident waves: a m n = { i n + 1 2 n + 1 2 γ − 1 , n ,             m =   − 1 i n + 1 2 n + 1 2 n ( n + 1 ) γ 1 , n ,         m = 1 0 ,             e l s e (29) b m n = { − a − 1 , n ,                 m =   − 1 a 1 , n ,               m = 1 0 ,                 e l s e (30)

For an arbitrary particle, e.g., the lth particle, we assume the exciting field consisting of two parts: the first part is the original incident wave, and the other is the scattered wave from other particles. The expression for the exciting fields of the lth particle is E → e x ( l ) = ∑ n ∑ m [ p m n ( l ) R g M → m n ( k r r l → ) + q m n ( l ) R g N → m n ( k r r l → ) ] (31) where p m n   and q m n are unknown exciting field coefficients. r r l → equals r → − r l → .

For the lth particle, the exciting field consists of the incidence wave, and the part that comes from the scattered field of other particles can be expressed as E → e x ( l ) ( r → ) =   E → i n c ( r → ) + ∑ j ≠ l E → s ( j ) ( r → ) (32)

The electromagnetic field scattered by the jth particle   can be calculated by multiplying the T-matrix by the total exciting field E → e x ( j ) and replacing the first kind of spherical Bessel function j _n with third kind of spherical Bessel function h _n, which will give rise to scattered outward spherical waves. E → s ( j ) = ∑ n ∑ m [ p m n ( l ) T n ( M ) ( j ) M → m n ( k r r → j ) + q m n ( l ) T n ( N ) ( j ) N → m n ( k r r → j ) ] (33) where T n ( M ) ( j ) and T n ( N ) ( j ) are T-matrices calculated in (Eqs 22–23). Substitute (Eqs 31, 33) into (Eq. 32), and we obtain, ∑ n ∑ m [ p m n ( l ) R g M → m n ( k r r l → ) + q m n ( l ) R g N → m n ( k r r l → ) ] =   ∑ v ∑ μ [ a μ ν R g M → μ ν ( k r → ) + b μ ν R g N → μ ν ( k r → ) ] + ∑ j ≠ l ∑ ν ∑ μ [ p μ ν ( l ) T ν ( M ) ( j ) M → μ ν ( k r r → j ) + q μ ν ( l ) T ν ( N ) ( j ) N → μ ν ( k r r → j ) ] (34)

In order to solve this equation, we need to transform all spherical waves on the right side to be centered at r → l , i.e., transform M → μ ν ( k r r → j ) and R g M → μ ν ( k r r → j ) into linear combinations of M → μ ν ( k r r → l ) , N → μ ν ( k r r → l ) , R g M → μ ν ( k r r → l ) , and R g N → μ ν ( k r r → l ) . The addition theorem for spherical harmonics can be expressed as follows (Cruzan, 1962; Bruning and Yuen, 1971; Tsang and Kong, 1983; Stein, 1997). Let r → =   r → 0 + r → ′ , and we have R g M → m n ( k r → ) = ∑ v ∑ μ { R g A μ ν m n ( k r → 0 ) R g M → μ ν ( k r ′ → ) + R g B μ ν m n ( k r → 0 ) R g N → μ ν ( k r → ′ ) } (35) R g N → m n ( k r → ) =   ∑ v ∑ μ { R g B μ ν m n ( k r → 0 ) R g M → μ ν ( k r → ′ ) + R g A μ ν m n ( k r → 0 ) R g N → μ ν ( k r → ′ ) } (36)

For r 0 > r ′ , M → m n ( k r → ) = ∑ v ∑ μ { A μ ν m n ( k r → 0 ) R g M → μ ν ( k r → ′ ) + R g B μ ν m n ( k r → 0 ) R g N → μ ν ( k r → ′ ) } (37) N → m n ( k r → ) =   ∑ v ∑ μ { B μ ν m n ( k r → 0 ) R g M → μ ν ( k r → ′ ) + R g A μ ν m n ( k r → 0 ) R g N → μ ν ( k r → ′ ) } (38)

For r 0 < r ′ , M → m n ( k r → ) =   ∑ v ∑ μ { A μ ν m n ( k r → 0 ) M → μ ν ( k r → ′ ) + R g B μ ν m n ( k r → 0 ) N → μ ν ( k r → ′ ) } (39) N → m n ( k r → ) = ∑ v ∑ μ   { B μ ν m n ( k r → 0 ) M → μ ν ( k r → ′ ) + R g A μ ν m n ( k r → 0 ) N → μ ν ( k r → ′ ) } (40)

The expressions of A μ ν m n ( k r → ) and B μ ν m n ( k r → ) were presented by Tsang and Kong (Tsang et al., 2001). We substitute these equations of translational additional theorem into (Eq. 34) by replacing terms of r r → j and r → with r → l . Balancing (Eq. 34) yields p m n ( l ) =   ∑ v ∑ μ [ R g A m n μ ν ( k r → l ) a μ ν + R g B m n μ ν ( k r → l ) b μ ν ] + ∑ j ≠ l ∑ v ∑ μ [ A m n μ ν ( k r l r j → ) T ν ( M ) ( j ) p m n ( j ) + B m n μ ν ( k r l r j → ) T ν ( N ) ( j ) q m n ( j ) ] (41) q m n ( l ) = ∑ v ∑ μ [ R g B m n μ ν ( k r → l ) a μ ν + R g A m n μ ν ( k r → l ) b μ ν ] + ∑ j ≠ l ∑ ν ∑ μ [ B m n μ ν ( k r l r j → ) T ν ( M ) ( j ) p m n ( j ) + A m n μ ν ( k r l r j → ) T ν ( N ) ( j ) q m n ( j ) ] (42) for l = 1, 2, 3, …N. These equations can be expressed in matrix form as follows: p → ( l ) = ∑ j ≠ l σ ̿ ( k r l r → j ) T ̿ ( j ) p → ( j ) + R g σ ̿ ( k r → l ) a → i n c (43) q → ( l ) =   ∑ j ≠ l σ ̿ ( k r l r j → ) T ̿ ( j ) q → ( j ) + R g σ ̿ ( k r → l ) b → i n c (44) where a → i n c and b → i n c can be calculated through (Eqs 29, 30), σ ̿ ( k r → ) = [ A ̿ ( k r → ) B ̿ ( k r → ) B ̿ ( k r → ) A ̿ ( k r → ) ] , T ̿ = [ T ̿ ( M ) 0 0 T ̿ ( N ) ] .

Considering the incident field as plane wave without loss of generality, the multiple scattering equations can be written as p → ( l ) =   ∑ j ≠ l σ ̿ ( k r l r j → ) T ̿ ( j ) p → ( j ) + e k → ⋅ r → l a → i n c (45) q → ( l ) =   ∑ j ≠ l σ ̿ ( k r l r → j ) T ̿ ( j ) q → ( j ) + e k → ⋅ r → l b → i n c (46)

The above equations can be solved by iteration. After solving for the exciting coefficients p → and q → , we can get the scattering coefficients and the scattered waves a → s ( j ) = T ̿ ( j ) p → ( j ) (47) b → s ( j ) = T ̿ ( j ) q → ( j ) (48) E → s ( j ) = ∑ n ∑ m [ a m n s ( j ) M → m n ( k r r → j ) + b m n s ( j ) N → m n ( k r r → j ) ]   (49)

After obtaining all scattering coefficients of all particles, the scattered spherical waves centered at different particles should be transformed to be centered at one particular particle. Here, we set the particular center to be r → l . By employing the addition theorems (35)–(40), we can obtain E → s = E → s ( l ) + ∑ j ≠ l ∑ n , m ∑ ν μ [ ( a m n s ( j ) A μ ν m n ( k r l r j → ) + b m n s ( j ) B μ ν m n ( k r l r → j ) ) M → μ ν ( k r r → j ) + ( a m n s ( j ) R g B μ ν m n ( k r l r → j ) + b m n s ( j ) R g A μ ν m n ( k r l r j → ) ) N → μ ν ( k r r → j ) ] = ∑ n , m [ a m n s M → m n ( k r r → j ) + b m n s N → m n ( k r r → j ) ] (50)

Then, the scattering cross section in charged particle clusters considering the multiple scattering effects is given by C s = 2 π k 2 2 ∑ n = 1 ∞ ( 2 n + 1 )   ( ∑ m | a m n s | 2 ∑ m | a m n s | 2 + ∑ m | b m n s | 2 ∑ m | b m n s | 2 ) (51)

For an incident EMW with frequency between 10 and 200 GHz with electric field polarized in the x direction and propagating in the z direction, EMW has the form E → i n = e i ( k z − ω t ) x ^ , where k is the wavenumber. For typical particles with relative permittivity of 3.2 and a radius of 1 mm, the ratio of scattering cross section and extinction cross section of a charged particle to those of uncharged particle are plotted in Figures 1A,B, respectively. C s c a and C e x t are scattering cross section and extinction cross section for a charged particle, respectively. C s c a 0 and C e x t 0 are scattering cross section and extinction cross section for an uncharged particle, respectively. The surface charge density varies from 1 to 100 mC/m², which is consistent with the general level of net charges carried by particles in a sand/dust storm (Zhou et al., 2018). We also consider the extreme case that the surface charge density equals infinity, which will cause the surface conductivity to be infinite. In this case, the scattering properties of the particle will be equivalent to those of an ideal metal particle. It is illustrated in Figure 1A that with the increase in surface charge density, the scattering cross section will be enhanced in the low-frequency range where Rayleigh scattering is applicable. As shown in Figure 1A, as charge density increases, the curve of the scattering cross section gradually approaches the curve of the metal particle and finally reaches a saturated level as charges accumulated. The maximum enhancement of C s c a for an ideal metal particle (i.e., particle with surface charge density to be infinite) is 6.983 times of C s c a 0 for an uncharged particle for frequencies below 10 GHz. On the contrary, the scattering cross section is reduced for charged particles in the higher-frequency range. For ideal metal particle or particles with surface charge density being infinite, C s c a / C s c a 0 reaches a minimum value 0.442 at f = 145.7   GHz . As the frequency gets much higher, that is, the wavelength is far less than the size of the particle, the influence of surface charges to the scattering cross section can be neglected. The reason for these totally different responses to surface charges is that the way free electrons oscillate in response to the applied electromagnetic field is highly dependent on the frequency of the incident waves. It is well known that metals are highly reflective, and electromagnetic waves are not allowed to penetrate through them in the low-frequency visible range, which can be attributed to the surface charges that screens the incident field. However, the free electrons have an intrinsic finite response time, above which they will not be fast enough to screen the penetration from incident fields, and therefore, the material is no longer strongly reflective and turns transparent. Here, the surface plasma frequency for dielectric charged particles, which are related to the electrostatic potential at the surface of the particle and the radius of the charged sphere, is much lower than that of metal since the carrier density for sand/dust particles are much lower compared with metal (Klačka and Kocifaj, 2007). However, a similar “reflective to transparent” transition can be inferred from the scattering and extinction cross sections presented in Figure 1. As illustrated in Figure 1B, the extinction cross section by a charged particle was enhanced. However, in the low-frequency range the enhancement of the extinction of EMWs did not increase monotonically with the increase in the surface charge density as for the scattering cross section. This is because the extinction by the particle is the sum of scattering and absorption. As the surface charge density increases from zero, at the first stage (corresponds to curves of 1 and 10   mC / m 2 ), the imaginary part of the effective complex permittivity, which accounts for energy loss, increases with the charge density, and therefore, the absorption by the particle increases. As the charge density further increases, at this stage (corresponds to curves from 10 to 100   mC / m 2 ), the charged particles were further “metalized” and exhibits highly reflective property, which shields the EMWs from penetrating through the particle. Consequently, the absorption cross section gets smaller as the surface charge density increases. For metallic particles, in the low-frequency range, the absorption is extremely low due to the screening effect to the incident field. This is the reason why the curve for the metallic particle (purple) shows the lowest extinction cross section in Figure 1B. In Figure 2, the electric field distributions (x component of the electric field) are plotted for a neutral particle, and a charged particle with surface current density equals 100 mC/m². As shown in Figures 2A,B at a lower frequency range ( f   =   1   GHz ), the field inside the charged particle is close to zero and is smaller compared with that of the neutral particle, which can be attributed to the screening effect of surface charges. As illustrated in Figures 2C,D, at a higher-frequency regime ( f   =   100   GHz ), for neutral and charged particles, the difference in the penetration of free electrons is not as pronounced as the screening effect in the lower-frequency range, which agrees with the smaller discrepancies in scattering cross sections at higher frequencies as depicted in Figure 1.

As shown in Figure 1B, compared with the scattering cross section, the extinction cross section is greatly enhanced by the surface charge in the low-frequency range, which satisfies the Rayleigh scattering limit. The enhancement increases when the frequency gets lower. For surface charge density equal to 10 mC/m², the enhancement of C e x t reaches 1.5 × 10⁸ times of that of C e x t 0 at 0.1 GHz. If we consider a higher surface charge density, the enhancement will decrease and will be down to 1.7 × 10⁷ for surface charge density equal to 100 mC/m² and only 6.983 times the neutral extinction for ideal metal particles. For higher frequency, the extinction cross section is decreased by surface charges, but the effect is not strong as that of scattering cross section. C e x t / C e x t 0 has a minimum value of 0.441 at f   =   114.7   GHz for ideal metal particles or particles with infinite surface charge density.

Next, the influence of the material of the particles on scattering and extinction properties is investigated. For particles with radius   a   = 1 mm and surface charge density η s   = 10 mC/m², the ratio of scattering cross section and extinction cross section of charged particles to those of the uncharged one are plotted in Figures 3A,B, respectively. The relative permittivities vary from 2 to 5. It can be seen that with the increase in the relative permittivity of particles, the enhancement effect by surface charges on the scattering cross section and extinction cross section is weakened in the low frequency, which satisfies the Rayleigh scattering limit. For a higher-frequency range above 50 GHz, as shown in the inset of Figure 3A, the decrease in the scattering cross section caused by the surface charges is enlarged by higher relative permittivity. In the higher frequency range, the extinction cross section remains unchanged for charged particles compared with uncharged particles.

The influence of the size of the particles on the scattering and extinction properties is also investigated. For a particle with relative permittivity equal to 3.2 and surface charge density equal to 10 mC/m², the ratio of scattering cross section and extinction cross section of charged particles to those of the uncharged one are plotted in Figures 4A,B, respectively. The radius of particles varies from 0.01 to 10 mm. It can be seen that for particles with smaller sizes, the enhancement effect by surface charges on the scattering cross section and extinction cross section is enlarged. This is because the enhancement effect is dependent on the bulk complex refractive index m = ε r , and the imaginary part of m represents the absorption of the particle. Since we mathematically assume the charge is confined to the surface layer of the particle, for larger particles, the surface-to-volume ratio is relatively low; therefore, the scattering/extinction caused by the surface charge is small compared with the distribution from the volume. On the other hand, as the frequency increases and goes beyond the Rayleigh scattering limit, the enhancement is reduced, and the ratio is approximate to unity.

For a multi-particle system, the dependence of the scattering cross section on the number density and surface charge density of the particles is studied. Under the classical assumption of independent scattering, the scattering cross section can be expressed as C s c a = N C s c a s i g , where C s c a s i g is the scattering cross section for a single particle, and N is the number of particles in the system. However, the practical scattering properties are different from the classical assumption due to the multiple scattering effect, and therefore, the scattering cross section is not linearly proportional to the number density of particles. First, we study the influence of the number density of the particles on the scattering cross section in a two-particle system. The number density is directly related to the distance between particles. A higher-number density is equivalent to a closer distance between two particles. Assume the incident EMW propagating in the z direction with frequency between 10 and 200 GHz, and the electric field is polarized in the x direction, which can be expressed as E → i n = e i ( k z − ω t ) x ^ . For particles with a relative permittivity of 3.2 and a radius of 1 mm, the scattering cross section C s c a of a two-particle system is calculated for different distances ( | r → 1 − r → 2 | = 2.2 and 3  mm ). The T-matrix calculation is based on the overall charged model described in the Single charged particle section. The scattering cross section of the two-particle system is calculated under three methods: the T-matrix-based multiple scattering model with the calculated scattering cross section denoted by C s c a , classical independent scattering model with results denoted by C s c a C l a , and CST simulation with results denoted by C s c a C S T . Scattering spectrum with no surface charges and surface charge densities equal to 1 and 10 mC/m² are plotted in Figures 5A,D for different distances between two particles. The enlarged spectrum illustrated in the inset of Figures 5A,D shows that when the surface charge density increases from 0 to 1 and 10 mC/m², the scattering cross section is slightly reduced, and the deviation is amplified with higher surface current densities. This result is consistent with Figure 1A for a single charged particle where the scattering cross section is shown to be suppressed by surface charges in a higher-frequency range above 60 GHz.

There is an obvious difference between the result of multiple scattering method and classical independent scattering method. High accuracy of multiple scattering model can be seen when compared with CST simulation. Using the CST simulation as a benchmark, we calculated the relative difference between the two methods and the simulation to evaluate the accuracy of the multiple scattering model and classical assumption in Figures 5B,E. It can be seen that the multiple scattering model whose result is almost the same as simulation has higher accuracy compared with the classical assumption. As the distance between particles decreases from 3 to 2.2 mm, which is equivalent to increasing the number density, the maximum deviation between the calculated result under classical assumption and the simulation result increase from 2.2 × 10 − 6 to 3.3 × 10 − 6   m 2 , and the mean absolute difference E [ | C s c a C l a − C s c a C S T | ] of the two results increases from 6.6 × 10 − 7 to 9.3 × 10 − 7   m 2 with a 40% increment. The mean relative differences are 1.26 × 10 − 13 and 1.28 × 10 − 13 , respectively, while as shown in Figures 5B,E the mean relative differences between the results calculated by the multiple scattering model and the CST results are below 4 × 10 − 14 for both separations. It reveals that a closer distance between particles or a higher-number density of particles will enhance the multiple scattering effect and make the result deviate further from the result calculated under the classical independent scattering assumption.

Next, to gain deeper insight into the electromagnetic response of clustered charged particles, we investigated the influence of surface charges on the multiple scattering effect. The ratio of the scattering cross section calculated under the multiple scattering model C s c a and classical assumption C s c a C l a is plotted in Figures 5C,F for the two separations. Surface charge densities are set to 0, 1, 10 mC/m², and infinity (the same scattering properties with ideal metal), respectively. There is no observable influence of surface charges to the multiple scattering effect for surface charge density under 10 mC/m². For extreme situations, where we set the surface charge density to infinity, as illustrated in Figure 5C when the inter-particle distance is equal to 1.1 mm, the enhancement of scattering cross section by multiple scattering effect is weakened, and the weakened scattering cross section reaches the minimum when the frequency equals 54 GHz. As particles are getting closer, the deviation of the scattering cross section of metal particles with that of neutral particles is amplified.

We then increase the number of particles and validate the T-matrix-based multiple scattering model in a cluster composed of four particles. The pattern of the particle cluster is shown in the inset in Figures 6A,D. Because the distances between two of the particles in the cluster may be different, the distance R from each center of the particle to the center of the cluster is used to indicate the inter-particle distances. As illustrated in Figures 6A,D the scattering cross sections for four-particle clusters with R equals 1.5 and 2 mm, and the particle radius equals 1 mm are almost two times larger compared with the scattering cross section of the two-particle system since the number of particles is doubled. In the inset in Figures 6A,D the scattering spectrum is demonstrated to be suppressed by surface charges, which is in accordance with the results shown in Figures 1A, 5A,D. In Figures 6B,E when R decreases from 2 to 1.5 mm, the mean absolute difference E [ | C s c a C l a − C s c a C S T | ] between the results under classical independent scattering assumption and the CST simulation results increase from 1.5 × 10 − 6 to 2.6 × 10 − 6   m 2 with an 80% increment, which confirm the conclusion that the multiple scattering effect is amplified as the particles are getting closer to each other. The mean relative differences are 1.45 × 10 − 6 and 3.36 × 10 − 5 , respectively. Compared with the results calculated under the multiple scattering model, as shown in Figures 6B,E the relative differences are several times larger, which validate the accuracy of the multiple scattering model.

Finally, in Figure 7, we depict the ratio of the scattering cross section calculated under the multiple scattering model to the result under the classical assumption C s c a / C s c a C l a as a function of the maximum order Na in the formula of addition theorem and the number of iterations in solving Eqs 45, 46. Here, the frequency is set to 100 GHz, and the particles are set to be uncharged. The convergence is demonstrated to be closely related to the distance between particles. For consistency, different inter-particle distances are indicated by the distance R from the center point of the particle cluster to the center of individual spheres. As illustrated in Figure 7, for two and four particles, the results for four and six iterations are nearly identical, which means four times of iteration is enough for convergence while keeping relatively low computational costs. Therefore, in the former calculations, we adopted four times of iterations to ensure the accuracy of the calculation. As for the maximum order Na applied in the formula of addition theorem, for a two-particle system shown in Figure 7A, the result is convergent when Na is greater than 6 and 8 for R equals 1.1 and 1.5 mm, respectively. For four-particle clusters, the convergence condition is met for Na greater than 8 and 10 for R equals 1.5 and 2.0 mm, respectively. As demonstrated in Figure 7, the slope of the curve is higher, and the convergence condition is met at smaller Na for larger inter-particle distances.

Based on the T-matrix approach, the problem of electromagnetic waves scattered by a single-charged particle can be solved, and the influence of surface charged density and relative permittivity of the particle is investigated. For sizes of particles much smaller than the wavelength of the incident EMWs, the enhancement of extinction cross section by surface charge is enlarged with the higher surface charge density but reduced by higher permittivity of the particle. The enhancement of scattering cross section by surface charge is slightly enlarged by a higher surface charge density and reduced by a higher permittivity of the particle. For a higher-frequency range, the influence of surface charge and relative permittivity on the extinction cross section can be neglected, but the scattering cross section is slightly reduced by higher surface charge and larger relative permittivity of the particle. With the increase in the size of the particles, the enhancement of scattering/extinction cross sections due to the surface charges is smaller.

Next, exploiting the T-matrix approach based on the addition theorem and the iterative method, the multiple scattering problem of clustered charged particles can be solved. The influence of inter-particle distances and surface charge density on the multiple scattering effect is investigated. Scattering cross sections of two-particle systems and four-particle clusters are calculated via the multiple scattering method with varying inter-particle distances. The deviation from the result calculated under the classical independent scattering assumption is obvious when inter-particle distance is smaller than three times that of the particle radius and is further enlarged as particles are getting closer. Surface charge density has little influence on the multiple scattering effect within several hundred mC/m².