API

Smuthi is a Python package with the following modules and sub-packages.

smuthi.coordinates module

smuthi.coordinates.angular_frequency(vacuum_wavelength)

Angular frequency \(\omega = 2\pi c / \lambda\)

Parameters:vacuum_wavelength (float) – Vacuum wavelength in length unit
Returns:Angular frequency in the units of c=1 (time units=length units). This is at the same time the vacuum wavenumber.
smuthi.coordinates.complex_contour(vacuum_wavelength, neff_waypoints, neff_resolution)
smuthi.coordinates.k_z(k_parallel=None, n_effective=None, k=None, omega=None, vacuum_wavelength=None, refractive_index=None)

z-component \(k_z=\sqrt{k^2-\kappa^2}\) of the wavevector. The branch cut is defined such that the imaginary part is not negative. Not all of the arguments need to be specified.

Parameters:
  • k_parallel (numpy ndarray) – In-plane wavenumber \(\kappa\) (inverse length)
  • n_effective (numpy ndarray) – Effective refractive index \(n_\mathrm{eff}\)
  • k (float) – Wavenumber (inverse length)
  • omega (float) – Angular frequency \(\omega\) or vacuum wavenumber (inverse length, c=1)
  • vacuum_wavelength (float) – Vacuum wavelength \(\lambda\) (length)
  • refractive_index (complex) – Refractive index \(n_i\) of material
Returns:

z-component \(k_z\) of wavenumber with non-negative imaginary part (inverse length)

smuthi.coordinates.set_default_k_parallel(vacuum_wavelength, neff_waypoints=None, neff_resolution=0.01, neff_max=None, neff_imag=0.05)

smuthi.field_expansion module

Classes and functions to manage the expansion of the electric field in plane wave and spherical wave basis sets.

class smuthi.field_expansion.FarField(polar_angles='default', azimuthal_angles='default', signal_type='intensity')

Represent the far field intensity of an electromagnetic field.

\[P = \sum_{j=1}^2 \iint \mathrm{d}^2 \Omega \, I_{\Omega,j}(\beta, \alpha),\]

where \(P\) is the radiative power, \(j\) indicates the polarization and \(\mathrm{d}^2 \Omega = \mathrm{d}\alpha \sin\beta \mathrm{d}\beta\) denotes the infinitesimal solid angle.

Parameters:
  • polar_angles (numpy.ndarray) – Polar angles (default: from 0 to 180 degree in steps of 1 degree)
  • azimuthal_angles (numpy.ndarray) – Azimuthal angles (default: from 0 to 360 degree in steps of 1 degree)
  • signal_type (str) – Type of the signal (e.g., ‘intensity’ for power flux far fields).
alpha_grid()
Returns:Meshgrid with \(\alpha\) values.
append(other)

Combine two FarField objects with disjoint angular ranges. The other far field is appended to this one.

Parameters:other (FarField) – far field to append to this one.
azimuthal_integral()

Far field as a function of polar angle only.

\[P = \sum_{j=1}^2 \int \mathrm{d} \beta \, I_{\beta,j}(\beta),\]

with

\[I_{\beta,j}(\beta) = \int \mathrm{d} \alpha \, \sin\beta I_j(\beta, \alpha),\]
Returns:\(I_{\beta,j}(\beta)\) as numpy ndarray. First index is polarization, second is polar angle.
beta_grid()
Returns:Meshgrid with \(\beta\) values.
bottom()

Split far field into top and bottom part.

Returns:FarField object with only the intensity for bottom hemisphere (\(\beta\geq\pi/2\))
export(output_directory='.', tag='far_field')

Export far field information to text file in ASCII format.

Parameters:
  • output_directory (str) – Path to folder where to store data.
  • tag (str) – Keyword to use in the naming of data files, allowing to assign them to this object.
integral()

Integrate intensity to obtain total power \(P\).

Returns:\(P_j\) as numpy 1D-array with length 2, the index referring to polarization.
top()

Split far field into top and bottom part.

Returns:FarField object with only the intensity for top hemisphere (\(\beta\leq\pi/2\))
class smuthi.field_expansion.FieldExpansion

Base class for field expansions.

diverging(x, y, z)

Test if points are in domain where expansion could diverge. Virtual method to be overwritten in child classes.

Parameters:
  • x (numpy.ndarray) – x-coordinates of query points
  • y (numpy.ndarray) – y-coordinates of query points
  • z (numpy.ndarray) – z-coordinates of query points
Returns:

numpy.ndarray of bool datatype indicating if points are inside divergence domain.

electric_field(x, y, z)

Evaluate electric field. Virtual method to be overwritten in child classes.

Parameters:
  • x (numpy.ndarray) – x-coordinates of query points
  • y (numpy.ndarray) – y-coordinates of query points
  • z (numpy.ndarray) – z-coordinates of query points
Returns:

Tuple of (E_x, E_y, E_z) numpy.ndarray objects with the Cartesian coordinates of complex electric field.

valid(x, y, z)

Test if points are in definition range of the expansion. Virtual method to be overwritten in child classes.

Parameters:
  • x (numpy.ndarray) – x-coordinates of query points
  • y (numpy.ndarray) – y-coordinates of query points
  • z (numpy.ndarray) – z-coordinates of query points
Returns:

numpy.ndarray of bool datatype indicating if points are inside definition domain.

class smuthi.field_expansion.PiecewiseFieldExpansion

Manage a field that is expanded in different ways for different domains, i.e., an expansion of the kind

\[\mathbf{E}(\mathbf{r}) = \sum_{i} \mathbf{E}_i(\mathbf{r}),\]

where

