Core

core

core contains the main simulation object and the scatter routine to calculate the self-consistent fields inside the nano-structure.

simulation description class

class pyGDM2.core.simulation(struct, efield, dtype='f')

Bases: object

Main GDM simulation container class

Defines a linear GDM simulation. Contains information on:
  • struct : structures.struct:
    • the geometry of the nano structure
    • its dielectric function
    • its environment
  • efield : fields.efield
    • the incident field and the wavelenghts
    • possibly further incident field configuration
Parameters:
  • struct (structures.struct) – structure object
  • efield (fields.efield) – fundamental field
  • dtype ((optional) str, default: ‘f’) – precision of simulation
__init__(struct, efield, dtype='f')

Initialization

scattered field inside particle

pyGDM2.core.scatter(sim, method='lu', cgKwargs={}, cg_recycle_pc=True, pc_method='ilu', pc_recalc_thresh=0.5, multithreaded=True, Nprocesses=None, verbose=True)

Perform a linear scattering GDM simulation

Calculate the electric field distribution in a nano-structure. Multithreaded parallel execution only using the “dyson” method for inversion of the problem.

Note: 2D in SI units, 3D in cgs(gaussian) units

Parameters:
  • sim (simulation) – simulation description
  • method (string, default: "lu") –
    inversion method. One of [“lu”, “numpyinv”, “scipyinv”, “pinv2”, “superlu”, “cg”, “pycg”, “dyson”]
    • ”lu” LU-decomposition (scipy.linalg.lu_factor)
    • ”numpyinv” numpy inversion (np.linalg.inv, if numpy compiled accordingly: LAPACK’s dgesv)
    • ”scipyinv” scipy default inversion (scipy.linalg.inv)
    • ”pinv2” is usually slower (req. scipy.linalg.pinv2)
    • ”superlu” often has a <N^3 efficieny but is very expensive in memory (scipy.sparse.linalg.splu)
    • ”cg” and “pycg” use iterative conjugate gradient solving, scale with N^2. Not appropriate for many field configs per wavelength! (req. scipy or pyamg)
    • ”dyson” uses a dyson sequence. Only solver that does not depend on external libraries
  • cgKwargs (dict, default: {}) – keywords passed to conjugate-gradient method, if used
  • cg_recycle_pc (bool, default: True) – Use same Preconditioner for all wavelengths (default: True) Remark: For large dispersion / large wavelength differences cg_recycle_pc=True might slow down convergence significantly. Buy mostly this will speed up the Calcuation.
  • pc_method (str, default 'ilu') – Preconditioning method. ‘ilu’: scipy incomplete LU, ‘lu’: scipy complete LU (exact, but speed-imporvements for nearby wavelengths), ‘amg’ pyAMG amg, ‘none’: No preconditioning
  • PCrecalcThresh (float, default 0.5) – Preconditioning recalculation time threshold. If speedup of PC-recycling falls below this threshold, recalculate PC on next wavelength.
  • multithreaded (bool, default: True) – Whether to use multi-threading for multi-incident field processing (at each wavelength). Requires “pathos”.
  • Nprocesses (int, default: None) – if “multithreaded”==True, specifies number of processes. if “None”: automatic
  • verbose (bool, default False) – Print timing info to stdout
Returns:

E – list consists of incident field config and scattered electric field in structure. Each element of the list is composed as follows:

elements :
[dict: kwargs of incident efield,

list: field as list of tuples [[Ex1,Ey1,Ez1], [Ex2,Ey2,Ez2], …] ]

Return type:

list of lists

Notes

pyGDM2.core.scatter_mpi(sim, verbose=False, multithreaded=False, **kwargs)

MPI wrapper to scatter() for embarrassingly parallel calculation of spectra

requires: mpi4py

run with “mpirun -n X python scriptname.py”, where X is the number of parallel processes (ideally an integer divisor of the number of wavelengths)

Parameters:
  • sim (simulation) – simulation description
  • multithreaded (bool, default: False) – Whether to use multi-threading for multi-incident field processing (at each wavelength). Requires “pathos”. warning: use with caution together with MPI!
  • **kwargs – all kwargs are passed to scatter()

Notes

  • On single machines it is usually easier to install scipy compiled with parallel BLAS (parallel LU / CG routines). Usually the parallel BLAS is will be already installed automatically. Try to not use both parallelisation techniques simultaneously unless properly configured for instance via a batch script (e.g. SLURM). Overloading the CPUs will usually result in decreased calculation speed.
  • see scatter() for documentation
  • by default, will pass multithreaded = False to scatter()

decay rate of dipole transition

pyGDM2.core.decay_rate(sim, method='lu', verbose=False)

Calculate the decay rate of a dipole emitter

Calculate the propagator necessary to derive the change of decay rate of an electric or magnetic dipolar emitter in the vicinity of a photonic nanostructure. The result of this routine can be used to calculate the actual decay rate for a dipole transition of any orientation / amplitude using linear.decay_eval().

