caution with current version of `-int igt_so`

[Maxim Yurkin]

Dear colleagues,

You may have noticed that the master branch of ADDA contains an option to integrate the Green's tensor over the voxel in a semi-analytical manner (IGT(SO), `-int igt_so`). By design, it should have all the benefits of full IGT (`-int igt`), but without any significant extra computational time. This was demonstrated already two years ago at a couple of conference talks:

Yurkin M.A. and Smunev D.A. Analytical integration of the Green’s tensor in the discrete dipole approximation, Bremen Zoom Workshop on Light Scattering 2023, 20–21 March 2023, Online, pp. 67–71.  (paper, slides)
Yurkin M.A. and Shershnev D.V. Robust DDA simulations for exciting sources inside a particle, The 20th Electromagnetic and Light Scattering Conference, 15–19 May 2023, Almunecar, Spain, p. O04.  (abstract, slides, video)

However, we can't yet claim IGT(SO) to be universally better that other formulations, since the involved approximations can be significant in some cases. SO (second order) in its name means that the interaction term G_ij (an integral of the Green's tensor over the volume of cuboid voxel) is computed with correct consideration of all terms of order (kd)^2, where k is the vacuum wave number and d is the dipole (voxel) size (largest dimension of the cuboid). That implies that error of G_ij is O[(kd)^4]. When summed over all cuboids this leads to (integration) error O[(kd)^2] or smaller. By its order, this is not larger than other sources of error in the DDA, which allows this formulation to reach quadratic convergence of the DDA in some cases (e.g., cubes with moderate refractive index). However, the prefactor of this scaling is not necessarily small, which means that for fixed standard kd (e.g., the rule-of-thumb leads to kd=pi/(5*|m|)) this term maybe larger than even linear term. Such errors has already been noticed in the above presentations (slides 23,24 and 13,14, respectively) for source placed inside the particle, but they may also become relevant for more mundane simulations, like far-field scattering by large optically soft particles. 

The main conclusions:

1) Polarizability calculation IGT(SO), `-pol igt_so`, is not affected by this issue. We are not aware of any more accurate way to calculate it, and the resulting error (when summed over all dipoles) has never been significant, since here only the diagonal of the interaction matrix is affected. In other words, this polarizability formulation is used for both SO and full IGT in the interaction term. Thus, when talking about IGT(SO) vs. IGT in the following, we only mean the interaction term.

2) The IGT(SO) discussed here (in the master branch) is significantly different from table-based implementation in v.1.4.0. The latter has limited applicability in terms of particle sizes and, thus, is not recommended for any use.

3) IGT(SO) is almost perfect for small kd, e.g., when simulating nanoparticles.

4) If you have any doubt about `-int igt_so`, compare the result for one or a few test cases with `-int igt`. The latter will require large initialization time (comparable to several hundreds of iterations), but should work reliably for all range of parameters. If you see significant differences, try `-int igt 3` or `-int igt 5` (truncated IGT). This may lead to almost all benefits of full IGT, but at time comparable to IGT(SO). In principle, truncated IGT does introduce some error, but it is usually noticeable only in convergence studies, when the grid is refined for the fixed problem, but not for practical simulations with moderate kd. In this sense, truncated IGT is the reverse of IGT(SO) in terms of pros and cons for various applications. If one is not that good, the other will most probably work good. Another modern formulation, FCD, can also be competitive, but we do not have general guidelines for that. 

5) If you encounter any cases, where IGT(SO) leads to significant errors (compared to full IGT), please let us know, e.g., by replying to this message. Meanwhile, we have a few ideas for further improving IGT(SO) and will provide an update, once any of them is implemented.

Maxim.