\[\begin{split}\mathbf{E}_i(\mathbf{r}) = \begin{cases} \tilde{\mathbf{E}}_i(\mathbf{r}) & \text{ if }\mathbf{r}\in D_i \\ 0 & \text{ else} \end{cases}\end{split}\]

and \(\tilde{\mathbf{E_i}}(\mathbf{r})\) is either a plane wave expansion or a spherical wave expansion, and \(D_i\) is its domain of validity.

compatible(other)

Returns always true, because any field expansion can be added to a piecewise field expansion.

diverging(x, y, z)

Test if points are in domain where expansion could diverge.

Parameters:
  • x (numpy.ndarray) – x-coordinates of query points
  • y (numpy.ndarray) – y-coordinates of query points
  • z (numpy.ndarray) – z-coordinates of query points
Returns:

numpy.ndarray of bool datatype indicating if points are inside divergence domain.

electric_field(x, y, z)

Evaluate electric field.

Parameters:
  • x (numpy.ndarray) – x-coordinates of query points
  • y (numpy.ndarray) – y-coordinates of query points
  • z (numpy.ndarray) – z-coordinates of query points
Returns:

Tuple of (E_x, E_y, E_z) numpy.ndarray objects with the Cartesian coordinates of complex electric field.

valid(x, y, z)

Test if points are in definition range of the expansion.

Parameters:
  • x (numpy.ndarray) – x-coordinates of query points
  • y (numpy.ndarray) – y-coordinates of query points
  • z (numpy.ndarray) – z-coordinates of query points
Returns:

numpy.ndarray of bool datatype indicating if points are inside definition domain.

class smuthi.field_expansion.PlaneWaveExpansion(k, k_parallel='default', azimuthal_angles='default', kind=None, reference_point=None, lower_z=-inf, upper_z=inf)

A class to manage plane wave expansions of the form

\[\mathbf{E}(\mathbf{r}) = \sum_{j=1}^2 \iint \mathrm{d}^2\mathbf{k}_\parallel \, g_{j}(\kappa, \alpha) \mathbf{\Phi}^\pm_j(\kappa, \alpha; \mathbf{r} - \mathbf{r}_i)\]

for \(\mathbf{r}\) located in a layer defined by \(z\in [z_{min}, z_{max}]\) and \(\mathrm{d}^2\mathbf{k}_\parallel = \kappa\,\mathrm{d}\alpha\,\mathrm{d}\kappa\).

The double integral runs over \(\alpha\in[0, 2\pi]\) and \(\kappa\in[0,\kappa_\mathrm{max}]\). Further, \(\mathbf{\Phi}^\pm_j\) are the PVWFs, see plane_vector_wave_function().

Internally, the expansion coefficients \(g_{ij}^\pm(\kappa, \alpha)\) are stored as a 3-dimensional numpy ndarray.

If the attributes k_parallel and azimuthal_angles have only a single entry, a discrete distribution is assumed:

\[g_{j}^\pm(\kappa, \alpha) \sim \delta^2(\mathbf{k}_\parallel - \mathbf{k}_{\parallel, 0})\]
Parameters:
  • k (float) – wavenumber in layer where expansion is valid
  • k_parallel (numpy ndarray) – array of in-plane wavenumbers (can be float or complex) If ‘default’, use smuthi.coordinates.default_k_parallel
  • azimuthal_angles (numpy ndarray) – \(\alpha\), from 0 to \(2\pi\) If ‘default’, use smuthi.coordinates.default_azimuthal_angles
  • kind (str) – ‘upgoing’ for \(g^+\) and ‘downgoing’ for \(g^-\) type expansions
  • reference_point (list or tuple) – [x, y, z]-coordinates of point relative to which the plane waves are defined.
  • lower_z (float) – the expansion is valid on and above that z-coordinate
  • upper_z (float) – the expansion is valid below that z-coordinate
coefficients

numpy ndarray – coefficients[j, k, l] contains \(g^\pm_{j}(\kappa_{k}, \alpha_{l})\)

azimuthal_angle_grid()

Meshgrid of azimuthal_angles with respect to n_effective

compatible(other)

Check if two plane wave expansions are compatible in the sense that they can be added coefficient-wise

Parameters:other (FieldExpansion) – expansion object to add to this object
Returns:bool (true if compatible, false else)
diverging(x, y, z)

Test if points are in domain where expansion could diverge.

Parameters:
  • x (numpy.ndarray) – x-coordinates of query points
  • y (numpy.ndarray) – y-coordinates of query points
  • z (numpy.ndarray) – z-coordinates of query points
Returns:

numpy.ndarray of bool datatype indicating if points are inside divergence domain.

electric_field(x, y, z)

Evaluate electric field.

Parameters:
  • x (numpy.ndarray) – x-coordinates of query points
  • y (numpy.ndarray) – y-coordinates of query points
  • z (numpy.ndarray) – z-coordinates of query points
Returns:

Tuple of (E_x, E_y, E_z) numpy.ndarray objects with the Cartesian coordinates of complex electric field.

k_parallel_grid()

Meshgrid of n_effective with respect to azimuthal_angles

k_z()
k_z_grid()
valid(x, y, z)

Test if points are in definition range of the expansion.

Parameters:
  • x (numpy.ndarray) – x-coordinates of query points
  • y (numpy.ndarray) – y-coordinates of query points
  • z (numpy.ndarray) – z-coordinates of query points
Returns:

numpy.ndarray of bool datatype indicating if points are inside definition domain.

class smuthi.field_expansion.SphericalWaveExpansion(k, l_max, m_max=None, kind=None, reference_point=None, lower_z=-inf, upper_z=inf, inner_r=0, outer_r=inf)

A class to manage spherical wave expansions of the form

