Processing

linear

linear contains evaluation functions for linear effects.

extinction/absorption/scattering cross-sections

pyGDM2.linear.extinct(sim, field_index)

Extinction, scattering and absorption corss-sections

Calculates extinction, scattering and absorption crosssections for each wavelength of the GDM simulation

Parameters:
Returns:

  • extinct (float) – extinction cross-section
  • scatter (float) – scattering cross-section
  • absorpt (float) – apsorption cross-section

Notes

For the calculation of the cross-sections from the complex nearfield, see e.g.: Draine, B. T. & Flatau, P. J. Discrete-dipole approximation for scattering calculations. Journal of the Optical Society of America, A 11, 1491 (1994).

nearfield around structure

pyGDM2.linear.nearfield(sim, field_index, r_probe)

Nearfield distribution in proximity of nanostructre

For a given incident field, calculate the electro-magnetic field in the proximity of the nanostructure (positions defined by MAP).

Parameters:
  • sim (core.simulation) – simulation description
  • field_index (int) – index of evaluated self-consistent field to use for calculation. Can be obtained for specific parameter-set using tools.get_closest_field_index()
  • r_probe (tuple (x,y,z) or list of 3-lists/-tuples) – (list of) coordinate(s) to evaluate nearfield on. Format: tuple (x,y,z) or list of 3 lists: [Xmap, Ymap, Zmap] (the latter can be generated e.g. using tools.generate_NF_map())
Returns:

  • 4 lists of 6-tuples, complex
    • scattered Efield
    • total Efield (inlcuding fundamental field)
    • scattered Bfield
    • total Bfield (inlcuding fundamental field)
  • the tuples are of shape (X,Y,Z, Ax,Ay,Az) with Ai the corresponding
  • complex field component

Notes

For details of the calculation of the scattered field outside the nano-object using the self-consistent field inside the particle, see e.g.: Girard, C. Near fields in nanostructures. Reports on Progress in Physics 68, 1883–1933 (2005).

farfield patterns

pyGDM2.linear.farfield(sim, field_index, r=10000.0, tetamin=0, tetamax=1.5707963267948966, Nteta=10, Nphi=36, polarizerangle='none', return_value='map')

spatially resolved and polarization-filtered far-field scattering

For a given incident field, calculate the electro-magnetic field (E-component) in the far-field around the nanostructure (on a sphere of radius r).

Parameters:
  • sim (core.simulation) – simulation description
  • field_index (int) – index of evaluated self-consistent field to use for calculation. Can be obtained for specific parameter-set using tools.get_closest_field_index()
  • r (float, default: 10000.) – radius of integration sphere (distance to coordinate origin in nm)
  • tetamax (tetamin,) – minimum and maximum polar angle for semi-sphere in radians (calculate from tetamin to tetamax)
  • Nphi (Nteta,) – number of polar and azimuthal angles on sphere to calculate, (angle ranges: polar = 0 - tetamax, azimuthal: 0 - 2*Pi)
  • polarizerangle (float or 'none', default: 'none') – optional polarization filter angle **in degrees**(!). If ‘none’ (default), the total field-intensity is calculated (= no polarization filter)
  • return_value (str, default: 'map') –
    Values to be returned. Either ‘map’ (default) or ‘integrated’.
    • ”map” : (default) return spatially resolved farfield intensity at each spherical coordinate (5 lists)
    • ”int_Es” : return the integrated scattered field (as float)
    • ”int_E0” : return the integrated fundamental field (as float)
    • ”int_Etot” : return the integrated total field (as float)
Returns:

  • return_value == “map” (5 arrays of shape (Nteta, Nphi) :) –
    • tetalist : teta angles
    • philist : phi angles
    • I_sc : intensity of scattered field,
    • I_tot : intensity of total field (I_tot=|E_sc+E_0|^2),
    • I0 : intensity of incident field
  • return_typ == “int_XX” (float) – integrated total intensity over specified solid angle

Notes

For details of the asymptotic (non-retarded) far-field propagators for a dipole above a substrate, see e.g.: Colas des Francs, G. & Girard, C. & Dereux, A. Theory of near-field optical imaging with a single molecule as light source. The Journal of Chemical Physics 117, 4659–4666 (2002)

heat deposition inside structure

pyGDM2.linear.heat(sim, field_index, power_scaling_e0=1.0, return_value='total', return_units='nw')

calculate the total induced heat in the nanostructure

Parameters:
  • sim (core.simulation) – simulation description
  • field_index (int) – index of evaluated self-consistent field to use for calculation. Can be obtained for specific parameter-set using tools.get_closest_field_index()
  • power_scaling_e0 (float, default: 1) – incident laser power scaling. power_scaling_e0 = 1 corresponds to 1 mW per micron^2. See [1].
  • return_value (str, default: 'total') –

    Values to be returned. Either ‘total’ (default) or ‘structure’.

    • ”total” : return the total deposited heat in nW (float)
    • ”structure” : return spatially resolved deposited heat at each meshpoint in nW (list of tuples [x,y,z,q])
  • return_units (str, default: "nw") – units of returned heat values, either “nw” or “uw” (nano or micro Watts)
Returns:

Q

  • `return_value`=”structure” : (return float) total deposited heat in nanowatt (optionally in microwatt).
  • `return_value`=”structure” : (return list of tuples [x,y,z,q]) The returned quantity q is the total deposited heat at mesh-cell position [x,y,z] in nanowatt. To get the heating power-density, please divide by the mesh-cell volume.

Return type:

float or list of tuples [x,y,z,q]

Notes

For details on heat/temperature calculations and raster-scan simulations, see:

