Open-source implementation of the discrete-dipole approximation for a scatterer in an absorbing host medium

Theoretical description of light scattering by single particles is a well-developed field, but most of it applies to particles located in vacuum or non-absorbing host medium. Although the case of absorbing host medium has also been discussed in literature, a complete description and unambiguous definition of scattering quantities are still lacking. Similar situation is for simulation methods – some computer codes exist, but their choice is very limited, compared to the case of vacuum. Here we describe the extension of the popular open-source code ADDA to support the absorbing host medium. It is based on the discrete dipole approximation and is, thus, applicable to particles with arbitrary shape and internal structure. We performed test simulations for spheres and compared them with that using the Lorenz-Mie theory. Moreover, we developed a unified description of the energy budget for scattering by a particle in a weakly absorbing host medium, relating all existing local (expressed as volume integrals over scatterer volume) and far-field scattering quantities.


Introduction
Light scattering by small particles is a ubiquitous phenomenon in nature and is widely used in science and technology. Theoretical description of this process is well developed (at least for single particles), allowing one to calculate any scattering quantities. For instance, the Lorenz-Mie theory, applicable to spheres, was developed more than 100 years ago [1]. Later, the T-matrix method [2] and the discretedipole approximation (DDA) [3,4], among other methods, enabled simulation of light scattering by arbitrarily shaped particles, even in inhomogeneous environment (for instance, placed onto a flat substrate).
Surprisingly, the case of a homogeneous but absorbing environment (host medium) is still often omitted both in theoretical descriptions and computational codes. Importantly, this case is relevant for many practical applications, e.g., for particles submerged in water and droplets in oil. For instance, the Lorenz-Mie theory was extended to the case of absorbing host medium only in 1996 [5], although some far-field results were obtained earlier [6]. The Maxwell's equations themselves are perfectly applicable for complex wavenumber , and many light-scattering methods allow straightforward extension (by analytic continuation of all formulae). However, the extension of the existing codes has technical issues, e.g., with using special mathematical functions of complex arguments. Nevertheless, publicly available Lorenz-Mie codes are available for that case [7,8].
The situation is similar for the DDA. The underlying volume-integral equation (VIE) is valid for complex [9] and the first DDA simulations for absorbing host medium has been recently reported [10]. However, the corresponding code is not publicly available. We have independently extended the  [11] to the case of an absorbing host medium. In this paper we report its validation by comparison with the Lorenz-Mie theory.
A more fundamental issue, often overlooked, is the definition of the (measureable) scattering quantities. In the case of vacuum (or non-absorbing host medium), the analysis of the energy budget [12] leads to the classical definition of the intrinsic particle quantities: various cross sections, Mueller and amplitude scattering matrices [13]. They need to be calculated once for a given particle and real h (refractive index of the host medium), after which the measurement of any detector in the far-field zone can be easily simulated, both for a single particle and for a collection of many well-separated particles. Unfortunately, such theoretical description is lacking for the case of complex h , despite existing efforts in this direction [14][15][16]. The main problem is that the transformation of the integration surface for the flow of the Poynting vector from the near-to far-field region leads to the additional term associated with medium absorption. Thus, our second goal is to develop such description for relatively weak absorption, including the definition of the extinction cross section based on the energy-budget arguments.

Volume-integral equation
We consider a non-magnetic medium with the dielectric permittivity 1 ≝ 1 The VIE follows directly from the Maxwell's equation and can be written as [9,12]: where ( ) ≝ [ 2 ( ) − 1 ] ( ) is the polarization density, ( ) is the internal field, inc ( ) is the incident field, is the frequency, 0 is the spherical exclusion volume centered at , and ̅ ( , ′ ) is the Green's dyadic: where = − ′, = | |, and ≝ 2 0 1 is the complex wavenumber: Eq. (1) can be discretized into a set of voxels (dipoles) which constitute a particle of interest, and their polarizations can be determined numerically. This is the basis of the DDA and its implementations, including the code ADDA. Importantly, it remains valid in the case of absorbing medium.