Parameters:
  • sim (simulation) – simulation description
  • method (string, default: "lu") –
    inversion method. One of [“lu”, “numpyinv”, “scipyinv”, “pinv2”, “superlu”, “dyson”]
    • ”lu” LU-decomposition (scipy.linalg.lu_factor)
    • ”numpyinv” numpy inversion (np.linalg.inv, if numpy compiled accordingly: LAPACK’s dgesv)
    • ”scipyinv” scipy default inversion (scipy.linalg.inv)
    • ”pinv2” is usually slower (req. scipy.linalg.pinv2)
    • ”superlu” often has a <N^3 efficieny but is very expensive in memory (scipy.sparse.linalg.splu)
    • ”dyson” uses a dyson sequence. Only solver that does not depend on external libraries
  • verbose (bool default=False) – print some info
Returns:

list of listsS_P again is a list of the 3x3 tensors S_P^EE (or S_P^BB) for each map-position (EE: electric dipole, BB: magnetic dipole)

Return type:

each element of the main list contains [wavelength, S_P]

Notes

sim must contain either an electric or magnetic dipole emitter as fundamental field. Its orientation and amplitude will however be ignored. These need to be specified later, when the actual evaluation is done using linear.decay_eval()

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. arXiv:1707.07006 (2017)

Other functions

pyGDM2.core.get_side_by_side(sim, wavelength, invertible=True, times_alpha=True, as_csc=True)

Build side-by-side matrix CA0 for the structure at given wavelength

Parameters:
  • sim (simulation) – simulation description
  • wavelength (float) – Wavelength at which to calculate susceptibility matrix (in nm)
  • invertible (bool, default: True) – return invertible matrix (I-CA0)
  • times_alpha (bool, default: True) – multiply by alpha. If True: CA0 = alpha.G; else CA0 = G
  • as_csc (bool, default: True) – return as csc sparse matrix format
Returns:

CA0 – return the side-by-side matrix CA0 (“S”). Inverse of CAM0=(I-CA0) is generalized propagator

Return type:

array-like

Notes

For the analytical expression of the asymptotic propagators with substrate, see e.g.: Girard, C. Near fields in nanostructures. Reports on Progress in Physics 68, 1883–1933 (2005).

pyGDM2.core.get_general_propagator(CAM0, method='lu', sim=None, wavelength=None, return_linear_operator=False, return_susceptibility=False)

Invert Matrix CAM0

Parameters:
  • CAM0 (array-like) – Matrix to invert (I-CA0). As returned by get_side_by_side() (“dyson” calculates CAM0 itself, may be set to None in that case)
  • method (string, default: "lu") –
    inversion method. One of [“lu”, “numpyinv”, “scipyinv”, “pinv2”, “superlu”, “dyson”]
    • ”lu” LU-decomposition (scipy.linalg.lu_factor)
    • ”numpyinv” numpy inversion (np.linalg.inv, if numpy compiled accordingly: LAPACK’s dgesv)
    • ”scipyinv” scipy default inversion (scipy.linalg.inv)
    • ”pinv2” is usually slower (scipy.linalg.pinv2)
    • ”superlu” often has a <N^3 efficieny but is very expensive in memory (scipy.spares.linalg.splu)
    • ”dyson” uses a dyson sequence. Does not depend on external libraries
  • simDict (dict) – simulation description, generated from ‘genSimDict’
  • Ilambda (int, default: None) – Wavelength, only required for method “dyson”
  • retLinearOperator (bool, default: True) – return K as LinearOperator Class
  • return_susceptibility (bool, default: False) – return the field susceptibility S inside the structure instead of the generalized propagator (K = 1 + chi*S)
Returns:

- K

Return type:

Inverse of CAM0. Generalized Propagator

Notes

For details on the concept of the generalized propagator, see e.g.: Martin, O. J. F. & Girard, C. & Dereux, A. Generalized Field Propagator for Electromagnetic Scattering and Light Confinement. Phys. Rev. Lett. 74, 526–529 (1995).

pyGDM2.core.get_efield_by_cg(CAM0, E0, pcTol=0.002, fill_factor=15, cgsTol=0.01, maxiter=200, PC=None, method='cg', pc_method='ilu', verbose=False)

Invert Matrix CAM0

Parameters:
  • CAM0 (array-like) – Matrix to invert (I-CA0). As returned by get_side_by_side.
  • E0 (array-like) – Incident electric field, 3N super-vector (you may use _fieldListToSuperVector())
  • cgsTol (pcTol,) – tolerances for preconditioning / conjugate gradients pcTol > 1: Don’t use preconditioning
  • maxiter (int, default: 500) – Max. iterations for CG
  • fill_factor (int, default: 15) – Preconditioning Parameter
  • PC (array-like, default: None) – Optional Preconditer Matrix. If not None, overrides preconditioning.
  • method (string, default: 'cg') – iterative solver to use. ‘cg’: scipy stab-cg; ‘pycg’ pyAMG implementation of stab-cg (threadsafe)
  • pc_method (string, default: 'ilu') –
    Preconditioning method.
    • ’ilu’: scipy incomplete LU,
    • ’lu’: scipy complete LU (exact, but speed-improvements for close wavelengths when recycling),
    • ’amg’ pyAMG amg,
    • ’none’: No preconditioning
Returns:

  • E (array-like)
  • M (array-like, STATUS : int)
  • - E (Solution of CAM0*E = E0 (3N supervector))
  • - M (Used Preconditioner Matrix (for possible re-use))
  • - STATUS (Status of conj. gradient. (0: converged, !=0: not converged))