open source FDTD solver with GPU support
FDTD is well suited for purposes of optical force calculation in different geometries used in light trapping of nanoparticles (e.g. in optical tweezers). Here we show a simple example of calculation of scattering force acting on dielectric particle under a plane wave illumination and its dependence on particle permittivity. The result is compared to Mie solution (Mie data courtesy of Martin Siler, ISI ASCR).
Image on the right shows a scheme of the computational volume used for the simulation. A normally incident plane wave with chosen wavelength (here 1 um) was entered into computational volume using Total/Scattered field (TSF) approach. A dielectric sphere was placed in the computational volume and a box for force calculation was placed over it. Whole computational volume was terminated by an absorbing boundary condition. To calculate dependence of the force on spherical particle permittivity or size, we performed single computation for every permittivity or size. To obtain such a dependence therefore needs significant amount of time. On the other hand, also particles of different shape than sphere could be introduced easily using the same approach which is a benefit over some simpler computational methods. If dependency on wavelength is calculated, it could be also in principle speeded up by broadband excitation.
A result of a calculation of force acting on a particle in the illumination direction for relative permittivity in the range of 2-4 is shown below, compared to the same quantities obtained using Mie-Lorenz theory for particle radius of 900 nm (top) and 1800 nm (bottom). For larger particles we can see many peaks related to different resonance phenomena inside particles.
As we can see, results are similar but not exactly the same. There are several effects affecting the FDTD performance in this case. First of all, it is discretisation and related staircasing. Result is very sensitive to real particle volume, which does not match exactly the requested value after discretisation. For example, resulting peaks for nominal sphere radius of 900 nm were slightly shifted. This was compensated by using radius value of only 896 nm (which is less than half percent deviation from the nominal value) and we can see that this tuning could be done even better. Even more important is the effect of number of computation steps. There are different resonance modes developing in larger particles, providing sharp peaks in the force dependence and obviously it takes some time to let these modes to develop in FDTD. A comparison of the force dependence for four different numbers of steps is shown below for a single resonance peak in larger particle dependency. Note that the particle volume was still slightly wrong so the peak position was shifted by 0.007 in permittivity for comparison purposes. We can clearly see that for longer computation we get sharper resonance peak.
As discussed above, the benefit of using FDTD is in variability of the parameters. The same setup can be used to calculate the force on a non-ideal particle, either dielectric or metallic as shown below (a simulation snapshot). Here the particle was created using Blender modeler and imported into GSvit via tetrahedral mesh created by Tetgen. In the image below (A) denotes the Blender model, (B) is a cross-section as seen in GSvit, (C) is a cross-section snapshot of the simulation assuming that particle is dielectric and (D) similar snapshot for simulation assuming that the particle is gold.