[1] Baffou, G., Quidant, R. & Girard, C.: Heat generation in plasmonic nanostructures: Influence of morphology Applied Physics Letters 94, 153109 (2009)

[2] Teulle, A. et al.: Scanning optical microscopy modeling in nanoplasmonics Journal of the Optical Society of America B 29, 2431 (2012).

temperature rise close to structure

pyGDM2.linear.temperature(sim, field_index, r_probe, kappa_env=0.6, kappa_subst=None, incident_power=1.0)

calculate the temperature rise at locations outside the nano-particle

Calculate the temperature increase close to a optically excited nanostructure using the approach described in [2] and [3] (Ref. [3] introduces a further correction term for particles lying on a substrate)

Parameters:
  • sim (core.simulation) – simulation description
  • field_index (int) – index of evaluated self-consistent field to use for calculation. Can be obtained for specific parameter-set using tools.get_closest_field_index()
  • r_probe (tuple (x,y,z) or list of 3-lists/-tuples) – (list of) coordinate(s) to evaluate nearfield on. Format: tuple (x,y,z) or list of 3 lists: [Xmap, Ymap, Zmap] (the latter can be generated e.g. using tools.generate_NF_map())
  • kappa_env (float, default: 0.6) – heat conductivity of environment. default: kappa_env = 0.6 (water). (air: kappa=0.024, for more material values see e.g. [4]). In W/mK.
  • kappa_subst (float, default: None) – heat conductivity of substrate. default: None –> same as substrate. Using the mirror-image technique described in [3]. (glass: kappa=0.8)
  • incident_power (float, default: 1.0) – incident laser power density in mW per micron^2. See also [1].
Returns:

  • if evaluating at a single position, D_T (float) – temperature increase in Kelvin at r_probe
  • if evaluating at a list of positions, list of tuples [x, y, z, D_T] – where D_T is the temperature increase in Kelvin at (x,y,z), which are the positions defined by r_probe

Notes

For details on heat/temperature calculations and raster-scan simulations, see:

[1] Baffou, G., Quidant, R. & Girard, C.: Heat generation in plasmonic nanostructures: Influence of morphology Applied Physics Letters 94, 153109 (2009)

[2] Baffou, G., Quidant, R. & Girard, C.: Thermoplasmonics modeling: A Green’s function approach Phys. Rev. B 82, 165424 (2010)

[3] Teulle, A. et al.: Scanning optical microscopy modeling in nanoplasmonics Journal of the Optical Society of America B 29, 2431 (2012).

For the thermal conductivity of common materials, see:

[4] Hugh D Young, Francis Weston Sears: University Physics, chapter 15, 8th. edition: Addison-Wesley, 1992 (see also e.g.: http://hyperphysics.phy-astr.gsu.edu/hbase/Tables/thrcn.html)

dipole decay rate evaluation

pyGDM2.linear.decay_eval(sim, SBB, mx, my, mz, verbose=False)

evaluate decay rate of electric or magnetic dipole transition

Evaluate the decay rate modification of a dipole with complex amplitude (mx,my,mz) using pre-calculated tensors (core.decay_rate()).

Parameters:
  • sim (core.simulation) – simulation description
  • SBB (int or list of lists) – index of wavelength in sim or tensor-list as returned by core.decay_rate()
  • mx,my,mz (float) – x/y/z amplitude of dipole transition vector
  • verbose (bool default=False) – print some runtime info
Returns:

gamma_map

Each element consists of:
  • x,y,z: coordinates of evaluation position
  • gamma: normalized decay-rate (real) gamma / gamma_0 at each map-position

Return type:

list of lists [x,y,z, gamma]

Notes

  • For detailed information about the underlying formalism, see: Wiecha, P. R., Girard, C., Cuche, A., Paillard, V. & Arbouet, A. Decay Rate of Magnetic Dipoles near Non-magnetic Nanostructures. Phys. Rev. B 97, 085411 (2018).
  • Requires scipy v0.17.0 or later to work also inside the volume of a nanostructure

nonlinear

nonlinear contains evaluation functions for non-linear effects.

TPL and SP-LDOS

pyGDM2.nonlinear.tpl_ldos(sim, field_index, nonlin_order=2.0, beta=100000.0, verbose=False)

calculate the TPL for nano-object under given illumination condition

Calculate the two-photon photolouminescence (TPL). Higher order nonlinear photoluminescence can also be calculated by changing the parameter nonlin_order.

  • Can be used to simulate TPL signals using `nonlin_order`=2 (default)
  • Using `nonlin_order`=1 together with an unphysically tight focused beam –> calculate the LDOS. WARNING: The focus waist must not be smaller than some times the stepsize!
  • Might be used for Raman intensities using `nonlin_order`=1.
Parameters:
  • sim (core.simulation) – simulation description
  • field_index (int) – index of evaluated self-consistent field to use for calculation. Can be obtained for specific parameter-set using tools.get_closest_field_index()
  • nonlin_order (float, default: 2.0) – order of nonlinear response. (I_npl ~ | |E|^2 | ^ nonlin_order).
  • beta (float, default: 1E5) – incident laser power scaling. Use 1E5 for correct LDOS values.
  • verbose (bool, default: False) – Enable some info printing
Returns:

I_tpl – TPL intensity in farfield from nano-object under given illumination

Return type:

float

Notes

For details on TPL/LDOS calculation via focused beam rasterscan simulations, see: Viarbitskaya, S. et al. Tailoring and imaging the plasmonic local density of states in crystalline nanoprisms, Nat. Mater. 12, 426–432 (2013)