Here are some examples of calculations you can do with
scuffscatter. Input files and
commandline runscripts for all these examples are included
in the share/scuffem/examples subdirectory of
the scuffem installation.
4a. Mie scattering
We first demonstrate how to use
scuffscatter
to solve the canonical textbook problem of Mie scattering  the
scattering of a plane wave from a dielectric sphere. The files for this
example are in the share/scuffem/examples/SolidSphere subdirectory
of the scuffem source distribution.
We begin by creating a
gmsh
geometry file for the sphere
(Sphere.geo ).
We turn this geometry file into a mesh file by running the following
command:
% gmsh 2 clscale 1.0 Sphere.geo
This produces a file named
Sphere.msh,
which looks like this:
You can adjust the fineness of the surface mesh by varying the
clscale parameter (which stands
for "characteristic length scale"); finer meshes will be
more accurate but will take longer to simulate.
Next we create a
scuffem geometry file
that will tell
scuffscatter about our geometry,
including both the surface mesh and the material properties
(dielectric function) of the sphere. As a first example,
we'll use a dielectric model for silicon carbide
that expresses the relative permittivity as a rational function
of ω; in this case we'll call the geometry file
SiCSphere.scuffgeo.
MATERIAL SiliconCarbide
EpsInf = 6.7;
a0 = 3.32377e28;
a1 = +8.93329e11;
b0 = 2.21677e28;
b1 = 8.93329e11;
Eps(w) = EpsInf * (a0 + i*a1*w + w*w) / ( b0 + i*b1*w + w*w);
ENDMATERIAL
OBJECT TheSphere
MESHFILE Sphere.msh
MATERIAL SiliconCarbide
ENDOBJECT
We create a simple file called
OmegaValues.dat
containing a list of angular frequencies at which to run the scattering problem:
0.010
0.013
...
10.0
(We pause to note one subtlety here: Angular frequencies specified using the
Omega or OmegaFile arguments are interpreted
in units of c / 1 μm = 3•10^{14} rad/sec. These
are natural frequency units to use for problems involving micronsized
objects; in particular, for Mie scattering from a sphere of radius
1 μm, as we are considering here, the numerical value of Omega
is just the quantity kR (wavenumber times radius) known as
the ``size parameter'' in the Mie scattering literature.
In contrast,
when specifying functions of angular frequency like Eps(w) in
MATERIAL...ENDMATERIAL sections of geometry files or in any other
scuffem material description,
the w variable is always interpreted in
units of 1 rad/sec , because these are the units in which
tabulated material properties and functional forms for model dielectric
functions are typically expressed.)
Finally, we'll create a little text file called Args that will
contain a list of commandline options for
scuffscatter;
these will include
(1) a specification of the geometry,
(2) the frequency list,
(3) the name of an output file for the power, force, and torque
(4) the name of a cache file for geometric data (this file doesn't
exist yet, but will be created by our first run of
scuffscatter), and
(5) a specification of the incident field.
geometry SiCSphere.scuffgeo
OmegaFile OmegaValues.dat
PFTFile SiCSphere.PFT
Cache Sphere.cache
pwDirection 0 0 1
pwPolarization 1 0 0
And now we just pipe this little file into the standard input of
scuffscatter:
% scuffscatter < Args
This produces the file SiCSphere.PFT ,
which contains one line per simulated frequency;
each line contains data on the scattered and total power,
the force, and the torque on the particle at that frequency.
(On my fairly standard workstation (8 Xeon E5420 cores, 2.5 GHz),
this calculation takes a few minutes to run. You can monitor its progress
by following the scuffscatter.log file.
Note that, during computationallyintensive operations such as the BEM
matrix assembly, the code should be using all available CPU cores
on your workstation; if you find that this is not the case (for example,
by monitoring CPU usage using
htop) you may need to
reconfigure and recompile with different openmp/pthreads
configuration options.)
Here's a comparison of the
scuffscatter results
with the analytical Mie series, as computed using
this Mathematica script.
[Like most Mie codes, this script computes the absorption and scattering
crosssections, which we multiply by the incoming beam flux
(1/2Z_{0} for a unitstrength plane wave) to get values
for the absorbed and scattered power.]
Now let's redo the calculation for a sphere made of gold
instead of silicon carbide. In this case we will name our
scuffem geometry file
GoldSphere.scuffgeo :
MATERIAL Gold
wp = 1.37e16;
gamma = 5.32e13;
Eps(w) = 1  wp^2 / (w * (w + i*gamma));
ENDMATERIAL
OBJECT TheSphere
MESHFILE Sphere.msh
MATERIAL Gold
ENDOBJECT
Since most of the commandline arguments to
scuffscatter
will be the same as before, we can reuse the
same Args file, with the options
that need to be given new values specified on the
command line:
% scuffscatter geometry GoldSphere.scuffgeo PFTFile GoldSphere.PFT < Args
(Note that we don't have to change the name of the cache file
specified with the Cache option; because we are using
the same surface mesh as before, and because
cached geometric data
in scuffem are independent of material
properties, we can take advantage of geometric data computed during
the earlier run for the silicon carbide sphere.)
Now our data look like this:
In some cases it's useful to look at how the induced surface currents
vary over the surface of the object. Let's rerun the SiC example,
now at just the single angular frequency of ω=0.1, and
ask for a surface current plot.
scuffscatter geometry SiCSphere.scuffgeo Omega 0.1 Cache Sphere.cache pwDirection 0 0 1 pwPolarization 1 0 0 PlotSurfaceCurrents
This produces a file named SiCSphere.0.1.pp , which we can
open in gmsh like this:
gmsh SiCSphere.0.1.pp
4b. Electrostatics of a spherical dielectric shell
For our next trick, we'll consider a spherical shell of dielectric
material illuminated by a plane wave of such low frequency that we
may think of the incident field as a spatially constant DC electric
field. In this case it is easy to obtain an
exact analytical solution of the scattering problem,
which we will reproduce numerically using
scuffscatter.
We will take the outer and inner radii to be
R_{out}=1 and
R_{in}=0.5.
The files for this example are in the
share/scuffem/examples/SphericalShell
subdirectory of the scuffem source distribution.
To represent a spherical shell in
scuffem, we need two surface
meshes, one each for the inner and outer spherical surfaces.
These are described by GMSH mesh files
Sphere_R1P0.msh and Sphere_R0P5.msh
We describe the shell as an inner vacuum sphere embedded in
the outer sphere; the geometry file for this situation is
SphericalShell.scuffgeo:
OBJECT OuterSphere
MESHFILE Sphere_R1P0.msh
MATERIAL CONST_EPS_10
ENDOBJECT
OBJECT InnerSphere
MESHFILE Sphere_R0P5.msh
MATERIAL Vacuum
ENDOBJECT


