Silicon nano-sphere

01/2021: updated to pyGDM v1.1+

Comparing pyGDM to Mie theory for a silicon nano-sphere (D=150nm).


from pyGDM2 import structures
from pyGDM2 import materials
from pyGDM2 import fields

from pyGDM2 import core
from pyGDM2 import propagators
from pyGDM2 import tools
from pyGDM2 import linear
from pyGDM2 import visu

import numpy as np
import matplotlib.pyplot as plt

## --- load pre-calculated Mie-data
wl_mie, qext_mie, qsca_mie = np.loadtxt("scat_mie_Si_D150nm.txt").T
qabs_mie = qext_mie - qsca_mie

Simulation setup

## --- Setup incident field
field_generator = fields.planewave
## log-interval spectrum (denser at low lambda):
wavelengths = np.exp(np.linspace(np.log(300), np.log(1000), 50))
kwargs = dict(theta = [0.0])
efield = fields.efield(field_generator, wavelengths=wavelengths,

scale_factor = 1.4
step = 18.75/scale_factor
radius = 4.*scale_factor
geometry = structures.sphere(step, R=radius, mesh='hex', ORIENTATION=2)
material = materials.silicon()
struct = structures.struct(step, geometry, material)

dyads = propagators.DyadsQuasistatic123(n1=1, n2=1)

sim = core.simulation(struct, efield, dyads)

print('(hex) ----- N_dipoles =', len(sim.struct.geometry), end='')
structure initialization - automatic mesh detection: hex
structure initialization - consistency check: 1159/1159 dipoles valid
/home/hans/.local/lib/python3.8/site-packages/pyGDM2/ UserWarning: 3D data. Falling back to XY projection...
  warnings.warn("3D data. Falling back to XY projection...")
(hex) ----- N_dipoles = 1159

Run the simulation

## main simulation
E = core.scatter(sim, method='lu', verbose=True)

## extinction spectrum
field_kwargs = tools.get_possible_field_params_spectra(sim)[0]
wl, spec = tools.calculate_spectrum(sim, field_kwargs, linear.extinct)
a_ext, a_sca, a_abs = spec.T
a_geo = tools.get_geometric_cross_section(sim)
/home/hans/.local/lib/python3.8/site-packages/numba/core/ UserWarning: Numba extension module 'numba_scipy' failed to load due to 'ValueError(No function '__pyx_fuse_0pdtr' found in __pyx_capi__ of 'scipy.special.cython_special')'.
timing for wl=300.00nm - setup: EE 8901.2ms, inv.: 5052.6ms, repropa.: 1758.7ms (1 field configs), tot: 15713.7ms
timing for wl=307.46nm - setup: EE 298.4ms, inv.: 2397.7ms, repropa.: 19.7ms (1 field configs), tot: 2716.5ms
timing for wl=315.11nm - setup: EE 289.7ms, inv.: 1979.9ms, repropa.: 18.3ms (1 field configs), tot: 2288.7ms
timing for wl=322.95nm - setup: EE 227.4ms, inv.: 1898.8ms, repropa.: 13.3ms (1 field configs), tot: 2140.1ms
timing for wl=330.98nm - setup: EE 253.0ms, inv.: 2455.0ms, repropa.: 27.5ms (1 field configs), tot: 2737.1ms
timing for wl=339.22nm - setup: EE 787.0ms, inv.: 4580.3ms, repropa.: 36.4ms (1 field configs), tot: 5405.1ms
timing for wl=347.65nm - setup: EE 769.3ms, inv.: 5904.0ms, repropa.: 31.0ms (1 field configs), tot: 6706.5ms
timing for wl=356.30nm - setup: EE 848.0ms, inv.: 3236.1ms, repropa.: 13.5ms (1 field configs), tot: 4098.4ms
timing for wl=365.17nm - setup: EE 223.5ms, inv.: 1788.3ms, repropa.: 16.5ms (1 field configs), tot: 2029.4ms
timing for wl=374.25nm - setup: EE 221.7ms, inv.: 1969.8ms, repropa.: 16.2ms (1 field configs), tot: 2208.9ms
timing for wl=383.56nm - setup: EE 226.6ms, inv.: 1710.9ms, repropa.: 14.0ms (1 field configs), tot: 1952.5ms
timing for wl=393.10nm - setup: EE 235.6ms, inv.: 1819.8ms, repropa.: 13.7ms (1 field configs), tot: 2069.8ms
timing for wl=402.88nm - setup: EE 222.6ms, inv.: 4402.3ms, repropa.: 34.6ms (1 field configs), tot: 4661.0ms
timing for wl=412.90nm - setup: EE 662.7ms, inv.: 5772.9ms, repropa.: 35.1ms (1 field configs), tot: 6472.0ms
timing for wl=423.17nm - setup: EE 718.9ms, inv.: 4231.8ms, repropa.: 14.4ms (1 field configs), tot: 4966.0ms
timing for wl=433.70nm - setup: EE 231.5ms, inv.: 1680.6ms, repropa.: 12.4ms (1 field configs), tot: 1925.0ms
timing for wl=444.49nm - setup: EE 248.0ms, inv.: 1742.7ms, repropa.: 15.4ms (1 field configs), tot: 2007.0ms
timing for wl=455.54nm - setup: EE 261.7ms, inv.: 1885.0ms, repropa.: 12.8ms (1 field configs), tot: 2160.4ms
timing for wl=466.87nm - setup: EE 282.5ms, inv.: 1791.1ms, repropa.: 13.7ms (1 field configs), tot: 2088.4ms
timing for wl=478.49nm - setup: EE 229.6ms, inv.: 3400.3ms, repropa.: 31.5ms (1 field configs), tot: 3662.8ms
timing for wl=490.39nm - setup: EE 864.8ms, inv.: 6044.9ms, repropa.: 26.7ms (1 field configs), tot: 6938.1ms
timing for wl=502.59nm - setup: EE 880.9ms, inv.: 5472.6ms, repropa.: 12.7ms (1 field configs), tot: 6367.5ms
timing for wl=515.09nm - setup: EE 268.6ms, inv.: 1883.5ms, repropa.: 17.1ms (1 field configs), tot: 2170.5ms
timing for wl=527.90nm - setup: EE 225.2ms, inv.: 1965.9ms, repropa.: 14.5ms (1 field configs), tot: 2206.6ms
timing for wl=541.03nm - setup: EE 211.4ms, inv.: 1949.2ms, repropa.: 14.7ms (1 field configs), tot: 2176.0ms
timing for wl=554.49nm - setup: EE 267.6ms, inv.: 2047.9ms, repropa.: 15.0ms (1 field configs), tot: 2331.1ms
timing for wl=568.29nm - setup: EE 280.8ms, inv.: 4503.2ms, repropa.: 27.1ms (1 field configs), tot: 4812.6ms
timing for wl=582.42nm - setup: EE 804.5ms, inv.: 4616.9ms, repropa.: 32.6ms (1 field configs), tot: 5455.6ms
timing for wl=596.91nm - setup: EE 773.2ms, inv.: 4505.9ms, repropa.: 13.6ms (1 field configs), tot: 5293.7ms
timing for wl=611.76nm - setup: EE 233.4ms, inv.: 1762.7ms, repropa.: 15.5ms (1 field configs), tot: 2012.2ms
timing for wl=626.98nm - setup: EE 228.3ms, inv.: 1785.8ms, repropa.: 14.1ms (1 field configs), tot: 2029.1ms
timing for wl=642.57nm - setup: EE 220.2ms, inv.: 1828.3ms, repropa.: 14.4ms (1 field configs), tot: 2063.6ms
timing for wl=658.56nm - setup: EE 214.5ms, inv.: 1875.8ms, repropa.: 15.3ms (1 field configs), tot: 2106.3ms
timing for wl=674.94nm - setup: EE 224.7ms, inv.: 3088.3ms, repropa.: 34.3ms (1 field configs), tot: 3348.8ms
timing for wl=691.73nm - setup: EE 764.3ms, inv.: 4357.0ms, repropa.: 25.3ms (1 field configs), tot: 5148.2ms
timing for wl=708.93nm - setup: EE 846.5ms, inv.: 4506.1ms, repropa.: 29.8ms (1 field configs), tot: 5384.3ms
timing for wl=726.57nm - setup: EE 630.3ms, inv.: 2681.8ms, repropa.: 13.3ms (1 field configs), tot: 3326.5ms
timing for wl=744.64nm - setup: EE 254.4ms, inv.: 1991.2ms, repropa.: 13.4ms (1 field configs), tot: 2259.8ms
timing for wl=763.17nm - setup: EE 268.8ms, inv.: 1725.8ms, repropa.: 16.0ms (1 field configs), tot: 2011.8ms
timing for wl=782.15nm - setup: EE 221.2ms, inv.: 1710.9ms, repropa.: 15.0ms (1 field configs), tot: 1948.1ms
timing for wl=801.61nm - setup: EE 257.9ms, inv.: 2069.1ms, repropa.: 25.9ms (1 field configs), tot: 2355.1ms
timing for wl=821.55nm - setup: EE 809.9ms, inv.: 4853.8ms, repropa.: 29.9ms (1 field configs), tot: 5695.0ms
timing for wl=841.98nm - setup: EE 685.0ms, inv.: 4615.1ms, repropa.: 25.1ms (1 field configs), tot: 5326.8ms
timing for wl=862.93nm - setup: EE 649.0ms, inv.: 3846.4ms, repropa.: 15.5ms (1 field configs), tot: 4511.6ms
timing for wl=884.39nm - setup: EE 247.5ms, inv.: 1745.4ms, repropa.: 14.2ms (1 field configs), tot: 2007.9ms
timing for wl=906.39nm - setup: EE 238.0ms, inv.: 1908.2ms, repropa.: 13.3ms (1 field configs), tot: 2160.7ms
timing for wl=928.94nm - setup: EE 272.5ms, inv.: 1685.5ms, repropa.: 18.8ms (1 field configs), tot: 1977.9ms
timing for wl=952.05nm - setup: EE 241.4ms, inv.: 1895.3ms, repropa.: 13.2ms (1 field configs), tot: 2151.0ms
timing for wl=975.73nm - setup: EE 261.6ms, inv.: 3039.0ms, repropa.: 25.4ms (1 field configs), tot: 3327.5ms
timing for wl=1000.00nm - setup: EE 756.9ms, inv.: 4784.2ms, repropa.: 25.8ms (1 field configs), tot: 5568.2ms

Plot the spectrum

plt.title("silicon sphere, D=150nm")

## --- Mie
plt.plot(wl_mie, qext_mie, 'b--', dashes=[2,1],label='ext.')
plt.plot(wl_mie, qabs_mie, 'g--', dashes=[2,1],label='abs.')
plt.plot(wl_mie, qsca_mie, 'r--', dashes=[2,1],label='scat.')

## --- pyGDM
plt.scatter(wl, a_ext/a_geo, marker='x', linewidth=1.5, color='b', label='')
plt.scatter(wl, a_abs/a_geo, marker='x', linewidth=1.5, color='g', label='')
plt.scatter(wl, a_sca/a_geo, marker='x', linewidth=1.5, color='r', label='')

## --- for legend only
plt.plot([0], [0], 'k--', dashes=[2,1], label='Mie')
plt.scatter([0], [0], marker='x', linewidth=1.5, color='k', label='pyGDM')
## -- legend
plt.legend(loc='best', fontsize=12)

plt.xlabel("wavelength (nm)")
plt.ylabel("ext. / scat. / abs. efficiency")
plt.xlim( [wl.min(), wl.max()] )
plt.ylim( [0, 10] )


The agreement with Mie theory is ok but not ideal. This can be easily improved by increasing the number of meshpoints (see our paper), which increases of course the simulation time, for which reason we stick to a coarser mesh for this demonstration.