\[\mathbf{E}(\mathbf{r}) = \sum_{\tau=1}^2 \sum_{l=1}^\infty \sum_{m=-l}^l a_{\tau l m} \mathbf{\Psi}^{(\nu)}_{\tau l m}(\mathbf{r} - \mathbf{r}_i)\]

for \(\mathbf{r}\) located in a layer defined by \(z\in [z_{min}, z_{max}]\) and where \(\mathbf{\Psi}^{(\nu)}_{\tau l m}\) are the SVWFs, see smuthi.vector_wave_functions.spherical_vector_wave_function().

Internally, the expansion coefficients \(a_{\tau l m}\) are stored as a 1-dimensional array running over a multi index \(n\) subsumming over the SVWF indices \((\tau,l,m)\). The mapping from the SVWF indices to the multi index is organized by the function multi_to_single_index().

Parameters:
  • k (float) – wavenumber in layer where expansion is valid
  • l_max (int) – maximal multipole degree \(l_\mathrm{max}\geq 1\) where to truncate the expansion.
  • m_max (int) – maximal multipole order \(0 \leq m_\mathrm{max} \leq l_\mathrm{max}\) where to truncate the expansion.
  • kind (str) – ‘regular’ for \(\nu=1\) or ‘outgoing’ for \(\nu=3\)
  • reference_point (list or tuple) – [x, y, z]-coordinates of point relative to which the spherical waves are considered (e.g., particle center).
  • lower_z (float) – the expansion is valid on and above that z-coordinate
  • upper_z (float) – the expansion is valid below that z-coordinate
  • inner_r (float) – radius inside which the expansion diverges (e.g. circumscribing sphere of particle)
  • outer_r (float) – radius outside which the expansion diverges
coefficients

numpy ndarray – expansion coefficients \(a_{\tau l m}\) ordered by multi index n

coefficients_tlm(tau, l, m)

SWE coefficient for given (tau, l, m)

Parameters:
  • tau (int) – SVWF polarization (0 for spherical TE, 1 for spherical TM)
  • l (int) – SVWF degree
  • m (int) – SVWF order
Returns:

SWE coefficient

compatible(other)

Check if two spherical wave expansions are compatible in the sense that they can be added coefficient-wise

Parameters:other (FieldExpansion) – expansion object to add to this object
Returns:bool (true if compatible, false else)
diverging(x, y, z)

Test if points are in domain where expansion could diverge.

Parameters:
  • x (numpy.ndarray) – x-coordinates of query points
  • y (numpy.ndarray) – y-coordinates of query points
  • z (numpy.ndarray) – z-coordinates of query points
Returns:

numpy.ndarray of bool datatype indicating if points are inside divergence domain.

electric_field(x, y, z)

Evaluate electric field.

Parameters:
  • x (numpy.ndarray) – x-coordinates of query points
  • y (numpy.ndarray) – y-coordinates of query points
  • z (numpy.ndarray) – z-coordinates of query points
Returns:

Tuple of (E_x, E_y, E_z) numpy.ndarray objects with the Cartesian coordinates of complex electric field.

valid(x, y, z)

Test if points are in definition range of the expansion.

Parameters:
  • x (numpy.ndarray) – x-coordinates of query points
  • y (numpy.ndarray) – y-coordinates of query points
  • z (numpy.ndarray) – z-coordinates of query points
Returns:

numpy.ndarray of bool datatype indicating if points are inside definition domain.

smuthi.field_expansion.blocksize(l_max, m_max)

Number of coefficients in outgoing or regular spherical wave expansion for a single particle.

Parameters:
  • l_max (int) – Maximal multipole degree
  • m_max (int) – Maximal multipole order
Returns:

Number of indices for one particle, which is the maximal index plus 1.

smuthi.field_expansion.multi_to_single_index(tau, l, m, l_max, m_max)

Unique single index for the totality of indices characterizing a svwf expansion coefficient.

The mapping follows the scheme:

single index spherical wave expansion indices
\(n\) \(\tau\) \(l\) \(m\)
1 1 1 -1
2 1 1 0
3 1 1 1
4 1 2 -2
5 1 2 -1
6 1 2 0
... ... ... ...
... 1 l_max m_max
... 2 1 -1
... ... ... ...
Parameters:
  • tau (int) – Polarization index :math:`tau`(0=spherical TE, 1=spherical TM)
  • l (int) – Degree \(l\) (1, ..., lmax)
  • m (int) – Order \(m\) (-min(l,mmax),...,min(l,mmax))
  • l_max (int) – Maximal multipole degree
  • m_max (int) – Maximal multipole order
Returns:

single index (int) subsuming \((\tau, l, m)\)

smuthi.field_expansion.pwe_to_ff_conversion(vacuum_wavelength, plane_wave_expansion)

Compute the far field of a plane wave expansion object.

Parameters:
Returns:

A smuthi.field_evaluation.FarField object containing the far field intensity.

smuthi.field_expansion.pwe_to_swe_conversion(pwe, l_max, m_max, reference_point)

Convert plane wave expansion object to a spherical wave expansion object.

Parameters:
  • pwe (PlaneWaveExpansion) – Plane wave expansion to be converted
  • l_max (int) – Maximal multipole degree of spherical wave expansion
  • m_max (int) – Maximal multipole order of spherical wave expansion
  • reference_point (list) – Coordinates of reference point in the format [x, y, z]
Returns:

SphericalWaveExpansion object.

smuthi.field_expansion.swe_to_pwe_conversion(swe, k_parallel='default', azimuthal_angles='default', layer_system=None, layer_number=None, layer_system_mediated=False)

Convert SphericalWaveExpansion object to a PlaneWaveExpansion object.