We will run two separate calculations. First, we will fix the
relative permittivity of the shell at ε^{r}=10
and look at the z component of the electric field at points
on the z axis ranging from the origin (the center of the
concentric spheres) to the exterior medium. We create a file called
LineOfPoints
which lists the Cartesian coordinates of each evaluation point:
0.0 0.0 0.000
0.0 0.0 0.025
0.0 0.0 0.050
...
0.0 0.0 2.000
We will pass this file to scuffscatter using
the EPFile option:
% scuffscatter EPFile LineOfPoints < Args
where the Args file looks like this:
geometry SphericalShell.scuffgeo
cache SphericalShell.cache
omega 0.001
pwDirection 1.0 0.0 0.0
pwPolarization 0.0 0.0 1.0
Note that we choose a frequency low enough to ensure we are well within the
electrostatic limit, and that the constant zdirected electrostatic
field described in the memo above becomes a linearly polarized
plane wave traveling in the x direction.
This run of the code produces files SphericalShell.scattered and
SphericalShell.total . Plotting the 8th vs. the 3rd column of the
latter (real part of E_{z} vs. z, as noted
here) yields
good agreement with the analytical calculation, modulo some funkiness at
points on or near the boundary surfaces which is to be expected in an
SIE/BEM calculation:
Next, we will vary the shell permittivity and look at the
electric field at the center of the shell. In this case
the analytical solution makes the interesting prediction
which we will try to verify numerically.
This calculation is slightly trickier than the last one, because
scuffscatter doesn't offer commandline
options for varying the dielectric constant of an object.
One way around this is to use the
python interface to
scuffem, as discussed
on this page.
Here we will pursue a different solution involving a shell
script that modifies the .scuffgeo file for each
different value of ε we want to simulate.
That script looks like this:
#!/bin/bash
cat EpsValues  while read EPS
do
# copy the .scuffgeo file with EPS_10 replaced by EPS_xx
sed "s/EPS_10/EPS_${EPS}/" SphericalShell.scuffgeo > temp.scuffgeo
# run scuffscatter to get Efield at origin
/bin/rm f CenterPoint.total
/bin/rm f CenterPoint.scattered
scuffscatter geometry temp.scuffgeo EPFile CenterPoint < Args
# extract the zcomponent of the field from the output file
EZ=`awk '{print $8}' CenterPoint.total`
echo "${EPS} ${EZ}" >> EzVsEps.dat
done
(Here EpsValues is a file containing the values of
epsilon that we want to simulate, and
CenterPoint is a file containing just the first line
of the file LineOfPoints for the cartesian coordinates
of the origin.)
The result of the calculation looks like this:
4c. Spatiallyresolved study of planewave transmission
through a infinitearea thin dielectric film
The previous examples dealt with compact scatterers. We'll
next consider an
extended geometry  namely, a thin dielectric film
of finite thickness in the z direction but infinitely
extended in the x and y directions. This is
the same geometry for which we used
scufftransmission to look at
the planewave transmission and reflection coefficients as a
function of frequency in
this example,
but here we'll do a different calculation  namely, we'll
pick a single frequency and look at how the electric
and magnetic fields vary in space, both inside and outside
the thin film.
(The files for this example may be found in the
share/scuffem/examples/ThinFilm
directory of the scuffem installation.)
The mesh file and .scuffgeo file for this geometry
are discussed in the
documentation for scufftransmission.
The geometry consists of a film of thickness T=1μm,
with relative dielectric constant ε^{r}=100,
illuminated from below by a plane wave at normal incidence. (We'll take
the incident field to be linearly polarized with E field
pointing in the x direction.) The lower and
upper surfaces of the film are at z=0 and z=T.
For this geometry it is easy to solve Maxwell's equation directly
to obtain the E and H fields directly at points
below, within, and above the film:
We will try to reproduce this behavior using
scuffscatter. First create
a little text file
(ThinFilm.EvalPoints)
containing the coordinates of a bunch of points on a straight line
passing from below the film to above the film. Then put the following
commandline arguments into a file called ThinFilm_58.args:
geometry ThinFilm_58.scuffgeo
cache ThinFilm_58.scuffcache
omega 1.0
EPFile ThinFilm.EvalPoints
pwDirection 0 0 1
pwPolarization 1 0 0
and pipe it into scuffscatter:
scuffscatter < scuffscatter.args
This produces files named ThinFilm.scattered and
ThinFilm.total .
Plotting the 4th vs the 3rd column
(Re E_{x} vs. z)
as well as the the 12th vs the 3rd column
(Re H_{y} vs. z)
of the .total file yields excellent agreement
with theory:
4d. Diffraction of a gaussian beam by a finite disc
As a final example, we shine a laser beam on a small metallic
disc and observe the resulting diffraction pattern on an
observation screen located behind the disc.
The files for this example are in the
share/scuffem/examples/DiffractionPatterns subdirectory
of the scuffem source distribution.
