Matrix Setup - Layers and Periodicity

The newest MSTM code (version 4) introduces periodic boundary conditions (in XY plane) and matrix layers (in Z direction). Periodic boundary can be setup with set_boundary method of SPR_v4 class, and layers - with set_layers. Both spectra and near field calculations are allowed.

Important: this works only in latest releases of MSTM. This wrapper was tested for 2023 year release (github.com/dmckwski/MSTM/tree/main/december2023). The compiled binaries are supplied in our repositary (https://github.com/lavakyan/mstm-spectrum/releases).

Warning

Near field calculations with MSTM v.4 sometimes give strange results. Investigation in the process.

Of course, no sphere allowed to intercept PBC and layer limits.

Periodic boundary conditions

Periodic boundary conditions (PBC) in MSTMv4 can be setup only in XY plane. The returned values are - reflectance (R), - absorptance (A), - transmittance (T) coefficients for parallel, orthogonal orientations (suffices ‘_par’ and ‘_ort’. Also stored the orientation-averaged values. The coefficients should be 0 .. 1 interval, but for closely placed spheres, with low gap between particles, T can become negative.

The extinction cross section can be calculated from transmittance as [Kreibig_book1995]:

\[\sigma_{ext} = - \frac{1}{n L} ln(T),\]

where n is a particle number density and L is an optical length. Note, that negative T will produce NaN values.

Since \(n = N / (L_x L_y L)\), optical length does not affect the cross section and

\[\sigma_{ext} = - \frac{L_x L_y}{N} ln(T),\]

where \(L_x\) and \(L_y\) are periodic cell parameters, and N is number of spheres.

Since other calculation modes returns efficiencies, not cross sections, the division by effective scattering area is performed in the code. However, the estimation of multple spheres area is tricky. Currently it is evaluated simply as \(\sum_i \pi R_i^2\) which may be wrong if spheres projections overlap.

The general recommendation is to use more reliable T, A and R values stored after caclulation as the fields of SPR object.

Example: chain of nanoparticles

The following script compares extinction spectra of - isolated particle, - finite group of particles in line, - infinite chain of particles (usinf PBC).

import numpy as np
import matplotlib.pyplot as plt
from mstm_studio.mstm_v4 import SPR_v4
from mstm_studio.mstm_spectrum import Material, Spheres, SingleSphere

# specify the path to gold dielectric function
mat1 = Material('nk/etaGold.txt')

# calculation parameters
wls = np.linspace(400, 700, 61)  # wavelength, nm
cell = 30  # distance between particle centers, nm
a = 10  # particle radius, nm

# create caclulator object
spr = SPR_v4(wls, mstm_path='~/bin/mstm2023.x')  # , temp_dir='./temp/')
spr.environment_material = 1.5  # glass

# 1. isolated particle
spheres = SingleSphere(0, 0, 0, a, mat1)
spr.set_spheres(spheres)
spr.set_incident_field(fixed=True, beta_angle=0.0, alpha_angle=0.0)
spr.simulate()
isol_ext = spr.extinction_par

# 2. N spheres chain
def ext_of_spheres_chain(N):
    ''' helper function to calc
    extinction for chain '''
    spheres = Spheres()  # empty
    for i in range(-N//2, N//2+1):
        spheres.append(
            SingleSphere(x=i*cell, y=0, z=0, a=a,
                         mat_filename=mat1))
    spr.set_spheres(spheres)
    spr.simulate()
    return spr.extinction_par

chain3_ext = ext_of_spheres_chain(N=3)
chain5_ext = ext_of_spheres_chain(N=5)
chain7_ext = ext_of_spheres_chain(N=7)

# use of PBC
spheres = SingleSphere(0, 0, 0, a, mat1)
spr.set_spheres(spheres)
spr.set_boundary(pbc=True,  # use periodic boundaries
                 cell_x=cell,  # chain direction
                 cell_y=100)  # some gap in off-chain direction
                              # tests show that effic. are indep. of it
spr.simulate()
pbc_ext = spr.extinction_par

# do plots:
plt.figure(figsize=(4,3))
plt.plot(wls, chain7_ext, linewidth=1, linestyle=':', color='green', label='chain of 7')
plt.plot(wls, chain5_ext, 'g-.', linewidth=1, label='chain of 5')
plt.plot(wls, chain3_ext, 'g-',  linewidth=1, label='chain of 3')
plt.plot(wls, pbc_ext, 'b-', linewidth=1.5, label='PBC')
plt.plot(wls, isol_ext, 'r--', linewidth=1.5, label='isolated')
plt.xlabel(r'Wavelength, nm')
plt.ylabel(r'Ext. effic.')
plt.legend()
plt.tight_layout()
plt.savefig('example_periodic_chain.png')
plt.show()

Output figure:

../_images/example_periodic_chain.png

Layered matrix

By default there is no layers and matrix spans for all space -inf < Z +inf. The material constants for matrix are set up in the usual way, by specifiing environment_material property. I.e. spr.environment_material = Material(1.5).

The single layer boundary is added by specifieng list of one material in set_layers method:

After that the negative Z half space (-inf < Z < 0) is described with environment_material and positive (0 < Z < +inf) with mat passed to this method.

The further increase of layer numbers required specification of list of thier depths. Note, that length of matrices list should be bigger by 1 than lengths of depths list. MSTMv4 has hard-coded maximal amount of layers.

Example: Silver particle in layer

The nearfield from 10 nm Ag particle in a water layer of 20 nm depth on a glass.

import numpy as np
import matplotlib.pyplot as plt
from mstm_studio.mstm_v4 import SPR_v4
from mstm_studio.mstm_spectrum import Material, Spheres, SingleSphere

# specify the path to gold dielectric function
mat1 = Material('nk/etaGold.txt')

# calculation parameters
wls = np.linspace(400, 700, 61)  # wavelength, nm
cell = 30  # distance between particle centers, nm
a = 10  # particle radius, nm

# create caclulator object
spr = SPR_v4(wls, mstm_path='~/bin/mstm2023.x')  # , temp_dir='./temp/')
spr.environment_material = 1.5  # glass

# 1. isolated particle
spheres = SingleSphere(0, 0, 0, a, mat1)
spr.set_spheres(spheres)
spr.set_incident_field(fixed=True, beta_angle=0.0, alpha_angle=0.0)
spr.simulate()
isol_ext = spr.extinction_par

# 2. N spheres chain
def ext_of_spheres_chain(N):
    ''' helper function to calc
    extinction for chain '''
    spheres = Spheres()  # empty
    for i in range(-N//2, N//2+1):
        spheres.append(
            SingleSphere(x=i*cell, y=0, z=0, a=a,
                         mat_filename=mat1))

Output figure:

../_images/example_silver_in_layer.png

MSTM v4 classes