Parameters:
  • swe (SphericalWaveExpansion) – Spherical wave expansion to be converted
  • k_parallel (numpy array or str) – In-plane wavenumbers for the pwe object. If ‘default’, use smuthi.coordinates.default_k_parallel
  • azimuthal_angles (numpy array or str) – Azimuthal angles for the pwe object If ‘default’, use smuthi.coordinates.default_azimuthal_angles
  • layer_system (smuthi.layers.LayerSystem) – Stratified medium in which the origin of the SWE is located
  • layer_number (int) – Layer number in which the PWE should be valid.
  • layer_system_mediated (bool) – If True, the PWE refers to the layer system response of the SWE, otherwise it is the direct transform.
Returns:

Tuple of two PlaneWaveExpansion objects, first upgoing, second downgoing.

smuthi.graphical_output module

smuthi.initial_field module

smuthi.layers module

smuthi.memoizing module

Provide functionality to store intermediate results in lookup tables (memoize)

class smuthi.memoizing.Memoize(fn)

To be used as a decorator for functions that are memoized.

smuthi.particle_coupling module

smuthi.particles module

Provide class for the representation of scattering particles.

class smuthi.particles.FiniteCylinder(position=None, euler_angles=None, refractive_index=(1+0j), cylinder_radius=1, cylinder_height=1, l_max=None, m_max=None, t_matrix_method=None)

Particle subclass for finite cylinders.

Parameters:
  • position (list) – Particle position in the format [x, y, z] (length unit)
  • refractive_index (complex) – Complex refractive index of particle
  • cylinder_radius (float) – Radius of cylinder (length unit)
  • cylinder_height (float) – Height of cylinder, in z-direction if not rotated (length unit)
  • l_max (int) – Maximal multipole degree used for the spherical wave expansion of incoming and scattered field
  • m_max (int) – Maximal multipole order used for the spherical wave expansion of incoming and scattered field
class smuthi.particles.Particle(position=None, euler_angles=None, refractive_index=(1+0j), l_max=None, m_max=None)

Base class for scattering particles.

Parameters:
  • position (list) – Particle position in the format [x, y, z] (length unit)
  • euler_angles (list) – Particle Euler angles in the format [alpha, beta, gamma]
  • refractive_index (complex) – Complex refractive index of particle
  • l_max (int) – Maximal multipole degree used for the spherical wave expansion of incoming and scattered field
  • m_max (int) – Maximal multipole order used for the spherical wave expansion of incoming and scattered field
class smuthi.particles.Sphere(position=None, refractive_index=(1+0j), radius=1, l_max=None, m_max=None)

Particle subclass for spheres.

Parameters:
  • position (list) – Particle position in the format [x, y, z] (length unit)
  • refractive_index (complex) – Complex refractive index of particle
  • radius (float) – Particle radius (length unit)
  • l_max (int) – Maximal multipole degree used for the spherical wave expansion of incoming and scattered field
  • m_max (int) – Maximal multipole order used for the spherical wave expansion of incoming and scattered field
  • t_matrix_method (dict) – Dictionary containing the parameters for the algorithm to compute the T-matrix
class smuthi.particles.Spheroid(position=None, euler_angles=None, refractive_index=(1+0j), semi_axis_c=1, semi_axis_a=1, l_max=None, m_max=None, t_matrix_method=None)

Particle subclass for spheroids.

Parameters:
  • position (list) – Particle position in the format [x, y, z] (length unit)
  • refractive_index (complex) – Complex refractive index of particle
  • semi_axis_c (float) – Spheroid half axis in direction of axis of revolution (z-axis if not rotated)
  • semi_axis_a (float) – Spheroid half axis in lateral direction (x- and y-axis if not rotated)
  • l_max (int) – Maximal multipole degree used for the spherical wave expansion of incoming and scattered field
  • m_max (int) – Maximal multipole order used for the spherical wave expansion of incoming and scattered field
  • t_matrix_method (dict) – Dictionary containing the parameters for the algorithm to compute the T-matrix

smuthi.post_processing module

smuthi.read_input module

smuthi.scattered_field module

Manage post processing steps to evaluate the scattered near and far field

smuthi.scattered_field.extinction_cross_section(initial_field, particle_list, layer_system)

Evaluate and display the differential scattering cross section as a function of solid angle.

Parameters:
  • initial_field (smuthi.initial_field.PlaneWave) – Plane wave object
  • particle_list (list) – List of smuthi.particles.Particle objects
  • layer_system (smuthi.layers.LayerSystem) – Representing the stratified medium
Returns:

Dictionary with following entries
  • ‘forward’: Extinction in the positinve z-direction (top layer)
  • ‘backward’: Extinction in the negative z-direction (bottom layer)

smuthi.scattered_field.scattered_far_field(vacuum_wavelength, particle_list, layer_system, polar_angles='default', azimuthal_angles='default')

Evaluate the scattered far field.

Parameters:
  • vacuum_wavelength (float) – in length units
  • particle_list (list) – list of smuthi.Particle objects
  • layer_system (smuthi.layers.LayerSystem) – represents the stratified medium
  • polar_angles (numpy.ndarray or str) – polar angles values (radian). if ‘default’, use smuthi.coordinates.default_polar_angles
  • azimuthal_angles (numpy.ndarray or str) – azimuthal angle values (radian) if ‘default’, use smuthi.coordinates.default_azimuthal_angles
Returns:

A FarField object of the scattered field.

smuthi.scattered_field.scattered_field_piecewise_expansion(vacuum_wavelength, particle_list, layer_system, k_parallel='default', azimuthal_angles='default')

Compute a piecewise field expansion of the scattered field.

