Scattering
from a multilayered sphere
The routine returns the total extinction, dipole extinction, quadrupole extinction,
(elastic) scattering, and absorption cross sections, albedo, and, optionally, the
T-matrix elements for the scattering off a single multilayered sphere.
The number LCS of concentric sphere layers can be arbitrary.
The program has built-in checks for convergence of the output values.
By changing the parameter NMAT value, you can switch between
plurality of different materials. Currently, the following options are
available: dispersionless dielectric, Drude metal, Ag, Au, ZnS, Cu, Al, Pt, and Si.
The option of a perfect conductor (the Dirichlet boundary condition) is activated
by the .TRUE. value of the logical parameter YNPERFCON.
The algorithm has been described in details in Sec. 3 of
my article
You can download the source code "sphere.f"
here.
The instructions how to run and compile it, together with
the definition of its parameters, are written
within the comment lines at the beginning of the source file.
You will also need to download material data
file Audat.dat.
(For security reasons, *.dat files cannot be posted; therefore
the additional *.txt extension. You have to rename the data file
back to Audat.dat after the download.)
You can also download the code limited Windows executable sphere, which calculates the scattering off a single SiO2@Au coated sphere in air. Refractive index of SiO2 in the sphere core is taken to be 1.45. Refractive indices of gold are taken to be gold bulk data without any adjustable parameter as read from the data file Audat.dat, which is also to be downloaded.
Optical trapping of a multilayered sphere
The program opttrap.f
calculates the axial trapping efficiency Q for
a multilayered sphere. The axial trapping efficiency Q, or the
dimensionless normalized axial force, when multiplied by the focused laser
beam power P, determines the resulting axial force F on a sphere
immersed in a homogeneous medium with refractive index n according to the formula
F=(n/c) PQ. The program was made according to the theory of
P. A. Maia Neto and H. M. Nussenzveig, Theory of optical tweezers,
Europhys. Lett. 50, 702-708 (2000) by extending it to the case of
(multi)coated spheres. For the results obtained by this program, see:
Radiation power of a dipole radiating inside
or outside a multilayered sphere
Dipole can be located either outside or
embedded in any layer of the multilayered sphere.
Provides both radiative and non-radiative decay rates. Avoids
the cumbersome algorithm of H. Chew, P. J. McNulty, and M. Kerker,
Raman and fluorescent scattering by molecules embedded
in concentric spheres, J. Opt. Soc. Am. 66, 440-444 (1976) and
provides a new transfer matrix alternative.
Important for the description of inelastic light-scattering (fluorescence or Raman)
spectroscopy for characterizing single micrometer sized particles
(chemical speciation) of both inorganic and organic compounds, aerosols
and particulates, in LIDAR applications for remote sensing of both molecular
and particulate constituents of atmosphere, for identification of biological
particles in fluorescence-activated flow of cytomeres, to monitor specific
cell functions, or in the cell identification and sorting systems, in the
investigation thermal radiation from spherical microparticles, engineering of the
radiative decay for biophysical and biomedical applications, imaging of buried
saturated fluorescent molecules and imaging of surfaces in near-field
optical microscopy, in the study of the effects of light absorption
and amplification on the stimulated transition rates of the electric-dipole
emission of atoms or molecules embedded in micro- or nano-structured spheres,
stimulated Raman scattering, the interplay
between lasing and stimulated Raman scattering, etc.
Theoretical details can be found in my article
Here you can download F77 codes CHEW and
CHEWFS. The first one, CHEW,
calculates the radiative and nonradiative (due to Ohmic losses) decay rates
(the latter are sometime also called energy transfer rates),
whereas CHEWFS calculates induced
frequency shifts and total decay rates directly from
Green´s function at coinciding spatial arguments.
The sphere options are identical
to that in scattering from a multilayered sphere.
Thus, you can run it for both homogeneous and (multi)coated spheres.
You can also download a corresponding limited Windows executables.
Chew.exe (download also material data
file Audat.dat, which contains gold bulk data
without any adjustable parameter) calculates the electric dipole
radiated power and the dipole power loss due to Ohmic losses
for a coated SiO2@X@SiO2 sphere in water. Refractive index
of SiO2 is taken to be 1.45, that of water 1.33, and that of X you
can supply yourself. For a Windows executable
generated from CHEWFS, see an
on-line Appendix A of my
Chem. Phys. 317(1), 1-15 (2005).
In case you are using a Linux or Unix machine, here are the respective makefiles
mkchew and mkchewfs.
Some of the referred routines are attached at the end of the source files.
Effective medium parameters for a composite of (coated) spheres in the long
wavelength limit
Even though all the constituents media are nonmagnetic,
the effective medium can nevertheless exhibit magnetic properties.
In fact, a recent article by
C. L. Holloway, E. F. Kuester, J. Baker-Jarvis, and P. Kabos,
A double negative (DNG) composite medium composed of magnetodielectric spherical
particles embedded in a matrix,
IEEE Trans. Antennas Propag. 51, 2596-2603 (2003) suggested that a
negative-refraction material can, in principle, be
made as a composite of (not necessarily magnetic) spheres.
Together with Vassilis Yannopapas [see our manuscript
J. Phys.: Condens. Matter. 17, 3717-3734 (2005)]
we have recently expanded on the ideas of
Holloway et al. and showed that even a composite of inherently
non-magnetic homogeneous spheres can provide
a negative refractive index metamaterial.
The absence of inherently
magnetic materials in our porposal is a key
difference which makes it possible to achieve a negative
refractive index band even in the deep infrared region.
Our results are explained in the context of the extended
Maxwell-Garnett theory
(see my F77 code EFFE2P) and reproduced by
the ab initio calculations
based on multiple scattering theory.
Here you can also download a limited Windows executable
effe (download also material data file
Audat.dat, which contains gold bulk data
without any adjustable parameter), which calculates an effective medium properties
for a composite of coated Au@SiO2 spheres in air with volume fraction of 34%.
The code and executable are to be used for wavelengths typically at least one order of
magnitude longer than the sphere radius.
Scattering from an axially symmetric particle
If you use Mishchenko´s code
at http://www.giss.nasa.gov/~crmim
you will often find that its non-NAG alternative for a nonabsorbing dielectric particle
and both NAG and non-NAG alternatives in a strongly absorbing
case fail
to converge, complain about singular matrices, and do not provide any output.
An improvement of the Mishchenko code is provided which
expands its range especially in the direction of strongly absorbing (e.g., metallic)
particles
(see A. Moroz, Improvement of Mishchenko's T-matrix code for absorbing
particles,
Appl. Opt. 44,
3604-3609 (2005)
(accompanying F77 code is available here).
Rather simple amendment described in the article allows one to extend
Mishchenko´s code, which has been used to generate beautiful
benchmarking results for dielectric particles of various axisymmetric shapes
such as [see free online book by
M. I. Mishchenko, L. D. Travis, and A. A. Lacis,
Scattering, Absorption, and Emission of Light by Small
Particles (Cambridge University Press, Cambridge, 2002)]:
Electric intensity around an axially symmetric particle
Indispensible tool for the SERS, near-field imaging, near-Field Raman Microscopy,
near-field scanning optical microscopy (NSOM),
design of near-field optical probes, etc.
The particle options are identical to that in
scattering from an axially symmetric particle.
The set of programs calculates a free-space Green´s function of the Helmholtz equation in a given space dimension d which satisfies a Bloch periodic boundary condition in a subspace of the d-dimensional space, including the d-dimensional embedding space itself. The programs employ exponentially convergent representations of lattice sums as derived for different cases by Ewald, Ham and Segall, Kambe, Ozorio de Almeida, and, in the case of d=2 and Bloch periodicity in one dimension, by
The total lattice sum D_L is conventionally written as a sum of three partial contributions
If you look at my J. Phys. A: Math. Gen. 39, 11247-11282 (2006), the D_L^1 term for different subperiodic cases is given by Eqs. 83,85,86. The only special function is the incomplete gamma function which is easily calculated numerically on using recurrences given by Eqs. 88-89. There is nothing mysterious or difficult in calculating incomplete gamma function compared to e.g. calculation of the Bessel functions.
The D_L^2 term for different subperiodic cases is given by Eqs. 102-103. There is no special function involved. An integral therein is calculated easily numerically again on using recurrence given by Eq. 106.
The D_L^3 term for different subperiodic cases is given by Eqs. 118-119. There is no special function involved and the sum is calculated numerically in a straightforward way.
A historical overview, relevant definitions and references can be
found in the latter article and also here.
The set of programs forms an integral part of my multiple-scattering programs
for scattering off a periodic arrangement of scatterers.
The numerical code OLA for a 1D periodicity in 2D is available
here.
For a 2D periodicity in 3D, see an Appendix of the monograph by
which provides a full listing of the Pendry's LEED code. The routines of Pendry's LEED code can be used even by a pragmatic engineer as black box routines. The Pendry's LEED code comprises among others also, up to minor amendmends, the program lines that I have arranged in the code dlsumf2in3.f, which performs the Kambe's lattice sums.My code gff2in3.f can be then used as a main code. It generates the lattice points, calls the routines that calculate spherical harmonics and spherical Bessel functions, and helps you to obtain the free-space Green's function given the total lattice sum D_L generated by dlsumf2in3.f. The code has been employed in my publication
An optimization program for a finite periodic stacks of parallel layers
Implemenets the well-known Abeles transfer matrix method (reproduced in Sec. 1.6
"Wave propagation in a stratified medium" and Chapter 14 of Born and Wolf). Can be used to
optimize a finite stack parameters with respect to reflectivity, transmittivity,
and absorption. I have used it to optimize the properties of metallo-dielectric
stacks in my publication
One-dimensional local density of states
The source code ldosl1d.f calculates the local density of
states for a toy model
of a one-dimensional photonic crystal in one space dimension, which
has been discussed in my article
Two-dimensional layered KKR method
Solves the problem of scattering in two dimensions off a finite stack of discs arranged
periodically along parallel lines. Used to model the properties of a two-dimensional
periodic arrangement of infinite-length circular cylinders in the plane normal
to the cylinder axis.
Calculates reflection, transmission, and absorption of both classical
(water, acoustic, and electromagnetic) and quantum (electron) waves.
It is unique in a
sense that it allows for an optimization of the stack parameters with
respect to a maximal band gap, minimal absorbtion, etc. On a PC, up to cca 15.000
different configurations can be scanned overnight with program selecting the top
configurations, whether you want those with minimal absorption, the largest
bandgap, or otherwise. For a comparison of theory and experiment, see
some preliminary results in
Two-dimensional bulk KKR method
Solves the problem of scattering in an infinite periodic
arrangement of discs filling in entire two-dimensional space.
Calculates the two-dimensional band structure of both classical
(water, acoustic, and electromagnetic) and quantum (electron) waves.
It has been developed in collaboration with Han van der Lem. For the results in the
electromagnetic case, see
Three-dimensional layered KKR method
Solves the problem of scattering in three dimensions off a finite stack of
axially symmetric scatterers arranged
periodically along parallel (two-dimensional) planes.
Similar to the two-dimensional layered KKR method, but in one
dimension higher.
Improves its predecessor, the computer code of V. Yannopapas, N. Stefanou,
and A. Modinos,
Heterostructures of photonic crystals: frequency bands and transmission
coefficients, Comput. Phys. Commun. 113, 49 (1998)
[The scalar predecessor of the latter program has been developed in the book by
J. B. Pendry, Low Energy Electron Diffraction
(Academic Press, London, 1974)], in the following aspects:
Three-dimensional bulk KKR method
Solves the problem of scattering in three dimensions in an infinite
three-dimensional lattice of axially symmetric scatterers.
Three-dimensional variant of the two-dimensional bulk KKR method with identical
material and multilayered sphere options.
Regarding the band-structure calculations, the bulk KKR method is much
more robust and stable than the three-dimensional layered KKR method.
The reason is that, unlike
in the latter method, no matrix inversion is involved in calculating the
relevant secular matrix. The program was heavily
employed to generate results in several
my articles, which representative selection is as follows:
PhD theses where my codes have been used
ILCS
Usually, only a single sphere layer comprises dispersive material. ILCS defines
the number of the layer (when numbered from the core outwards). Hence ILCS=1 if
the core material is dispersive, ILCS=LCS if the outermost shell material
is dispersive.
LCS
LCS defines the total number of concentric sphere layers. Sphere
core is also considered to be a layer. Therefore, LCS=1 for a homogeneous spheres,
LCS=2 for a spheres with a single shell, etc. The number LCS can be arbitrary.
LMAX
LMAX is a convergence parameter. It defines the cut-off on the value of the
orbital angular momentum of spherical waves taken into
account in the calculation.
NFIN
The parameter defines the dimension of a corresponding material data file.
NMAT
By changing the parameter NMAT value, you can switch between
plurality of different materials. Currently, the following options are
available: dispersionless dielectric, Drude metal, Ag, Au, ZnS, Cu, Al, and Si.
The option of a perfect conductor (the Dirichlet boundary condition) is activated
by the .TRUE. value of the logical parameter YNPERFCON.
The dielectric constant of metals for wavelengths up to
cca 2 micrometers is determined by a linear interpolation between
experimental data according to Palik's Handbook of Optical Constants of Solids.
For longer wavelengths a Drude fit is used according to parameters by
M. A. Ordal et al, Optical properties of the metals
Al, Co, Cu, Au, Fe, Pb, Ni, Pd, Pt, Ag, Ti, V, and W
in the infrared and far infrared, Appl. Opt. 22, 1099-1119 (1983).
Regarding other materials, see
refractiveindex.info
and
Luxpop Index of refraction values.
NOUT
Defines the output unit.
NP
The parameter specifies the shape of particle from a predefined list.
NSTACK
By changing the value of NSTACK you can switch between different stacking sequences.
Currently 12 in three space dimensions and 2 in two space dimensions.
YNPERFCON
The option of a perfect conductor (the Dirichlet boundary condition) is activated
by the .TRUE. value of the logical parameter YNPERFCON. In all other cases, the value
of YNPERFCON should be .FALSE.
In writing my programs, I have occasionally made use of routines from
Numerical Recipes. For copyright reasons,
these routines are not distributed here. (I have not
received any reply yet from Numerical Recipes Software on my request
for a permission to distribute the routines.)
In the case of the free-space one-dimensional periodic Green function in two space dimensions [see A. Moroz, Exponentially convergent lattice sums, Opt. Lett. 26, 1119-1121 (2001) [pdf]], some additional instructions can be found here.
I would highly appreciate any comments and/or bug reports and informing me of any problems encountered with any of the above codes. If you choose to use them, please send email to
registering as a user; registered users will be notified when updates to the code(s) are made. If you publish results obtained using the above codes, please acknowledge the source of the code.