Energy budget and cross sections
The proper definition of power integrals and cross sections (absorption, scattering and extinction) in an absorbing medium is widely discussed in the literature [14][15][16]. In this case the flow of the incident Poynting vector does not vanish and the flow of any Poynting vector depends on the specific position of the integration surface. In particular, the integrals over the particle's surface differ from those taken in the far-field, the difference being the energy absorbed in the medium. Still, the power absorbed inside the particle has a clear physical meaning and can be calculated once ( ) or ( ) is known: where int is the particle volume. Certain scattering and extinction powers can potentially be defined similarly (integrating corresponding components of the Poynting vector over the particle surface), however the physical meaning of such integrals is ambiguous. A more practical approach is to define these quantities through the far-field measurements [16][17][18]. For sufficiently large distances, the scattered wave is the outgoing spherical one: where is the scattering amplitude. While Eq. (5) is valid for any 1 ″ , it is practically relevant only for relatively small 1 ″ so that there exist , such that it is large enough (for far-field limit) but the scattered field is still measureable. Obviously, the damping factor exp(− ″ ) affects the scattered power passing through the far-field sphere with radius , but it can be factored out to define the following normalized scattering power: If we assume that the incident field is a plane wave with amplitude 0 (in the origin) propagating in the direction 0 , one can similarly define the normalized extinction cross-section [16]: However, this definition is meaningful (related to far-field power flow) only if exp(2 ″ ) ≪ ′ , i.e. the absorption is weak and is not too large. Both Eqs. (6) and (7) are convenient because they are independent of and can be further used to assess the energy budget for arbitrary far-field sphere. For a plane incident wave it is convenient to define cross sections dividing the powers by the irradiance of the incident wave in the origin 0 = ′ | 0 | 2 (2 0 ) ⁄ [16,17]. The corresponding cross sections are computed by existing Lorenz-Mie codes, and we have implemented the same in the ADDA code.
The above theory is self-contained but has two important limitations. First, the far-field integration surface is limited to be a sphere (otherwise the normalized powers cannot be used), which is inconvenient for consideration of, e.g., extinction of a plane wave by a slab with many well-separated scatterers (the classical extinction problem for a non-absorbing host medium). Second, sca n requires integration of the scattered field, which is much slower in the framework of the VIE than a single volume integral to be computed for either ext n or abs . We have managed to address both these issues. In short, using the generalized Jones lemma [16] we proved that the integral of the extinction (cross) part of the Poynting vector is the only part, which is relatively independent of the integration surface. Its value is determined by the interference of the incident and scattered fields in the forward direction, hence the surface can be freely deformed as long as it remains smooth and fully inside the intermediate region of , described above, and the forwarddirection part stays intact (at a fixed distance from the scatterer). Next, for weakly absorbing host medium we expressed sca n as an integral over int similarly to that in vacuum [12]. The details of these derivations will be presented at the conference.

Numerical results
We modified the ADDA code to support arbitrary h . The internal code modifications are mostly limited to changing a number of variables and functions from real to complex type. We also had to modify functions calculating the cross sections according to the developed theory. From user's point of view, the main change is the new command-line option -mhost setting h (any complex number). Previously the simulations with ADDA were possible only for real h through the manual scaling of wavelength and particle refractive index (the scattering problem depends only on relative refractive index and wavenumber, except for trivial scaling of some scattering quantities [19]). We verified that the enhanced code obeys this scaling rule (data not shown).
To validate ADDA for absorbing host medium (complex h ), we compared the simulation results with the Lorenz-Mie theory using the publicly available Fortran code [7]. The results for a sphere with size parameter 10 and refractive index 1.5 + 0.1i are shown in Fig. 1. We tested three values of h (1 + 0i, 1 + 0.05i, and 1 + 0.1i) and two discretization levels (15 and 30 dipoles per the wavelength, corresponding to 24 and 48 dipoles per the particle diameter). One can see that ADDA results converge to the exact solution with increasing discretization, both for angle-resolved scattering intensity and cross sections. Moreover, the default discretization of 15 dipoles per the wavelength (for this particle refractive index) already gives reasonable results.  Figure 1 Comparison of ADDA simulations with the Lorenz-Mie theory for a sphere with size parameter 10 and refractive index 1.5 + 0.1i embedded in the medium with refractive index 1 + 0i (A), 1 + 0.05i (B), or 1 + 0.1i (C). Values of cross sections are shown in the insets.

Conclusion
We have implemented the support for the absorbing host medium in the popular open-source code ADDA and validated it by comparison with the exiting Lorenz-Mie code for spheres. The extended code is currently available at the separate branch: https://github.com/Sunmosk/adda-complex_freq, and will be incorporated into future releases of ADDA. Moreover, we extended the energy-budget analysis of the scattering problem to the case of weakly absorbing host medium, relating local (volume integrals over scatterer volume) and far-field scattering quantities.