Parameters:
  • vacuum_wavelength (float) – vacuum wavelength
  • particle_list (list) – list of smuthi.particles.Particle objects
  • layer_system (smuthi.layers.LayerSystem) – stratified medium
  • k_parallel (numpy.ndarray or str) – in-plane wavenumbers array. if ‘default’, use smuthi.coordinates.default_k_parallel
  • azimuthal_angles (numpy.ndarray or str) – azimuthal angles array if ‘default’, use smuthi.coordinates.default_azimuthal_angles
Returns:

scattered field as smuthi.field_expansion.PiecewiseFieldExpansion object

smuthi.scattered_field.scattered_field_pwe(vacuum_wavelength, particle_list, layer_system, layer_number, k_parallel='default', azimuthal_angles='default', include_direct=True, include_layer_response=True)

Calculate the plane wave expansion of the scattered field of a set of particles.

Parameters:
  • vacuum_wavelength (float) – Vacuum wavelength (length unit)
  • particle_list (list) – List of Particle objects
  • layer_system (smuthi.layers.LayerSystem) – Stratified medium
  • layer_number (int) – Layer number in which the plane wave expansion should be valid
  • k_parallel (numpy.ndarray or str) – in-plane wavenumbers array. if ‘default’, use smuthi.coordinates.default_k_parallel
  • azimuthal_angles (numpy.ndarray or str) – azimuthal angles array if ‘default’, use smuthi.coordinates.default_azimuthal_angles
  • include_direct (bool) – If True, include the direct scattered field
  • include_layer_response (bool) – If True, include the layer system response
Returns:

A tuple of PlaneWaveExpansion objects for upgoing and downgoing waves.

smuthi.scattered_field.scattering_cross_section(initial_field, particle_list, layer_system, polar_angles='default', azimuthal_angles='default')

Evaluate and display the differential scattering cross section as a function of solid angle.

Parameters:
  • initial_field (smuthi.initial.PlaneWave) – Initial Plane wave
  • particle_list (list) – scattering particles
  • layer_system (smuthi.layers.LayerSystem) – stratified medium
  • polar_angles (numpy.ndarray or str) – polar angles values (radian). if ‘default’, use smuthi.coordinates.default_polar_angles
  • azimuthal_angles (numpy.ndarray or str) – azimuthal angle values (radian) if ‘default’, use smuthi.coordinates.default_azimuthal_angles
Returns:

A tuple of FarField objects, one for forward scattering (i.e., into the top hemisphere) and one for backward scattering (bottom hemisphere).

smuthi.scattered_field.total_far_field(initial_field, particle_list, layer_system, polar_angles='default', azimuthal_angles='default')

Evaluate the total far field, the initial far field and the scattered far field. Cannot be used if initial field is a plane wave.

Parameters:
  • initial_field (smuthi.initial_field.InitialField) – represents the initial field
  • particle_list (list) – list of smuthi.Particle objects
  • layer_system (smuthi.layers.LayerSystem) – represents the stratified medium
  • polar_angles (numpy.ndarray or str) – polar angles values (radian). if ‘default’, use smuthi.coordinates.default_polar_angles
  • azimuthal_angles (numpy.ndarray or str) – azimuthal angle values (radian) if ‘default’, use smuthi.coordinates.default_azimuthal_angles
Returns:

A tuple of three FarField objects for total, initial and scattered far field. Mind that the scattered far field has no physical meaning and is for illustration purposes only.

smuthi.simulation module

smuthi.spherical_functions module

smuthi.spherical_functions.double_factorial(n)

Return double factorial.

Parameters:n (int) – Argument (non-negative)
Returns:Double factorial of n
smuthi.spherical_functions.dx_xh(n, x)

Derivative of \(x h_n(x)\), where \(h_n(x)\) is the spherical Hankel function.

Parameters:
  • n (int) – (n>0): Order of spherical Bessel function
  • x (array, complex or float) – Argument for spherical Hankel function
Returns:

Derivative \(\partial_x(x h_n(x))\) as array.

smuthi.spherical_functions.dx_xj(n, x)

Derivative of \(x j_n(x)\), where \(j_n(x)\) is the spherical Bessel function.

Parameters:
  • n (int) – (n>0): Order of spherical Bessel function
  • x (array, complex or float) – Argument for spherical Bessel function
Returns:

Derivative \(\partial_x(x j_n(x))\) as array.

smuthi.spherical_functions.factorial(n)

Return factorial.

Parameters:n (int) – Argument (non-negative)
Returns:Factorial of n
smuthi.spherical_functions.legendre_normalized(ct, st, lmax)

Return the normalized associated Legendre function \(P_l^m(\cos\theta)\) and the angular functions \(\pi_l^m(\cos \theta)\) and \(\tau_l^m(\cos \theta)\), as defined in A. Doicu, T. Wriedt, and Y. A. Eremin: “Light Scattering by Systems of Particles”, Springer-Verlag, 2006. Two arguments (ct and st) are passed such that the function is valid for general complex arguments, while the branch cuts are defined by the user already in the definition of st.

Parameters:
  • ct (array) – cosine of theta (or kz/k)
  • st (array) – sine of theta (or kp/k), need to have same dimension as ct, and st**2+ct**2=1 is assumed
  • lmax (int) – maximal multipole order
Returns:

  • list plm[l][m] contains \(P_l^m(\cos \theta)\). The entries of the list have same dimension as ct (and st)
  • list pilm[l][m] contains \(\pi_l^m(\cos \theta)\).
  • list taulm[l][m] contains \(\tau_l^m(\cos \theta)\).

smuthi.spherical_functions.spherical_bessel(n, x)

Spherical Bessel function. This is a wrapper for scipy.special.sph_jn to make it operate on numpy arrays.

As soon as some bug for complex arguments is resolved, this can be replaced by scipy.special.spherical_jn. https://github.com/ContinuumIO/anaconda-issues/issues/1415

Parameters:
  • n (int) – Order of spherical Bessel function
  • x (array, complex or float) – Argument for Bessel function
Returns:

Spherical Bessel function as array.

smuthi.spherical_functions.spherical_hankel(n, x)

Spherical Hankel function of first kind.

Parameters:
  • n (int) – Order of spherical Bessel function
  • x (array, complex or float) – Argument for Hankel function
Returns:

Spherical Hankel function as array.

smuthi.t_matrix module

smuthi.t_matrix.mie_coefficient(tau, l, k_medium, k_particle, radius)

Return the Mie coefficients of a sphere.

Parameters:
  • integer (l) – spherical polarization, 0 for spherical TE and 1 for spherical TM
  • integer – l=1,... multipole degree (polar quantum number)
  • float or complex (k_particle) – wavenumber in surrounding medium (inverse length unit)
  • float or complex – wavenumber inside sphere (inverse length unit)
  • float (radius) – radius of sphere (length unit)
Returns:

Mie coefficients as complex

smuthi.t_matrix.rotate_t_matrix(t, euler_angles)

Placeholder for a proper T-matrix rotation routine

smuthi.t_matrix.t_matrix(vacuum_wavelength, n_medium, particle)

Return the T-matrix of a particle.

Parameters:
  • vacuum_wavelength (float) –
  • n_medium (float or complex) – Refractive index of surrounding medium
  • particle (smuthi.particles.Particle) – Particle object
Returns:

T-matrix as ndarray

smuthi.t_matrix.t_matrix_sphere(k_medium, k_particle, radius, l_max, m_max)

T-matrix of a spherical scattering object.

Parameters:
  • k_medium (float or complex) – Wavenumber in surrounding medium (inverse length unit)
  • k_particle (float or complex) – Wavenumber inside sphere (inverse length unit)
  • radius (float) – Radius of sphere (length unit)
  • l_max (int) – Maximal multipole degree
  • m_max (int) – Maximal multipole order
  • blocksize (int) – Total number of index combinations
  • multi_to_single_index_map (function) – A function that maps the SVWF indices (tau, l, m) to a single index
Returns:

T-matrix as ndarray

smuthi.vector_wave_functions module

smuthi.vector_wave_functions.ab5_coefficients(l1, m1, l2, m2, p, symbolic=False)

a5 and b5 are the coefficients used in the evaluation of the SVWF translation operator. Their computation is based on the sympy.physics.wigner package and is performed with symbolic numbers.

Parameters:
  • l1 (int) – l=1,...: Original wave’s SVWF multipole degree
  • m1 (int) – m=-l,...,l: Original wave’s SVWF multipole order
  • l2 (int) – l=1,...: Partial wave’s SVWF multipole degree
  • m2 (int) – m=-l,...,l: Partial wave’s SVWF multipole order
  • p (int) – p parameter
  • symbolic (bool) – If True, symbolic numbers are returned. Otherwise, complex.
Returns:

A tuple (a5, b5) where a5 and b5 are symbolic or complex.

smuthi.vector_wave_functions.plane_vector_wave_function(x, y, z, kp, alpha, kz, pol)

Electric field components of plane wave (PVWF).

\[\mathbf{\Phi}_j = \exp ( \mathrm{i} \mathbf{k} \cdot \mathbf{r} ) \hat{ \mathbf{e} }_j\]

with \(\hat{\mathbf{e}}_0\) denoting the unit vector in azimuthal direction (‘TE’ or ‘s’ polarization), and \(\hat{\mathbf{e}}_1\) denoting the unit vector in polar direction (‘TM’ or ‘p’ polarization).

The input arrays should have one of the following dimensions:

  • x,y,z: (N x 1) matrix
  • kp,alpha,kz: (1 x M) matrix
  • Ex, Ey, Ez: (M x N) matrix

or

  • x,y,z: (M x N) matrix
  • kp,alpha,kz: scalar
  • Ex, Ey, Ez: (M x N) matrix
Parameters:
  • x (numpy.ndarray) – x-coordinate of position where to test the field (length unit)
  • y (numpy.ndarray) – y-coordinate of position where to test the field
  • z (numpy.ndarray) – z-coordinate of position where to test the field
  • kp (numpy.ndarray) – parallel component of k-vector (inverse length unit)
  • alpha (numpy.ndarray) – azimthal angle of k-vector (rad)
  • kz (numpy.ndarray) – z-component of k-vector (inverse length unit)
  • pol (int) – Polarization (0=TE, 1=TM)
Returns:

  • x-coordinate of PVWF electric field (numpy.ndarray)
  • y-coordinate of PVWF electric field (numpy.ndarray)
  • z-coordinate of PVWF electric field (numpy.ndarray)

smuthi.vector_wave_functions.spherical_vector_wave_function(x, y, z, k, nu, tau, l, m)

Electric field components of spherical vector wave function (SVWF). The conventions are chosen according to A. Doicu, T. Wriedt, and Y. A. Eremin: “Light Scattering by Systems of Particles”, Springer-Verlag, 2006

Parameters:
  • x (numpy.ndarray) – x-coordinate of position where to test the field (length unit)
  • y (numpy.ndarray) – y-coordinate of position where to test the field
  • z (numpy.ndarray) – z-coordinate of position where to test the field
  • k (float or complex) – wavenumber (inverse length unit)
  • nu (int) – 1 for regular waves, 3 for outgoing waves
  • tau (int) – spherical polarization, 0 for spherical TE and 1 for spherical TM
  • l (int) – l=1,... multipole degree (polar quantum number)
  • m (int) – m=-l,...,l multipole order (azimuthal quantum number)
Returns:

  • x-coordinate of SVWF electric field (numpy.ndarray)
  • y-coordinate of SVWF electric field (numpy.ndarray)
  • z-coordinate of SVWF electric field (numpy.ndarray)

smuthi.vector_wave_functions.transformation_coefficients_vwf(tau, l, m, pol, kp=None, kz=None, pilm_list=None, taulm_list=None, dagger=False)

Transformation coefficients B to expand SVWF in PVWF and vice versa:

\[B_{\tau l m,j}(x) = -\frac{1}{\mathrm{i}^{l+1}} \frac{1}{\sqrt{2l(l+1)}} (\mathrm{i} \delta_{j1} + \delta_{j2}) (\delta_{\tau j} \tau_l^{|m|}(x) + (1-\delta_{\tau j} m \pi_l^{|m|}(x))\]

For the definition of the \(\tau_l^m\) and \(\pi_l^m\) functions, see A. Doicu, T. Wriedt, and Y. A. Eremin: “Light Scattering by Systems of Particles”, Springer-Verlag, 2006

Parameters:
  • tau (int) – SVWF polarization, 0 for spherical TE, 1 for spherical TM
  • l (int) – l=1,... SVWF multipole degree
  • m (int) – m=-l,...,l SVWF multipole order
  • pol (int) – PVWF polarization, 0 for TE, 1 for TM
  • kp (numpy array) – PVWF in-plane wavenumbers
  • kz (numpy array) – complex numpy-array: PVWF out-of-plane wavenumbers
  • pilm_list (list) – 2D list numpy-arrays: alternatively to kp and kz, pilm and taulm as generated with legendre_normalized can directly be handed
  • taulm_list (list) – 2D list numpy-arrays: alternatively to kp and kz, pilm and taulm as generated with legendre_normalized can directly be handed
  • dagger (bool) – switch on when expanding PVWF in SVWF and off when expanding SVWF in PVWF
Returns:

Transformation coefficient as array (size like kp).

smuthi.vector_wave_functions.translation_coefficients_svwf(tau1, l1, m1, tau2, l2, m2, k, d, sph_hankel=None, legendre=None, exp_immphi=None)

Coefficients of the translation operator for the expansion of an outgoing spherical wave in terms of regular spherical waves with respect to a different origin:

\[\mathbf{\Psi}_{\tau l m}^{(3)}(\mathbf{r} + \mathbf{d} = \sum_{\tau'} \sum_{l'} \sum_{m'} A_{\tau l m, \tau' l' m'} (\mathbf{d}) \mathbf{\Psi}_{\tau' l' m'}^{(1)}(\mathbf{r})\]

for \(|\mathbf{r}|<|\mathbf{d}|\).

Parameters:
  • tau1 (int) – tau1=0,1: Original wave’s spherical polarization
  • l1 (int) – l=1,...: Original wave’s SVWF multipole degree
  • m1 (int) – m=-l,...,l: Original wave’s SVWF multipole order
  • tau2 (int) – tau2=0,1: Partial wave’s spherical polarization
  • l2 (int) – l=1,...: Partial wave’s SVWF multipole degree
  • m2 (int) – m=-l,...,l: Partial wave’s SVWF multipole order
  • k (float or complex) – wavenumber (inverse length unit)
  • d (list) – translation vectors in format [dx, dy, dz] (length unit) dx, dy, dz can be scalars or ndarrays
  • sph_hankel (list) – Optional. sph_hankel[i] contains the spherical hankel funciton of degree i, evaluated at k*d where d is the norm of the distance vector(s)
  • legendre (list) – Optional. legendre[l][m] contains the legendre function of order l and degree m, evaluated at cos(theta) where theta is the polar angle(s) of the distance vector(s)
Returns:

translation coefficient A (complex)

smuthi.vector_wave_functions.translation_coefficients_svwf_out_to_out(tau1, l1, m1, tau2, l2, m2, k, d, sph_bessel=None, legendre=None, exp_immphi=None)

Coefficients of the translation operator for the expansion of an outgoing spherical wave in terms of outgoing spherical waves with respect to a different origin:

\[\mathbf{\Psi}_{\tau l m}^{(3)}(\mathbf{r} + \mathbf{d} = \sum_{\tau'} \sum_{l'} \sum_{m'} A_{\tau l m, \tau' l' m'} (\mathbf{d}) \mathbf{\Psi}_{\tau' l' m'}^{(3)}(\mathbf{r})\]

for \(|\mathbf{r}|>|\mathbf{d}|\).

Parameters:
  • tau1 (int) – tau1=0,1: Original wave’s spherical polarization
  • l1 (int) – l=1,...: Original wave’s SVWF multipole degree
  • m1 (int) – m=-l,...,l: Original wave’s SVWF multipole order
  • tau2 (int) – tau2=0,1: Partial wave’s spherical polarization
  • l2 (int) – l=1,...: Partial wave’s SVWF multipole degree
  • m2 (int) – m=-l,...,l: Partial wave’s SVWF multipole order
  • k (float or complex) – wavenumber (inverse length unit)
  • d (list) – translation vectors in format [dx, dy, dz] (length unit) dx, dy, dz can be scalars or ndarrays
  • sph_bessel (list) – Optional. sph_bessel[i] contains the spherical Bessel funciton of degree i, evaluated at k*d where d is the norm of the distance vector(s)
  • legendre (list) – Optional. legendre[l][m] contains the legendre function of order l and degree m, evaluated at cos(theta) where theta is the polar angle(s) of the distance vector(s)
Returns:

translation coefficient A (complex)

smuthi.nfmds package

Check if NFM-DS is installed and otherwise do so.

smuthi.nfmds.install_nfmds()

smuthi.nfmds.t_matrix_axsym module

Functions to call NFM-DS (null field method with discrete sources) Fortran code by Doicu et al. for the generation of T-matrices for non-spherical particles. The Fortran code comes with the book A. Doicu, T. Wriedt, and Y. A. Eremin: Light Scattering by Systems of Particles, 1st ed. Berlin, Heidelberg: Springer-Verlag, 2006. and can also be downloaded from https://scattport.org/index.php/programs-menu/t-matrix-codes-menu/239-nfm-ds

smuthi.nfmds.t_matrix_axsym.taxsym_read_tmatrix(filename, l_max, m_max)

Export TAXSYM.f90 output to SMUTHI T-matrix.

Parameters:
  • filename (str) – Name of the file containing the T-matrix output of TAXSYM.f90
  • l_max (int) – Maximal multipole degree
  • m_max (int) – Maximal multipole order
Returns:

T-matrix as numpy.ndarray

smuthi.nfmds.t_matrix_axsym.taxsym_run()

Call TAXSYM.f90 routine.

smuthi.nfmds.t_matrix_axsym.taxsym_write_input_cylinder(vacuum_wavelength=None, layer_refractive_index=None, particle_refractive_index=None, cylinder_height=None, cylinder_radius=None, use_ds=True, nint=None, nrank=None, filename='T_matrix_cylinder.dat')

Generate input file for the TAXSYM.f90 routine for the simulation of a finite cylinder.

Parameters:
  • vacuum_wavelength (float) –
  • layer_refractive_index (float) – Real refractive index of layer (complex values are not allowed)
  • particle_refractive_index (float or complex) – Complex refractive index of cylinder
  • cylinder_height (float) – Height of cylinder (length unit)
  • cylinder_radius (float) – Radius of cylinder (length unit)
  • use_ds (bool) – Flag to switch the use of discrete sources on (True) and off (False)
  • nint (int) – Nint parameter for internal use of NFM-DS (number of points along integral). Higher value is more accurate and takes longer
  • nrank (int) – l_max used internally in NFM-DS
  • filename (str) – Name of the file in which the T-matrix is stored
smuthi.nfmds.t_matrix_axsym.taxsym_write_input_spheroid(vacuum_wavelength=None, layer_refractive_index=None, particle_refractive_index=None, semi_axis_c=None, semi_axis_a=None, use_ds=True, nint=None, nrank=None, filename='T_matrix_spheroid.dat')

Generate input file for the TAXSYM.f90 routine for the simulation of a spheroid.

Parameters:
  • vacuum_wavelength (float) –
  • layer_refractive_index (float) – Real refractive index of layer (complex values are not allowed)
  • particle_refractive_index (float or complex) – Complex refractive index of spheroid
  • semi_axis_c (float) – Semi axis of spheroid along rotation axis
  • semi_axis_a (float) – Semi axis of spheroid along lateral direction
  • use_ds (bool) – Flag to switch the use of discrete sources on (True) and off (False)
  • nint (int) – Nint parameter for internal use of NFM-DS (number of points along integral). Higher value is more accurate and takes longer
  • nrank (int) – l_max used internally in NFM-DS
  • filename (str) – Name of the file in which the T-matrix is stored
smuthi.nfmds.t_matrix_axsym.tmatrix_cylinder(vacuum_wavelength=None, layer_refractive_index=None, particle_refractive_index=None, cylinder_height=None, cylinder_radius=None, l_max=None, m_max=None, use_ds=True, nint=None, nrank=None)

Return T-matrix for finite cylinder, using the TAXSYM.f90 routine from the NFM-DS.

Parameters:
  • vacuum_wavelength (float) –
  • layer_refractive_index (float) – Real refractive index of layer (complex values are not allowed).
  • particle_refractive_index (float or complex) – Complex refractive index of spheroid
  • cylinder_height (float) – Semi axis of spheroid along rotation axis
  • cylinder_radius (float) – Semi axis of spheroid along lateral direction
  • l_max (int) – Maximal multipole degree
  • m_max (int) – Maximal multipole order
  • use_ds (bool) – Flag to switch the use of discrete sources on (True) and off (False)
  • nint (int) – Nint parameter for internal use of NFM-DS (number of points along integral). Higher value is more accurate and takes longer
  • nrank (int) – l_max used internally in NFM-DS
Returns:

T-matrix as numpy.ndarray

smuthi.nfmds.t_matrix_axsym.tmatrix_spheroid(vacuum_wavelength=None, layer_refractive_index=None, particle_refractive_index=None, semi_axis_c=None, semi_axis_a=None, l_max=None, m_max=None, use_ds=True, nint=None, nrank=None)

T-matrix for spheroid, using the TAXSYM.f90 routine from the NFM-DS.

Parameters:
  • vacuum_wavelength (float) –
  • layer_refractive_index (float) – Real refractive index of layer (complex values are not allowed).
  • particle_refractive_index (float or complex) – Complex refractive index of spheroid
  • semi_axis_c (float) – Semi axis of spheroid along rotation axis
  • semi_axis_a (float) – Semi axis of spheroid along lateral direction
  • l_max (int) – Maximal multipole degree
  • m_max (int) – Maximal multipole order
  • use_ds (bool) – Flag to switch the use of discrete sources on (True) and off (False)
  • nint (int) – Nint parameter for internal use of NFM-DS (number of points along integral). Higher value is more accurate and takes longer
  • nrank (int) – l_max used internally in NFM-DS
Returns:

T-matrix as numpy.ndarray