9
Localized surface plasmons

Figures in this chapter that were generated by Mathematica routines contained in the Code section are labelled "Mathematica simulation" in the figure caption.  The reader can regenerate all of these graphs automatically by selecting "Evaluate Notebook" under the "Evaluation" menu bar in Mathematica.  If Mathematica asks "Do you want to automatically evaluate all the initialization cells in the notebook?" it is sufficient to answer "No."  If the reader is interested in generating similar graphs for different parameter ranges or different materials, it is generally straightforward to do this as well.  The Code section should first be expanded by double clicking on the down arrow of the rightmost bracket.  This will expand the notebook to display "Setup" and "Code for generating figures" subsubsections.  If the notebook has not yet been evaluated, the "Setup" subsubsection should be executed first by selecting (clicking once on) its rightmost bracket and entering "Shift-Return" on a standard keyboard or just "Enter" on an extended keyboard.  Then the "Code for generating figures" subsection should be expanded by double clicking on the downward arrow on the rightmost bracket of that subsubsection.  This subsubsection can be further expanded, if necessary, by again double clicking on the right bracket with the downward pointing arrow.  After expanding the subsubsections, eventually a list of figure captions is obtained.  By double-clicking on the up arrow of the figure caption, a text section and a code section for that figure are displayed.  The parameters for generating the graph are described in the text section.  The reader simply needs to modify these parameters in the code section, which is right below the text section, and then replot the figure by selecting (clicking once on) the rightmost bracket that encloses the code section for the figure and hitting "Shift-Return" on the keyboard or just "Enter" on an extended keyboard.  The figure should then be replotted using the new parameters.  The third line in the Setup subsection of the Code section is normally "option := optC" for generating color figures, but the line can be changed to "option := optB" for generating black and white figures instead if desired.
This notebook also contains a section with all of the figures in the text, a section with the exercises and code for generating relevant figures, and a section with all the chapter references.  Each section can, as usual, be opened by double clicking on the down arrow of the right bracket.

Code

Setup

ClearAll["Global`*"];
Needs["PlotLegends`"];
option :=optC;   (*  change this option to "optC" for color plots or "optB" for black and white  *)

Plotting Examples

Graphics Styles

fs=12;(*FontSize for LabelStyle*)
is={350,200};(*ImageSize for one figure*)
gr={400,150};(*ImageSize for GraphicsRow*)
gr2={480,180};(*ImageSize for GraphicsRow*)
gg={400,300};(*ImageSize for GraphicsGrid*)
gg2={480,300};(*ImageSize for GraphicsGrid*)
gs=-30; (*GraphicsSpacing for two figures*)
ggs={-20,20};(*GraphicsSpacing for four figures*)
thinthick=Thickness[0.005];
thick=Thickness[0.007]; (*Thickness of lines*)

Graphics color option: optC

optC={
LabelStyle->{FontSize -> fs},
GridLines->Automatic,
Frame->True,
ImageSize->is,
PlotRange->All,
PlotStyle->{
{Red,Thick},
{Green,Dashing[{0.06,0.02}],Thick},{Blue,Dashing[{0.02,0.01}],Thick},
{Magenta,DotDashed,Thick},{Purple,Thick,Dashing[{0.04,0.01}]},{Orange,Thick,Dashing[{0.01,0.01}]}}};

Graphics black option: optB

optB={
GridLines -> Automatic,
Frame -> True,
ImageSize->is,
LabelStyle -> {FontSize -> fs},
PlotRange->{All,All},
PlotStyle->{
{Black,thick},
{Black,thick,Dashing[{0.06,0.02}]},
{Black,thick,Dashing[{0.02,0.01}]},
{Black,thick,DotDashed},
{Black, thick, Dashing[{0.04, 0.01}]},
{Black,thick,Dashing[{0.04, 0.03}]}}};

General two plots using optX

xmin=-π;xmax=π;
ymin=-1;ymax=1;
Fig01=Plot[
{Sin[x],Sin[2x],Sin[3x],Sin[4x]},
{x,xmin,xmax},
PlotRange->{ymin,ymax},
FrameLabel -> {{"x",""},{"sinx", ""}},
Evaluate[optC]]

xmin=-π;xmax=π;
ymin=-1;ymax=1;
Fig02=Plot[
{Sin[x],Sin[2x]},
{x,xmin,xmax},Filling -> Axis,
FillingStyle -> Opacity[0.2],
PlotRange->{ymin,ymax},
FrameLabel -> {{"x",""},{"sinx", ""}},
Evaluate[optC]]

GraphicsRow[{Fig01,Fig02},
ImageSize->gr,
Frame->False,
Spacings->gs]

General four plots using optX

GraphicsGrid[{
  {Fig01,Fig02},
  {Fig01,Fig02}},
 ImageSize -> gg,
 Frame -> False,
 Spacings -> gs]

Manipulate plots

amin=1;
amax=4;
Manipulate[
Plot[
 Sin[a θ],
{θ,0,6},
Frame->True,
FrameLabel->{"θ (deg)","sin(a θ)","Function",StringJoin["a = ",ToString[.01Ceiling[100.a]],""]},ImageSize->.75is,
PlotRange->{-1,1},
GridLines->Automatic,
PlotStyle->{
{Red,thinthick},
{Green,thinthick,Dashing[{0.06,0.02}]},{Blue,thinthick,Dashing[{0.02,0.01}]},
{Magenta,thinthick,DotDashed},
{Purple,thinthick,Dashing[{0.04,0.01}]}}],
{{a,amin,"Phase"},amin,amax}]

Tabulated metallic refractive indices

A note on refractive indices:  There are several sources of tabulated refractive indices for silver, gold and other metals in the literature and on the web.  The Handbook of Optical Constants of Solids edited by Palik [Academic, San Diego, 1998] is probably the most frequently cited source in recent years and is one of the more comprehensive sources.  For the noble metals, the data of Johnson and Christy [Phys. Rev. B, 6 (1972) 4370-4379] are also frequently cited, especially in literature prior to the publication of the Handbook.  Indices from both of these sources have been included for use in simulations in this chapter.  However, the experimental precision in the Johnson and Christy data is somewhat limited and the data in the Handbook are frequently stitched together from several sources over an extended wavelength range.  The result of this is that calculations of extinction coefficients and peak field intensities often exhibit anomalous peaks at certain wavelengths corresponding to small changes in the refractive index from the limited experimental precision or a data stitch.  In order to avoid these unphysical effects in the modeling, it is generally much better to use an analytical model of the refractive indices that is smooth and continuous and has no precision problems.  Therefore, we have also included the anayltical models developed by Rakic et al. [Appl. Opt. 37 (1998) 5271-5283] and generally recommend their use for all modeling.

Refractive indices for gold are taken from D. W. Lynch and W. R. Hunter, "Gold" in Handbook of Optical Constants of Solids, ed. E. D. Palik (Academic, San Diego, 1998) 293-295.  Refractive indices for silver are taken from D. W. Lynch and W. R. Hunter, "Silver" in Handbook of Optical Constants of Solids, ed. E. D. Palik (Academic, San Diego, 1998) 355-357.  Refractive indices for copper are taken from D. W. Lynch and W. R. Hunter, "Copper" in Handbook of Optical Constants of Solids, ed. E. D. Palik (Academic, San Diego, 1998) 284-285.  Refractive indices for aluminum are taken from D. Y. Smith, E. Shiles and M. Inokuti, "The optical properties of metallic aluminum" in Handbook of Optical Constants of Solids, ed. E. D. Palik (Academic, San Diego, 1998) 397-400.

GoldIndices = {{0.2883,1.742,1.9},{0.2952,1.776,1.918},{0.3024,1.812,1.92},{0.310, 1.83, 1.916}, {0.3179, 1.84, 1.904}, {0.3263, 1.824, 1.878}, {0.3351, 1.798, 1.86}, {0.3444, 1.766, 1.846}, {0.3542, 1.74, 1.848}, {0.3647, 1.716, 1.862}, {0.3757, 1.696, 1.906}, {0.3875, 1.674, 1.936}, {0.400, 1.658, 1.956}, {0.4133, 1.636, 1.958}, {0.4275, 1.616, 1.94}, {0.4428, 1.562, 1.904}, {0.4592, 1.426, 1.846}, {0.4769, 1.242, 1.796}, {0.4959, 0.916, 1.84}, {0.5166, 0.608, 2.12}, {0.5391, 0.402, 2.54}, {0.5636, 0.306, 2.88}, {0.5904, 0.236, 2.96}, {0.6199, 0.194, 3.05}, {0.6526, 0.166, 3.15}, {0.6888, 0.16, 3.8}, {0.7293, 0.164, 4.35},{0.7749, 0.174, 4.86}, {0.8266, 0.188, 5.39}, {0.8856, 0.21, 5.88}, {0.9537, 0.236, 6.47}, {1.033, 0.272, 7.07}, {1.127, 0.312, 7.93}, {1.240, 0.372, 8.77}, {1.265, 0.389, 8.09}, {1.291, 0.403, 8.25}, {1.319, 0.419, 8.42}, {1.378, 0.454, 8.77}, {1.409, 0.473, 8.96}, {1.442, 0.493, 9.15}, {1.476, 0.515, 9.36}, {1.512, 0.537, 9.58}, {1.550, 0.559, 9.81}, {1.590, 0.583, 10.1}, {1.631, 0.609, 10.3}, {1.675, 0.636, 10.6}, {1.722, 0.665, 10.9}, {1.771, 0.696, 11.2}, {1.823, 0.73, 11.5}, {1.879, 0.767, 11.9}, {1.937, 0.807, 12.2}, {2.000, 0.85, 12.6}, {2.066, 0.896, 13}, {2.138, 0.947, 13.4}, {2.214, 1.002, 13.9}, {2.296, 1.063, 14.4}, {2.384, 1.13, 14.9}, {2.480, 1.205, 15.5}, {2.583, 1.287, 16.1}, {2.695, 1.379, 16.8}, {2.818, 1.482, 17.5}, {2.952, 1.598, 18.3}};

SilverIndices = {{0.200,1.072,1.24},{0.20327869,1.098,1.26},{0.20666667,1.125,1.27},{0.2137931,1.173,1.29},{0.22142857,1.208,1.3},{0.22962963,1.238,1.31},{0.23846154,1.265,1.33},{0.248,1.298,1.35},{0.25306122,1.32,1.35},{0.25833333,1.343,1.35},{0.26382979,1.372,1.35},{0.26956522,1.404,1.33},{0.27555556,1.441,1.31},{0.28181818,1.476,1.26},{0.28837209,1.502,1.19},{0.2952381,1.519,1.08},{0.29879518,1.522,0.992},{0.30243902,1.496,0.882},{0.30617284,1.432,0.766},{0.310,1.323,0.647},{0.31155779,1.246,0.586},{0.31392405,1.149,0.54},{0.31552163,1.044,0.514},{0.31794872,0.932,0.504},{0.31958763,0.815,0.526},{0.32207792,0.708,0.565},{0.32375979,0.616,0.609},{0.32631579,0.526,0.663},{0.33066667,0.371,0.813},{0.33243968,0.321,0.902},{0.33513514,0.294,0.986},{0.33972603,0.259,1.12},{0.34444444,0.238,1.24},{0.35428571,0.209,1.44},{0.36470588,0.186,1.61},{0.37575758,0.2,1.67},{0.3875,0.192,1.81},{0.400,0.173,1.95},{0.41333333,0.173,2.11},{0.42758621,0.16,2.26},{0.44285714,0.157,2.4},{0.45925926,0.144,2.56},{0.47692308,0.132,2.72},{0.496,0.13,2.88},{0.51666667,0.13,3.07},{0.53913043,0.129,3.25},{0.56363636,0.12,3.45},{0.59047619,0.121,3.66},{0.620,0.131,3.88},{0.65263158,0.14,4.15},{0.68888889,0.14,4.44},{0.72941176,0.148,4.74},{0.775,0.143,5.09},{0.82666667,0.145,5.5},{0.88571429,0.163,5.95},{0.95384615,0.198,6.43},{1.03333333,0.226,6.99},{1.12727273,0.251,7.67},{1.240,0.329,8.49},{1.26530612,0.375,7.78},{1.29166667,0.383,7.92},{1.31914894,0.392,8.06},{1.34782609,0.401,8.21},{1.37777778,0.411,8.37},{1.40909091,0.421,8.37},{1.44186047,0.431,8.7},{1.47619048,0.442,8.88},{1.51219512,0.455,9.08},{1.550,0.469,9.32},{1.58974359,0.485,9.57},{1.63157895,0.501,9.84},{1.67567568,0.519,10.1},{1.72222222,0.537,10.4},{1.77142857,0.557,10.7},{1.87878788,0.6,11.4},{1.9375,0.624,11.8},{2.000,0.65,12.2}};

CopperIndices={{0.2952,1.42,1.64},{0.3099,1.34,1.72},{0.3263,1.34,1.81},{0.3444,1.31,1.87},{0.3646,1.27,1.95},{0.3874,1.18,2.21},{0.4133,1.18,2.21},{0.4428,1.17,2.36},{0.4768,1.15,2.5},{0.5166,1.12,2.6},{0.539,1.04,2.59},{0.5635,0.826,2.6},{0.5904,0.468,2.81},{0.6199,0.272,3.24},{0.6525,0.214,3.67},{0.6702,0.215,3.86},{0.688,0.213,4.05},{0.7084,0.214,4.24},{0.7293,0.223,4.43},{0.8265,0.26,5.26},{1.240,0.433,8.46},{1.265,0.496,6.78},{1.291,0.505,6.92},{1.319,0.515,7.06},{1.348,0.525,7.21},{1.378,0.536,7.36},{1.409,0.547,7.53},{1.442,0.559,7.7},{1.476,0.572,7.88},{1.512,0.586,8.06},{1.550,0.606,8.26},{1.589,0.627,8.47},{1.631,0.649,8.69},{1.675,0.672,8.92},{1.722,0.697,9.16},{1.771,0.723,9.41},{1.823,0.752,9.68},{1.878,0.782,9.97},{1.937,0.815,10.3},{2.000,0.85,10.6},{2.066,0.89,11},{2.138,0.933,11.3},{2.214,0.98,11.8},{2.296,1.03,12.2},{2.384,1.09,12.7},{2.480,1.15,13.2},{2.583,1.22,13.7}};

AluminumIndices={{0.300,0.276,3.61},{0.310,0.294,3.74},{0.3179,0.31,3.84},{0.3263,0.326,3.95},{0.3351,0.344,4.06},{0.3444,0.364,4.17},{0.350,0.375,4.24},{0.3542,0.385,4.3},{0.3647,0.407,4.43},{0.3757,0.432,4.56},{0.3875,0.46,4.71},{0.400,0.49,4.86},{0.4133,0.523,5.02},{0.4275,0.558,5.2},{0.4428,0.598,5.38},{0.450,0.618,5.47},{0.4592,0.644,5.58},{0.4769,0.695,5.8},{0.4959,0.755,6.03},{0.500,0.769,6.08},{0.5061,0.789,6.15},{0.5166,0.826,6.28},{0.5276,0.867,6.42},{0.5391,0.912,6.55},{0.550,0.958,6.69},{0.551,0.963,6.7},{0.5636,1.02,6.85},{0.5767,1.08,7},{0.5904,1.15,7.15},{0.600,1.2,7.26},{0.6048,1.22,7.31},{0.6199,1.3,7.48},{0.6358,1.39,7.65},{0.650,1.47,7.79},{0.6526,1.49,7.82},{0.6702,1.6,8.01},{0.6888,1.74,8.21},{0.700,1.83,8.31},{0.7085,1.91,8.39},{0.7293,2.14,8.57},{0.750,2.4,8.62},{0.7514,2.41,8.62},{0.7749,2.63,8.6},{0.7999,2.8,8.45},{0.825,2.75,8.31},{0.8266,2.74,8.31},{0.850,2.61,8.22},{0.8551,2.58,8.21},{0.875,2.38,8.18},{0.8856,2.24,8.21},{0.900,2.06,8.3},{0.9184,1.86,8.44},{0.925,1.77,8.49},{0.950,1.49,8.88},{0.9537,1.47,8.95},{0.9919,1.37,9.49},{1.000,1.35,9.58},{1.033,1.26,10},{1.078,1.21,10.6},{1.127,1.2,11.2},{1.181,1.21,11.8},{1.200,1.21,12},{1.240,1.21,12.5},{1.305,1.23,13.2},{1.378,1.26,14},{1.459,1.33,14.9},{1.500,1.38,15.4},{1.550,1.44,16},{1.653,1.59,17.1},{1.771,1.77,18.3},{1.907,1.99,19.8},{2.000,2.15,20.7}};

GoldIndices2 and SilverIndices2 are from Johnson and Christy, "Optical constants of the noble metals," Phys. Rev. B, 6 (1972) 4370-4379.  Note that at short wavelengths there is a substantial difference between the refractive indices of Ag in this list and that from the Handbook of Optical Constants of Solids.  When comparing to examples in the literature it is important to use the same source for the optical constants!

GoldIndices2 = {{0.18788, 1.28, 1.188}, {0.19165, 1.32, 1.203}, {0.19528, 1.34, 1.226}, {0.19936, 1.33, 1.251}, {0.20328, 1.33, 1.277}, {0.20736, 1.3, 1.304}, {0.21197, 1.3, 1.35}, {0.2164, 1.3, 1.387}, {0.22143, 1.3, 1.427}, {0.22628, 1.31, 1.46}, {0.23134, 1.3, 1.497}, {0.23709, 1.32, 1.536}, {0.24266, 1.32, 1.577}, {0.2490, 1.33, 1.631}, {0.25514, 1.33, 1.688}, {0.2616, 1.35, 1.749}, {0.26898, 1.38, 1.803}, {0.27617, 1.43, 1.847}, {0.2844, 1.47, 1.869}, {0.29245, 1.49, 1.878}, {0.30097, 1.53, 1.889}, {0.31078, 1.53, 1.893}, {0.32041, 1.54, 1.898}, {0.33155, 1.48, 1.883}, {0.34254, 1.48, 1.871}, {0.35429, 1.5, 1.866}, {0.36795, 1.48, 1.895}, {0.38154, 1.46, 1.933}, {0.39744, 1.47, 1.952}, {0.41333, 1.46, 1.958}, {0.43056, 1.45, 1.948}, {0.45091, 1.38, 1.914}, {0.47148, 1.31, 1.849}, {0.4960, 1.04, 1.833}, {0.52101, 0.62, 2.081}, {0.54867, 0.43, 2.455}, {0.58216, 0.29, 2.863}, {0.61692, 0.21, 3.272}, {0.65957, 0.14, 3.697}, {0.70455, 0.13, 4.103}, {0.7561, 0.14, 4.542}, {0.82119, 0.16, 5.083}, {0.89209, 0.17, 5.663}, {0.98413, 0.22, 6.35}, {1.08772, 0.27, 7.15}, {1.21569, 0.35, 8.145}, {1.39326, 0.43, 9.519}, {1.61039, 0.56, 11.21}, {1.9375, 0.92, 13.78}};

SilverIndices2={{0.18788,1.07,1.212},{0.191.65,1.1,1.232},{0.19528,1.12,1.255},{0.19936,1.14,1.277},{0.20328,1.15,1.296},{0.20736,1.18,1.312},{0.21197,1.2,1.325},{0.21640,1.22,1.336},{0.22143,1.25,1.342},{0.22628,1.26,1.344},{0.23134,1.28,1.357},{0.23709,1.28,1.367},{0.24266,1.3,1.378},{0.24900,1.31,1.389},{0.25514,1.33,1.393},{0.26160,1.35,1.387},{0.26898,1.38,1.372},{0.27617,1.41,1.331},{0.28440,1.41,1.264},{0.29245,1.39,1.161},{0.30097,1.34,0.964},{0.31078,1.13,0.616},{0.32041,0.81,0.392},{0.33155,0.17,0.829},{0.34254,0.14,1.142},{0.35429,0.1,1.419},{0.36795,0.07,1.657},{0.38154,0.05,1.864},{0.39744,0.05,2.07},{0.41333,0.05,2.275},{0.43056,0.04,2.462},{0.45091,0.04,2.657},{0.47148,0.05,2.869},{0.49600,0.05,3.093},{0.52101,0.05,3.324},{0.54867,0.06,3.586},{0.58216,0.05,3.858},{0.61692,0.06,4.152},{0.65957,0.05,4.483},{0.70455,0.04,4.838},{0.75610,0.03,5.242},{0.82119,0.04,5.727},{0.89209,0.04,6.312},{0.98413,0.04,6.992},{1.08772,0.04,7.795},{1.21569,0.09,8.828},{1.39326,0.13,10.1},{1.61039,0.15,11.85},{1.9375,0.24,14.08}};

Routines to compute refractive indices of various metals

Drude[λ_]:=Module[{n,ω,ϵ,c,h,ωp,γ,ϵ∞},
c=2.99792456×(10^14);          (*  speed of light in μm/s  *)
h=4.1359×(10^(-15));                (*  Plank's constant in eV s  *)
ϵ∞=4;                                     (*  background dielectric constant  *)
ωp=15;                                       (*  ωp in eV  *)
γ = 0.1;                               (*  damping coefficient in eV  *)
ω = 7.7905/λ;                     (*  ang. frequency in eV when λ is in μm  *)
ϵ=ϵ∞-(ωp^2)/(ω^2 + I γ ω);
n = (Sqrt[ϵ]);
(*  Print["λ=",λ," n=",n];  *)
Return[n]];

(* from Rakic et al., Appl. Opt. vol. 37 (1998) 5271   *)
LorentzDrudeAg[λ_]:=Module[{n,ω,f,ωω,ϵ,c,ωp,γ,ϵ∞},
f={0.065,0.124,0.011,0.840,5.646};
γ={3.886,0.452,0.065,0.916,2.419};
ωω={0.816,4.481,8.185,9.083,20.29};
ωp=9.01;
ω = 1.24/λ;                 (*  ang. frequency in eV when λ is in μm  *)
Do[If[i==0,ϵ=1-(0.845(ωp^2)/(ω^2+I 0.048ω )),
ϵ=ϵ+f[[i]] ((ωp)^2)/((ωω[[i]])^2-ω^2-I ω γ[[i]])],{i,0,5}];
n = (Sqrt[ϵ]);
Return[n]];


(* from Rakic et al., Appl. Opt. vol. 37 (1998) 5271   *)
LorentzDrudeAu[λ_]:=Module[{n,ω,f,ωω,ϵ,c,ωp,γ,ϵ∞},
f={0.024,0.010,0.071,0.601,4.384};
γ={0.241,0.345,0.870,2.494,2.214};
ωω={0.415,0.830,2.969,4.304,13.32};
ωp=9.03;
ω = 1.24/λ;                 (*  ang. frequency in eV when λ is in μm  *)
Do[If[i==0,ϵ=1-(0.760(ωp^2)/(ω^2+I 0.053ω )),
ϵ=ϵ+f[[i]] ((ωp)^2)/((ωω[[i]])^2-ω^2-I ω γ[[i]])],{i,0,5}];
n = (Sqrt[ϵ]);
Return[n]];

(* from Rakic et al., Appl. Opt. vol. 37 (1998) 5271   *)
LorentzDrudeCu[λ_]:=Module[{n,ω,f,ωω,ϵ,c,ωp,γ,ϵ∞},
f={0.061,0.104,0.723,0.638};
γ={0.378,1.056,3.213,4.305};
ωω={0.291,2.957,5.300,11.18};
ωp=10.83;
ω = 1.24/λ;                 (*  ang. frequency in eV when λ is in μm  *)
Do[If[i==0,ϵ=1-(0.575(ωp^2)/(ω^2+I 0.030ω )),
ϵ=ϵ+f[[i]] ((ωp)^2)/((ωω[[i]])^2-ω^2-I ω γ[[i]])],{i,0,4}];
n = (Sqrt[ϵ]);
Return[n]];

(* from Rakic et al., Appl. Opt. vol. 37 (1998) 5271   *)
LorentzDrudeAl[λ_]:=Module[{n,ω,f,ωω,ϵ,c,ωp,γ,ϵ∞},
f={0.227,0.050,0.166,0.030};
γ={0.333,0.312,1.351,3.382};
ωω={0.162,1.544,1.808,3.473};
ωp=14.98;
ω = 1.24/λ;                 (*  ang. frequency in eV when λ is in μm  *)
Do[If[i==0,ϵ=1-(0.523(ωp^2)/(ω^2+I 0.047ω )),
ϵ=ϵ+f[[i]] ((ωp)^2)/((ωω[[i]])^2-ω^2-I ω γ[[i]])],{i,0,4}];
n = (Sqrt[ϵ]);
Return[n]];

Gold[λ_]:=Module[{n,n1,n2,K1,K2,λ1,λ2},
If[λ<GoldIndices[[1,1]],Print[" Input wavelength is too short."];Abort[]];
If[λ>GoldIndices[[Length[GoldIndices],1]],Print[" Input wavelength is too long."];Abort[]];
Catch[Do[λ2=GoldIndices[[i,1]];
n2=GoldIndices[[i,2]];
K2=GoldIndices[[i,3]];
If[λ==λ2,n=n2+I K2;Throw[n],
If[λ<λ2,λ1=GoldIndices[[i-1,1]];
n1=GoldIndices[[i-1,2]];
K1=GoldIndices[[i-1,3]];
n=n1 + (n2-n1)(λ-λ1)/(λ2-λ1) + I(K1 + (K2-K1)(λ-λ1)/(λ2-λ1));
Throw[n]]],{i,Length[GoldIndices]}]];
Return[n]];

Gold2[λ_]:=Module[{n,n1,n2,K1,K2,λ1,λ2},
If[λ<GoldIndices2[[1,1]],Print[" Input wavelength is too short."];Abort[]];
If[λ>GoldIndices2[[Length[GoldIndices2],1]],Print[" Input wavelength is too long."];Abort[]];
Catch[Do[λ2=GoldIndices2[[i,1]];
n2=GoldIndices2[[i,2]];
K2=GoldIndices2[[i,3]];
If[λ==λ2,n=n2+I K2;Throw[n],
If[λ<λ2,λ1=GoldIndices2[[i-1,1]];
n1=GoldIndices2[[i-1,2]];
K1=GoldIndices2[[i-1,3]];
n=n1 + (n2-n1)(λ-λ1)/(λ2-λ1) + I(K1 + (K2-K1)(λ-λ1)/(λ2-λ1));
Throw[n]]],{i,Length[GoldIndices2]}]];
Return[n]];

Silver[λ_]:=Module[{n,n1,n2,K1,K2,λ1,λ2},
If[λ<SilverIndices[[1,1]],Print[" Input wavelength is too short."];Abort[]];
If[λ>SilverIndices[[Length[SilverIndices],1]],Print[" Input wavelength is too long."];Abort[]];
Catch[Do[λ2=SilverIndices[[i,1]];
n2=SilverIndices[[i,2]];
K2=SilverIndices[[i,3]];
If[λ==λ2,n=n2+I K2;Throw[n],
If[λ<λ2,λ1=SilverIndices[[i-1,1]];
n1=SilverIndices[[i-1,2]];
K1=SilverIndices[[i-1,3]];
n=n1 + (n2-n1)(λ-λ1)/(λ2-λ1) + I(K1 + (K2-K1)(λ-λ1)/(λ2-λ1));
Throw[n]]],{i,Length[SilverIndices]}]];
Return[n]];

Silver2[λ_]:=Module[{n,n1,n2,K1,K2,λ1,λ2},
If[λ<SilverIndices2[[1,1]],Print[" Input wavelength is too short."];Abort[]];
If[λ>SilverIndices2[[Length[SilverIndices2],1]],Print[" Input wavelength is too long."];Abort[]];
Catch[Do[λ2=SilverIndices2[[i,1]];
n2=SilverIndices2[[i,2]];
K2=SilverIndices2[[i,3]];
If[λ==λ2,n=n2+I K2;Throw[n],
If[λ<λ2,λ1=SilverIndices2[[i-1,1]];
n1=SilverIndices2[[i-1,2]];
K1=SilverIndices2[[i-1,3]];
n=n1 + (n2-n1)(λ-λ1)/(λ2-λ1) + I(K1 + (K2-K1)(λ-λ1)/(λ2-λ1));
Throw[n]]],{i,Length[SilverIndices2]}]];
Return[n]];

Copper[λ_] := Module[{n, n1, n2, K1, K2, λ1, λ2},
   If[λ < CopperIndices[[1, 1]], Print[" Input wavelength is too short."]; Abort[]];
   If[λ > CopperIndices[[Length[CopperIndices], 1]], Print[" Input wavelength is too long."]; Abort[]];
   Catch[Do[λ2 = CopperIndices[[i, 1]];
     n2 = CopperIndices[[i, 2]];
     K2 = CopperIndices[[i, 3]];
     If[λ == λ2, n = n2 + I K2; Throw[n],
      If[λ < λ2, λ1 = CopperIndices[[i - 1, 1]];
       n1 = CopperIndices[[i - 1, 2]];
       K1 = CopperIndices[[i - 1, 3]];
       n = n1 + (n2 - n1) (λ - λ1)/(λ2 - λ1) + I (K1 + (K2 - K1) (λ - λ1)/(λ2 - λ1));
       Throw[n]]], {i, Length[CopperIndices]}]];
   Return[n]];

Aluminum[λ_] := Module[{n, n1, n2, K1, K2, λ1, λ2},
   If[λ < AluminumIndices[[1, 1]], Print[" Input wavelength is too short."]; Abort[]];
   If[λ > AluminumIndices[[Length[AluminumIndices], 1]], Print[" Input wavelength is too long."]; Abort[]];
   Catch[Do[λ2 = AluminumIndices[[i, 1]];
     n2 = AluminumIndices[[i, 2]];
     K2 = AluminumIndices[[i, 3]];
     If[λ == λ2, n = n2 + I K2; Throw[n],
      If[λ < λ2, λ1 = AluminumIndices[[i - 1, 1]];
       n1 = AluminumIndices[[i - 1, 2]];
       K1 = AluminumIndices[[i - 1, 3]];
       n = n1 + (n2 - n1) (λ - λ1)/(λ2 - λ1) + I (K1 + (K2 - K1) (λ - λ1)/(λ2 - λ1));
       Throw[n]]], {i, Length[AluminumIndices]}]];
   Return[n]];

Mie scattering from a sphere

This program computes the extinction and scattering coefficients for a sphere from Mie theory.  It has been adapted from the Fortran code in C. F. Bohren and D. R. Huffman, Absorption and Scattering of Light by Small Particles (Wiley, New York, 1983) Appendix A.  The SphereScattering function can be used to compute extinction, scattering, backscattering, and near field efficiencies as a function of wavelength or sphere radius.  It can also be used to compute polar diagrams of the scattered light intensity.  For calculations as a function of wavelength or sphere radius, the function returns four lists of (x, y) coordinates for plotting.  For each list, x is the wavelength in microns. y is the scattering efficiency in the first list, the extinction efficiency in the second list, the backscattering efficiency in the third list, and the near field efficiency in the fourth list.  For polar scattering diagrams, the function also returns four lists.  Let (Subscript[i, ⊥]) be the scattered irradiance per unit incident irradiance when the incident light is polarized perpendicular to the scattering plane, and let (Subscript[i, ∥]) be the scattered irradiance per unit incident irradiance when the incident light is polarized parallel to the scattering plane.  The first list is the set of coordinates (θ, (Subscript[i, ⊥])) as θ goes from 0° to 180°.  The second list is the same set of data expressed in a Cartesian frame, ((Subscript[i, ⊥]) cos(θ), (Subscript[i, ⊥])sin(θ)).  The third list is the set of coordinates (θ, (Subscript[i, ∥])) as θ goes from 0° to 180°, and the fourth list is the same set of data expressed in a Cartesian frame, ((Subscript[i, ∥]) cos(θ), (Subscript[i, ∥])sin(θ)).

The parameters that are set by the user before a calculation are as follows:
	(Subscript[λ, 1]):  initial wavelength in microns	
	(Subscript[λ, 2]):  final wavelength in microns
	(Subscript[λ, inc]): increment of wavelength in microns
	(Subscript[θ, inc]): increment of scattering angle in degrees
	(Subscript[r, 1]): initial radius of sphere in microns
	(Subscript[r, 2]): final radius of sphere in microns
	(Subscript[r, inc]): increment of radius of sphere in microns
	(Subscript[n, 1]): refractive index of outer medium
	(Subscript[n, 2]): refractive index of sphere

If (Subscript[λ, 1])≠(Subscript[λ, 2]) then the function computes efficiencies vs. wavelength.  If (Subscript[λ, 1])=(Subscript[λ, 2]), but (Subscript[r, 1])≠(Subscript[r, 2]), then the function computes efficiencies vs. sphere radius.  Finally, if (Subscript[λ, 1])=(Subscript[λ, 2]) and (Subscript[r, 1])=(Subscript[r, 2]), the function computes the polar diagrams wiht an angle increment of (Subscript[θ, inc]).
	The user can either set the refractive index of the sphere, (Subscript[n, 2]), to a specific (complex) value, or the user can set (Subscript[n, 2]) to a negative integer which tells the function to insert the wavelength-dependent refractive index for a particular metal as follows: 
	-1: complex index for gold interpolated from values in Palik 
	-2: complex index for silver interpolated from values in Palik  
	-3: complex index for copper interpolated from values in Palik 
	-4: complex index for aluminum interpolated from values in Palik
	-5: complex index of a Drude metal
	-6: complex index for gold from Drude-Lorentz fit from Rakik et al.
	-7: complex index for silver from Drude-Lorentz fit from Rakik et al.
	-8: complex index for copper from Drude-Lorentz fit from Rakik et al.
	-9: complex index for aluminum from Drude-Lorentz fit from Rakik et al.

mie[λ_,r_,n1_,n2_,nang_ ]:=Module[{s1,s2,Qsct,Qext,Qbck,Qnf,nsph,x,y,nstop,nang1,nmx,amu,d,rn,π0,π1,ψ0,ψ1,ψ2,ψ3,χ0,χ1,χ2,χ3,ξ0,ξ1,ξ2,ξ3,h0,h1,h2,h3,fn,jj,pi,τ,p,t,an,bn,m},
If[nang==1,nang1=2,nang1=nang];
nsph=Which[n2==-1,Gold[λ],n2==-2,Silver[λ],n2==-3,Copper[λ],n2==-4,Aluminum[λ],n2==-5,Drude[λ],n2==-6,LorentzDrudeAu[λ],n2==-7,LorentzDrudeAg[λ],n2==-8,LorentzDrudeCu[λ],n2==-9,LorentzDrudeAl[λ],n2==-10,Gold2[λ],n2==-11,Silver2[λ],True,n2];
m=nsph/n1;
x=2π r n1/λ;
y=m x;

nstop=Round[x+4 (x^(1/3))+2];
nmx=Round[Max[nstop,Abs[y]]]+15;
amu=Table[Cos[π(j-1)/(2(nang1-1))],{j,nang1}];
d[nmx]=0;
Table[rn=nmx-n+1;d[nmx-n]=rn/y - 1/(d[nmx-n+1] + rn/y),{n,nmx-1}];
Do[π0[j]=0;π1[j]=1,{j,nang1}];
Do[s1[j]=0;s2[j]=0,{j,2 nang1-1}];
ψ0=Cos[x];ψ1=Sin[x];ψ2=ψ1/x - ψ0;
χ0=-Sin[x];χ1=Cos[x];χ2=χ1/x - χ0;
ξ0=ψ0-I χ0;ξ1=ψ1-I χ1;ξ2=ψ2-I χ2;          (*spherical Hankel fns of 1st kind*)
h0=(ψ0+I χ0)/x;h1=(ψ1+I χ1)/x;h2=(ψ2+I χ2)/x;  (*Hankel fns of 2nd kind*)
Qsct=0;
Qnf=0;
Do[fn=(2 n+1)/(n(n+1));
ψ3=(2 (n+1)-1) ψ2/x - ψ1;
χ3=(2 (n+1)-1) χ2/x - χ1;
ξ3=ψ3-I χ3;
h3=(ψ3+I χ3)/x;
an=((d[n]/m + n/x)ψ2 - ψ1)/((d[n]/m + n/x)ξ2 - ξ1);
bn=((m d[n] + n/x)ψ2 - ψ1)/((m d[n] + n/x)ξ2 - ξ1);
Qsct=Qsct+(2 n+1)(Abs[an]^2+Abs[bn]^2);
Qnf=Qnf+2(Abs[an]^2 ((n+1) Abs[h1]^2 + n Abs[h3]^2) + (2n+1)Abs[bn]^2 Abs[h2]^2);
Do[jj=Round[2 nang1-j];
pi[j]=π1[j];
τ[j]=n amu[[j]] pi[j] - (n+1)π0[j];
p=(-1)^(n-1);
s1[j]=s1[j] + fn (an pi[j] + bn τ[j]);
t=(-1)^n;
s2[j]=s2[j]+fn(an τ[j]+ bn pi[j]);
If[j!=jj,s1[jj]=s1[jj] + fn (an pi[j] p + bn τ[j] t);
s2[jj]=s2[jj] + fn (an τ[j] t + bn pi[j] p)],{j,nang1}];
ψ0=ψ1;ψ1=ψ2;ψ2=ψ3;
χ0=χ1;χ1=χ2;χ2=χ3;
ξ0=ξ1;ξ1=ξ2;ξ2=ξ3;
h0=h1;h1=h2;h2=h3;
Do[π1[j]=((2 n+1)/n) amu[[j]] pi[j]-(n+1) π0[j]/n;
π0[j]=pi[j],{j,nang1}],{n,nstop}];
Qsct=2Qsct/x^2;
Qext=4 Re[s1[1]]/x^2;
Qbck=4(Abs[s1[2 nang1-1]]^2/x^2);
Return[{Qsct,Qext,Qbck,Qnf,s1,s2}]];

SphereScattering[λ1_,λ2_,λinc_,r1_,r2_,rinc_,n1_,n2_,θinc_]:=Module[{Qs,Qe,Qb,Qn,Qsct,Qbck,Qext,Qnf,p1,pp1,p2,pp2,p3,pp3,rinc1,nr,r,nang,s1,s2,s11norm,θ,s11,s12,s33,s34,perp,perp1,par,par1,pol,result,λinc1,nλ,λ},
If[λ2≠λ1,            (*  find wavelength dependence  *)
If[λ2>λ1,λinc1=Abs[λinc],λinc1=-Abs[λinc] ];
nλ=Round[Abs[(λ2-λ1)/λinc1]+1];
result=Table[λ=λ1+(j-1)λinc1;
mie[λ,r1,n1,n2,1 ],{j,nλ}];
Qs=Table[{λ1+(j-1)λinc1,result[[j,1]]},{j,nλ}];
Qe=Table[{λ1+(j-1)λinc1,result[[j,2]]},{j,nλ}];
Qb=Table[{λ1+(j-1)λinc1,result[[j,3]]},{j,nλ}];
Qn=Table[{λ1+(j-1)λinc1,result[[j,4]]},{j,nλ}];
Return[{Qs,Qe,Qb,Qn}],

If[r1≠r2,            (*  find radius dependence  *)
If[r2>r1,rinc1=Abs[rinc],rinc1=-Abs[rinc] ];
nr = Round[Abs[(r2-r1)/rinc1]+1];
result=Table[r=r1+(j-1)rinc1;
mie[λ1,r,n1,n2,1 ] ,{j,nr}];
Qs=Table[{r1+(j-1)rinc1,result[[j,1]]},{j,nr}];
Qe=Table[{r1+(j-1)rinc1,result[[j,2]]},{j,nr}];
Qb=Table[{r1+(j-1)rinc1,result[[j,3]]},{j,nr}];
Qn=Table[{r1+(j-1)rinc1,result[[j,4]]},{j,nr}];
Return[{Qs,Qe,Qb,Qn}],

nang = Round[Abs[90/θinc]+1];          (*  find angle dependence  *)
{Qsct,Qext,Qbck,Qnf,s1,s2}=mie[λ1,r1,n1,n2,nang ] ;
Print["λ = ",λ1," Qsct=",Qsct, " Qext=",Qext, " Qbck=",Qbck];
s11norm=(((Abs[s1[1]])^2)+ ((Abs[s2[1]])^2) )/2;
Do[θ=θinc(j-1);
s11=(Abs[s1[j]]^2+Abs[s2[j]]^2)/2;
s12=(Abs[s2[j]]^2-Abs[s1[j]]^2)/2;
perp[j]={θ,s11-s12};
perp1[j]={(s11-s12)Cos[π θ/180.0],(s11-s12)Sin[π θ/180.0]};
par[j]={θ,s11+s12};
par1[j]={(s11+s12)Cos[π θ/180.0],(s11+s12)Sin[π θ/180.0]};
pol=-s12/s11;
s33=Re[s2[j] Conjugate[s1[j]]]/s11;
s34=Im[s2[j] Conjugate[s1[j]]]/s11;
s11=s11/s11norm,{j,2 nang-1}] ;
p1=Table[perp[j],{j,2 nang-1}];   (* perp scattering in polar coords *)
pp1=Table[perp1[j],{j,2 nang-1}]; (* perp scattering in x-y coords  *)
p2=Table[par[j],{j,2 nang-1}];      (* par scattering in polar coords *)
pp2=Table[par1[j],{j,2 nang-1}];  (* par scattering in x-y coords  *)
p3=Table[0.5(perp[j]+par[j]),{j,2 nang-1}];      (* perp + par scattering in polar coords *)
pp3=Table[0.5(perp1[j]+par1[j]),{j,2 nang-1}];  (* perp + par scattering in x-y coords  *)
Return[{p1,pp1,p2,pp2,p3,pp3}]]]];

Near field of a sphere

The MieNearField function computes the magnitude of the electric field from Mie theory in the vicinity of a sphere for an incident plane wave polarized along the x axis with unit amplitude.  The function returns five lists of (x, y) coordinates.  The x coordinate in all lists is either the wavelength or one of the x, y, or z coordinates of the computed field amplitude in units of microns.  The y coordinate in the first list is|(Subscript[E, sx])|, the magnitude of the x component of the scattered field.   The y coordinate in the second list is|(Subscript[E, sy])| and the y coordinate in the third list is|(Subscript[E, sz])|.  The y coordinate in the fourth list is the magnitude of the scattered field, |(Subscript[E, s])|=(Sqrt[|(Subscript[E, sx])((|)^2)+|(Subscript[E, sy])((|)^2)+|(Subscript[E, sz])((|)^2)]) .  Finally, the y coordinate in the fifth list is the magnitude of the total field including incident plus scattered fields, |(Subscript[E, t])|=(Sqrt[|(Subscript[E, sx])+(Subscript[E, i])((|)^2)+|(Subscript[E, sy])((|)^2)+|(Subscript[E, sz])((|)^2)]).   The parameters that are set by the user before a calculation are as follows:

	(Subscript[λ, 1]): initial wavelength in microns
	(Subscript[λ, 2]): final wavelength in microns
	(Subscript[λ, inc]): increment of wavelength in microns	
	r: radius of sphere in microns
	(Subscript[n, 1]): refractive index of outer medium
	(Subscript[n, 2]): refractive index of the sphere
	(Subscript[x, 1]): initial x position in microns
	(Subscript[x, 2]): final x position in microns
	(Subscript[x, inc]): increment of x position in microns
	(Subscript[y, 1]): initial y position in microns
	(Subscript[y, 2]): final y position in microns
	(Subscript[y, inc]): increment of y position in microns	
	(Subscript[z, 1]): initial z position in microns
	(Subscript[z, 2]): final z position in microns
	(Subscript[z, inc]): increment of z position in microns

	If (Subscript[λ, 1])≠(Subscript[λ, 2]), the function computes the field magnitudes vs. wavelength at the point ((Subscript[x, 1]),(Subscript[y, 1]),(Subscript[z, 1])).  If (Subscript[λ, 1])=(Subscript[λ, 2]), but (Subscript[x, 1])≠(Subscript[x, 2]), the function computes the field magnitudes at wavelength (Subscript[λ, 1]) and at point (x,(Subscript[y, 1]),(Subscript[z, 1])) as x varies from (Subscript[x, 1]) to (Subscript[x, 2]) in increments of (Subscript[x, inc]).  If (Subscript[λ, 1])=(Subscript[λ, 2]) and (Subscript[x, 1])=(Subscript[x, 2]), but (Subscript[y, 1])≠(Subscript[y, 2]), the function computes the field magnitudes at wavelength (Subscript[λ, 1]) and at point ((Subscript[x, 1]),y,(Subscript[z, 1])) as y varies from (Subscript[y, 1]) to (Subscript[y, 2]) in increments of (Subscript[y, inc]).  Finally, if (Subscript[λ, 1])=(Subscript[λ, 2]), (Subscript[x, 1])=(Subscript[x, 2]) and (Subscript[y, 1])=(Subscript[y, 2]), but (Subscript[z, 1])≠(Subscript[z, 2]), the function computes the field magnitudes at wavelength (Subscript[λ, 1]) and at point ((Subscript[x, 1]),(Subscript[y, 1]),z) as z varies from (Subscript[z, 1]) to (Subscript[z, 2]) in increments of (Subscript[z, inc]).
	As with the SphereScattering function, the user can either set the refractive index of the sphere, (Subscript[n, 2]), to a specific (complex) value, or the user can set (Subscript[n, 2]) to a negative integer which tells the function to insert the wavelength-dependent refractive index for a particular metal as follows: 
	-1: complex index for gold interpolated from values in Palik 
	-2: complex index for silver interpolated from values in Palik  
	-3: complex index for copper interpolated from values in Palik 
	-4: complex index for aluminum interpolated from values in Palik
	-5: complex index of a Drude metal
	-6: complex index for gold from Drude-Lorentz fit from Rakik et al.
	-7: complex index for silver from Drude-Lorentz fit from Rakik et al.
	-8: complex index for copper from Drude-Lorentz fit from Rakik et al.
	-9: complex index for aluminum from Drude-Lorentz fit from Rakik et al.

ψ[n_,p_]:=Module[{L,x},
L=n+1/2;
x=BesselJ[L,p];
x Sqrt[π p/2]];

ψd[n_,p_]:=Module[{L,x1,x2},
L=n+1/2;
x1=BesselJ[L,p];
x2=BesselJ[L+1,p];
(x1(n+1)/p-x2)Sqrt[π p/2] ];

χ[n_,p_]:=Module[{L,x},
L=n+1/2;
x=BesselY[L,p];
-x Sqrt[π p/2]];

χd[n_,p_]:=Module[{L,x1,x2},
L=n+1/2;
x1=BesselY[L,p];
x2=BesselY[L+1,p];
-(x1(n+1)/p-x2)Sqrt[π p/2]];

ζ[L_,p_]:=ψ[L,p]-I χ[L,p];

ζd[L_,p_]:=ψd[L,p]-I χd[L,p];

yPLmd[L_,m_,x_]:=-(L+m)(L-m+1)LegendreP[L,m-1,x]-m x LegendreP[L,m,x]/Sqrt[1-(x^2)];

Compute[λ_,a_,nmed_,nsph_,r_,θ_,ϕ_]:=Module[{Qsct,Qext,Qbck,Er,Eθ,Eϕ,Ex,Ey,Ez,Einc,MaxIter,k0,k1,k2,k1r,k2r,d,c,h,j,x,y,sum1,sum2,sum3,aL,bL,cL,dL,eAL,eBL,mAL,mBL,aa,bb,cc,ncore},
ncore=Which[nsph==-1,Gold[λ],nsph==-2,Silver[λ],nsph==-3,Copper[λ],nsph==-4,Aluminum[λ],nsph==-5,Drude[λ],nsph==-6,LorentzDrudeAu[λ],nsph==-7,LorentzDrudeAg[λ],nsph==-8,LorentzDrudeCu[λ],nsph==-9,LorentzDrudeAl[λ],nsph==-10,Gold2[λ],nsph==-11,Silver[λ],True,nsph];
MaxIter=Round[Re[a k1]+4(Re[a k1])^(1/3) +2]+15;
k0=2 π/λ;k1=nmed k0;k2=ncore k0;k1r=k1 r;k2r=k2 r;
d=Abs[((ncore^2)-(nmed^2))/((ncore^2)+2 (nmed^2))];
c=((k1^2))((a^3))d;
h=a k1;j=a k2;
x=Cos[π θ/180];
If[x≥1,x=0.9999999,x];
If[x==0,x=0.0000001,x];
If[x≤-1,x=-0.9999999,x];
y=Sqrt[1-(x^2)];
sum1=0;sum2=0;sum3=0;
Qsct=0;Qext=0;Qbck=0;
Do[aL=(k2 ψd[L,h]ψ[L,j]-k1 ψd[L,j]ψ[L,h])/(k2 ζd[L,h]ψ[L,j]-k1 ψd[L,j]ζ[L,h]);
bL=(k2 ψ[L,h]ψd[L,j]-k1 ψd[L,h]ψ[L,j])/(k2 ζ[L,h]ψd[L,j]-k1 ζd[L,h]ψ[L,j]);
cL=(k2 ζd[L,h]ψ[L,h]-k2 ψd[L,h]ζ[L,h])/(k2 ζd[L,h]ψ[L,j]-k1 ψd[L,j]ζ[L,h]);
dL=(k2 ζ[L,h]ψd[L,h]-k2 ζd[L,h]ψ[L,h])/(k2 ζ[L,h]ψd[L,j]-k1 ζd[L,h]ψ[L,j]);
eBL=((I^(L+1)) (2L+1)/(L(L+1))) aL;
mBL=((I^(L+1)) (2L+1)/(L(L+1))) bL;
eAL=((I^(L+1)) (2L+1)/(L(L+1))) cL;
mAL=((I^(L+1)) (2L+1)/(L(L+1))) dL;
Qsct=Qsct+(2L+1)(Abs[aL]^2 + Abs[bL]^2);
Qext=Qext+(2L+1) Re[aL+bL];
Qbck=Qbck+(2L+1)(-1^L)(aL-bL);
If[Abs[r]≥a,   
(*  outside the sphere  *)
aa=+eBL ζd[L,k1r] yPLmd[L,1,x]+I mBL ζ[L,k1r] LegendreP[L,1,x]/y;
bb=-eBL ζd[L,k1r] LegendreP[L,1,x]/y-I mBL ζ[L,k1r] yPLmd[L,1,x];
cc=-L(L+1)eBL ζ[L,k1r] LegendreP[L,1,x];
aa=aa/k1r;
bb=bb/k1r;
cc=cc/((k1r^2)),
(*  inside the sphere  *)
aa=-eAL ψd[L,k2r] yPLmd[L,1,x]-I mAL ψ[L,k2r] LegendreP[L,1,x]/y;
bb=+eAL ψd[L,k2r] LegendreP[L,1,x]/y+I mAL ψ[L,k2r] yPLmd[L,1,x];
cc=+L(L+1)eAL ψ[L,k2r] LegendreP[L,1,x];
aa=aa/k2r;
bb=bb/k2r;
cc=cc/((k2r^2))];
sum1=sum1+aa;sum2=sum2+bb;sum3=sum3+cc,{L,MaxIter}];
Qsct=2Qsct/(a k1)^2;
Qext=2Qext/(a k1)^2;
Qbck=(Abs[Qbck]/(a k1))^2;
If[r≥a,Einc=Exp[I k1 r Cos[π θ/180]],Einc=0];
Er=+Cos[π ϕ/180]sum3;
Eθ=-Cos[π ϕ/180]sum1;
Eϕ=-Sin[π ϕ/180]sum2;
Ex=Er Sin[π θ/180]Cos[π ϕ/180]+Eθ Cos[π θ/180]Cos[π ϕ/180]-Eϕ Sin[π ϕ/180];
Ey=Er Sin[π θ/180]Sin[π ϕ/180]+Eθ Cos[π θ/180]Sin[π ϕ/180]+Eϕ Cos[π ϕ/180];
Ez=Er Cos[π θ/180]-Eθ Sin[π θ/180];
Return[{Einc,Ex,Ey,Ez,Er,Eθ,Eϕ,Qsct,Qext,Qbck}]];

MieNearField[λ1_,λ2_,λinc_,a_,nmed_,nsph_,x1_,x2_,xinc_,y1_,y2_,yinc_,z1_,z2_,zinc_]:=Module[{ExMag,EyMag,EzMag,EsMag,EtMag,θ,ϕ,z,r,r1,λ,λλinc,inc,rinc,result},
If[λ1≠λ2,    (*  vary wavelength  *)
θ=(180/π) ArcTan[z1,Sqrt[x1^2 + y1^2]];
ϕ=(180/π)ArcTan[x1,y1];
r=Sqrt[x1^2 + y1^2 + z1^2];            
inc=Round[Abs[(λ2-λ1)/λinc]+1];
If[λ2>λ1,λλinc=Abs[λinc],λλinc=-Abs[λinc]];
result=Table[λ=λ1+(i-1)λλinc;
Compute[λ,a,nmed,nsph,r,θ,ϕ],{i,inc}];
ExMag=Table[{λ1+(i-1)λλinc,Abs[result[[i,2]]]},{i,inc}];
EyMag=Table[{λ1+(i-1)λλinc,Abs[result[[i,3]]]},{i,inc}];
EzMag=Table[{λ1+(i-1)λλinc,Abs[result[[i,4]]]},{i,inc}];
EsMag=Table[{λ1+(i-1)λλinc,Sqrt[Abs[result[[i,2]]]^2+Abs[result[[i,3]]]^2+Abs[result[[i,4]]]^2]},{i,inc}];
EtMag=Table[{λ1+(i-1)λλinc,Sqrt[Abs[result[[i,1]]+result[[i,2]]]^2+Abs[result[[i,3]]]^2+Abs[result[[i,4]]]^2]},{i,inc}];
Return[{ExMag,EyMag,EzMag,EsMag,EtMag}],

If[x1≠x2,       (*  vary position on x axis  *)
λ=λ1;θ=90;r1=x1;            
inc=Round[Abs[(x2-x1)/xinc]+1];
rinc=xinc;
If[x2>x1,rinc=Abs[rinc],rinc=-Abs[rinc]];
result=Table[r=r1+(i-1)rinc;
If[r≥0,ϕ=0,ϕ=180];
r=Abs[r];
If[r==0,r=0.0001];
Compute[λ,a,nmed,nsph,r,θ,ϕ],{i,inc}];
ExMag=Table[{r1+(i-1)rinc,Abs[result[[i,2]]]},{i,inc}];
EyMag=Table[{r1+(i-1)rinc,Abs[result[[i,3]]]},{i,inc}];
EzMag=Table[{r1+(i-1)rinc,Abs[result[[i,4]]]},{i,inc}];
EsMag=Table[{r1+(i-1)rinc,Sqrt[Abs[result[[i,2]]]^2+Abs[result[[i,3]]]^2+Abs[result[[i,4]]]^2]},{i,inc}];
EtMag=Table[{r1+(i-1)rinc,Sqrt[Abs[result[[i,1]]+result[[i,2]]]^2+Abs[result[[i,3]]]^2+Abs[result[[i,4]]]^2]},{i,inc}];
Return[{ExMag,EyMag,EzMag,EsMag,EtMag}],

If[y1≠y2,        (*  vary position on y axis  *)
λ=λ1;θ=90;r1=y1;            
inc=Round[Abs[(y2-y1)/yinc]+1];
rinc=yinc;
If[y2>y1,rinc=Abs[rinc],rinc=-Abs[rinc]];
result=Table[r=r1+(i-1)rinc;If[r≥0,ϕ=90,ϕ=270;r=Abs[r]];
If[r==0,r=0.0001];
Compute[λ,a,nmed,nsph,r,θ,ϕ],{i,inc}];
ExMag=Table[{r1+(i-1)rinc,Abs[result[[i,2]]]},{i,inc}];
EyMag=Table[{r1+(i-1)rinc,Abs[result[[i,3]]]},{i,inc}];
EzMag=Table[{r1+(i-1)rinc,Abs[result[[i,4]]]},{i,inc}];
EsMag=Table[{r1+(i-1)rinc,Sqrt[Abs[result[[i,2]]]^2+Abs[result[[i,3]]]^2+Abs[result[[i,4]]]^2]},{i,inc}];
EtMag=Table[{r1+(i-1)rinc,Sqrt[Abs[result[[i,1]]+result[[i,2]]]^2+Abs[result[[i,3]]]^2+Abs[result[[i,4]]]^2]},{i,inc}];
Return[{ExMag,EyMag,EzMag,EsMag,EtMag}],

λ=λ1;ϕ=0;r1=z1;            (*  vary position on z axis  *) 
inc=Round[Abs[(z2-z1)/zinc]+1];
rinc=zinc;
If[z2>z1,rinc=Abs[rinc],rinc=-Abs[rinc]];
result=Table[r=r1+(i-1)rinc;If[r≥0,θ=0,θ=180;r=Abs[r]];
If[r==0,r=0.0001];
Compute[λ,a,nmed,nsph,r,θ,ϕ],{i,inc}];
ExMag=Table[{r1+(i-1)rinc,Abs[result[[i,2]]]},{i,inc}];
EyMag=Table[{r1+(i-1)rinc,Abs[result[[i,3]]]},{i,inc}];
EzMag=Table[{r1+(i-1)rinc,Abs[result[[i,4]]]},{i,inc}];
EsMag=Table[{r1+(i-1)rinc,Sqrt[Abs[result[[i,2]]]^2+Abs[result[[i,3]]]^2+Abs[result[[i,4]]]^2]},{i,inc}];
EtMag=Table[{r1+(i-1)rinc,Sqrt[Abs[result[[i,1]]+result[[i,2]]]^2+Abs[result[[i,3]]]^2+Abs[result[[i,4]]]^2]},{i,inc}];
Return[{ExMag,EyMag,EzMag,EsMag,EtMag}]]]]];

Mie scattering from a spherical shell on a spherical core

This program computes the extinction and scattering coefficients for a sphere with a shell.  It has been adapted from the Fortran code in C. F. Bohren and D. R. Huffman, Absorption and Scattering of Light by Small Particles (Wiley, New York, 1983) Appendix B.  The ShellScattering function computes the extinction coefficient, the scattering coefficient, and the backscattering coefficient for nanoshells.   The function returns three lists of (x, y) coordinates.  The x coordinate is either wavelength, shell thickness, or inner sphere radius for all three lists in units of microns.  The y coordinate for the first list is the scattering coefficient.  The y coordinate for the second list is the extinction coefficient, and the y coordinate for the third list is the backscattering coefficient.  
	For polar scattering diagrams, the function also returns four lists.  Let (Subscript[i, ⊥]) be the scattered irradiance per unit incident irradiance when the incident light is polarized perpendicular to the scattering plane, and let (Subscript[i, ∥]) be the scattered irradiance per unit incident irradiance when the incident light is polarized parallel to the scattering plane.  The first list is the set of coordinates (θ, (Subscript[i, ⊥])) as θ goes from 0° to 180°.  The second list is the same set of data expressed in a Cartesian frame, ((Subscript[i, ⊥]) cos(θ), (Subscript[i, ⊥])sin(θ)).  The third list is the set of coordinates (θ, (Subscript[i, ∥])) as θ goes from 0° to 180°, and the fourth list is the same set of data expressed in a Cartesian frame, ((Subscript[i, ∥]) cos(θ), (Subscript[i, ∥])sin(θ)).
	The parameters that are set by the user before a calculation are as follows:

		(Subscript[λ, 1]):  initial wavelength	
		(Subscript[λ, 2]):  final wavelength
		(Subscript[λ, inc]): increment of wavelength
		(Subscript[t, 1]): initial thickness of shell
		(Subscript[t, 2]): final thicknes of shell
		(Subscript[t, inc]): increment of shell thickness
		(Subscript[r, 1]): initial radius of sphere in microns
		(Subscript[r, 2]): final radius of sphere in microns
		(Subscript[r, inc]): increment of radius in microns
		(Subscript[n, 1]): refractive index of outer medium
		(Subscript[n, 2]): refractive index of coating
		(Subscript[n, 3]): refractive index of inner sphere
		(Subscript[θ, inc]): increment of angle for polar scattering diagram

	If (Subscript[λ, 1])≠(Subscript[λ, 2]), the function computes the coefficients vs. wavelength for shell thickness (Subscript[t, 1]) and shell radius (Subscript[r, 1]).  In this case, the function also returns three additional lists which are identical to the first three lists except that the x coordinate is energy in eV rather than wavelength.  If (Subscript[λ, 1])=(Subscript[λ, 2]), but (Subscript[t, 1])≠(Subscript[t, 2]), the function computes the coefficients vs. shell thickness at wavelength (Subscript[λ, 1]) and at radius (Subscript[r, 1]).  If (Subscript[λ, 1])=(Subscript[λ, 2]) and (Subscript[t, 1])=(Subscript[t, 2]), but (Subscript[r, 1])≠(Subscript[r, 2]) the function computes the coefficients vs. inner sphere radius at wavelength (Subscript[λ, 1]) and shell thickness (Subscript[t, 1]).  Finally, if (Subscript[λ, 1])=(Subscript[λ, 2]), (Subscript[t, 1])=(Subscript[t, 2]), and (Subscript[r, 1])=(Subscript[r, 2]) the function returns the four lists of polar scattering data as θ goes from 0° to 180° in increments of (Subscript[θ, inc]).
	As with the SphereScattering function, the user can either set the refractive index of the shell or sphere to a specific (complex) value, or the user can set (Subscript[n, 2]) or (Subscript[n, 3]) to a negative integer which tells the function to insert the wavelength-dependent refractive index for a particular metal as follows: 

	-1: complex index for gold interpolated from values in Palik 
	-2: complex index for silver interpolated from values in Palik  
	-3: complex index for copper interpolated from values in Palik 
	-4: complex index for aluminum interpolated from values in Palik
	-5: complex index of a Drude metal
	-6: complex index for gold from Drude-Lorentz fit from Rakik et al.
	-7: complex index for silver from Drude-Lorentz fit from Rakik et al.
	-8: complex index for copper from Drude-Lorentz fit from Rakik et al.
	-9: complex index for aluminum from Drude-Lorentz fit from Rakik et al.

f1[y_, n_, n1_, n2_, D1x1_, D1x2_, D1y2_, brack_, crack_, χy2_, χpy2_, ψy_, ψ1y_, ξy_, ξ1y_, χy_, χ1y_, χx2_, χy2_, χ1x2_, χ1y2_, Qsct_, Xbck_, Qext_,nang1_,π0_,π1_,amu_,s1_,s2_] := Module[{Qs , Xb, Qe, ψ0y, ψ1yy, χ0y, ξ1yy, χ0x2, χ0y2, χ1x22 , χ1yy, χ1y22, D0x1, D0x2, D0y2,dnbar,gnbar,an,bn,jj,pi,τ,p,fn,t,ss1,ss2},
ss1=s1;
ss2=s2;
   dnbar = (D1y2 - brack χpy2)/(1 - brack χy2); gnbar = (D1y2 - crack χpy2)/(1 - crack χy2);
   an = (((n1/n2) dnbar + n/y) ψy - ψ1y)/(((n1/n2) dnbar + n/y) ξy - ξ1y);
   bn = (((n2/n1) gnbar + n/y) ψy - ψ1y)/(((n2/n1) gnbar + n/y) ξy - ξ1y);
If[nang1≠1,                             (*   polar scattering diagram   *)
fn=(2 n+1)/(n(n+1));
p=(-1)^(n-1);
t=(-1)^n;
Do[jj=Round[2 nang1-j];
pi[j]=π1[j];
τ[j]=n amu[[j]] pi[j] - (n+1)π0[j];
ss1[j]=ss1[j] + fn (an pi[j] + bn τ[j]);
ss2[j]=ss2[j]+fn(an τ[j]+ bn pi[j]);
If[j!=jj,ss1[jj]=ss1[jj] + fn (an pi[j] p + bn τ[j] t);
ss2[jj]=ss2[jj] + fn (an τ[j] t + bn pi[j] p)],{j,nang1}]];
   Qs = Qsct + (2 n + 1) ((Abs[an])^2 + (Abs[bn])^2);
   Xb = Xbck + (2 n + 1) (-1)^n  (an - bn);
   Qe = Qext + (2 n + 1) (Re[an] + Re[bn]);
   ψ0y = ψ1y; ψ1yy = ψy;
   χ0y = χ1y; χ1yy = χy;
   ξ1yy = ψy - I χy;
   χ0x2 = χ1x2; χ1x22 = χx2;
   χ0y2 = χ1y2; χ1y22 = χy2;
   D0x1 = D1x1;
   D0x2 = D1x2;
   D0y2 = D1y2;
   Return[{Qs , Xb, Qe, ψ0y, ψ1yy, χ0y, χ1yy, ξ1yy, χ0x2, χ1x22, χ0y2 , χ1y22, D0x1, D0x2, D0y2,ss1,ss2,pi}]];

f2[λ_,r1_,r2_,n1_,n2_,n3_,nang1_] := Module[{Qs , Xb, Qe, ψ0y, ψ1yy, χ0y, χ1yy, ξ1yy, χ0x2, χ1x22, χ0y2 , χ1y22, D0x1, D0x2, D0y2,ncore,nshell,x,y,x1,x2,y2,Nstop,ψ1y,χ1y,ξ0y,ξ1y,χ1x2,χ1y2,Qsct,Qext,Qbck,Xbck,iflag,ψy,χy,ξy,D1x1,D1x2,D1y2,brack,crack,χy2,χpy2,χx2,χpx2,s1,s2,π0,π1,amu,pi,ss1,ss2},
ncore=Which[n3==-1,Gold[λ],n3==-2,Silver[λ],n3==-3,Copper[λ],n3==-4,Aluminum[λ],n3==-5,Drude[λ],n3==-6,LorentzDrudeAu[λ],n3==-7,LorentzDrudeAg[λ],n3==-8,LorentzDrudeCu[λ],n3==-9,LorentzDrudeAl[λ],n3==-10,Gold2[λ],n3==-11,Silver[λ],True,n3];
nshell=Which[n2==-1,Gold[λ],n2==-2,Silver[λ],n2==-3,Copper[λ],n2==-4,Aluminum[λ],n2==-5,Drude[λ],n2==-6,LorentzDrudeAu[λ],n2==-7,LorentzDrudeAg[λ],n2==-8,LorentzDrudeCu[λ],n2==-9,LorentzDrudeAl[λ],n2==-10,Gold2[λ],n2==-11,Silver[λ],True,n2];
x = 2 π r2 n1/λ;
y = 2 π r1 n1/λ;
x1= x(ncore/n1);
x2= x(nshell/n1);
y2= y(nshell/n1);
Nstop=Round[y+4 (y^(1/3))+2];    (* # terms in expansion *)
D0x1=Cot[x1];
D0x2=Cot[x2];
D0y2=Cot[y2];
ψ0y=Cos[y];ψ1y=Sin[y];
χ0y=-Sin[y];χ1y=Cos[y];
ξ0y=ψ0y-I χ0y;ξ1y=ψ1y-I χ1y;
χ0x2=-Sin[x2];χ1x2=Cos[x2];
χ0y2=-Sin[y2];χ1y2=Cos[y2];
Qsct=0;
Qext=0;
Xbck=0;
iflag=0;
If[nang1≠1,
Do[π0[j]=0;π1[j]=1,{j,nang1}];
Do[s1[j]=0;s2[j]=0,{j,2 nang1-1}];
amu=Table[Cos[π(j-1)/(2(nang1-1))],{j,nang1}],
π0=0;π1=0;s1=0;s2=0;amu=0];
Do[ψy=(2n-1)ψ1y/y - ψ0y;
χy=(2 n-1) χ1y/y - χ0y;
ξy=ψy - I χy;
D1y2=1/(n/y2 - D0y2) - n/y2;
If[iflag==1,{Qsct ,Xbck,Qext,ψ0y,ψ1y,χ0y,χ1y,ξ1y,χ0x2,χ1x2,χ0y2 ,χ1y2,D0x1,D0x2,D0y2,ss1,ss2,pi}=f1[y,n,n1,nshell,D1x1,D1x2,D1y2,brack,crack,χy2,χpy2,ψy,ψ1y,ξy,ξ1y,χy,χ1y,χx2,χy2,χ1x2,χ1y2,Qsct,Xbck,Qext,nang1,π0,π1,amu,s1,s2];
If[nang1≠1,
s1=ss1;
s2=ss2;
Do[π1[j]=((2 n+1)/n) amu[[j]] pi[j]-(n+1) π0[j]/n;
π0[j]=pi[j],{j,nang1}]],
D1x1=1/(n/x1 - D0x1) - n/x1;
D1x2=1/(n/x2 - D0x2) - n/x2;
χx2 =(2n - 1)χ1x2/x2 - χ0x2;
χy2 =(2 n - 1)χ1y2/y2 - χ0y2;
χpx2=χ1x2 - n χx2/x2;
χpy2=χ1y2 - n χy2/y2;
brack=((nshell/ncore) D1x1 - D1x2)(χy2 D1y2 - χpy2)/(((nshell/ncore)D1x1 χx2 - χpx2)(χx2 D1x2 - χpx2));
crack=((nshell/ncore) D1x2 - D1x1)(χy2 D1y2 - χpy2)/(((nshell/ncore)χpx2 - D1x1 χx2)(χx2 D1x2 - χpx2));
If[(Abs[brack χpy2]<10^-8 Abs[D1y2])&&(Abs[brack χy2]<10^-8)&&(Abs[crack χpy2]<10^-8 Abs[D1y2])&&(Abs[crack χy2]<10^-8),brack = 0; crack = 0; iflag = 1];
{Qsct ,Xbck,Qext,ψ0y,ψ1y,χ0y,χ1y,ξ1y,χ0x2,χ1x2,χ0y2 ,χ1y2,D0x1,D0x2,D0y2,ss1,ss2,pi}=f1[y,n,n1,nshell,D1x1,D1x2,D1y2,brack,crack,χy2,χpy2,ψy,ψ1y,ξy,ξ1y,χy,χ1y,χx2,χy2,χ1x2,χ1y2,Qsct,Xbck,Qext,nang1,π0,π1,amu,s1,s2]];
If[nang1≠1,
s1=ss1;
s2=ss2;
Do[π1[j]=((2 n+1)/n) amu[[j]] pi[j]-(n+1) π0[j]/n;
π0[j]=pi[j],{j,nang1}]],{n,Nstop}];
Qsct=2*Qsct/((y^2));
Qext=2*Qext/((y^2));
Qbck=(Abs[Xbck])^2/((y^2));
If[nang1==1,Return[{Qsct,Qext,Qbck}],Return[{s1,s2}]]];

ShellScattering[λ1_,λ2_,λinc_, t1_,t2_,tinc_,r1_, r2_,rinc_, n1_, n2_, n3_,θinc_]:=Module[{Qsct,Qbck,Qext,tinc1,nt,rr2,nr,rinc1,rr1,nλ,λ,result,Qs,Qe,Qb,QsE,QeE,QbE,s1,s2,nang,λinc1,s11norm,s11,s12,perp,perp1,par,par1,pol,s33,s34,p1,pp1,p2,pp2,θ},
If[λ2≠λ1,            (*  find wavelength dependence  *)
If[λ2>λ1,λinc1=Abs[λinc],λinc1=-Abs[λinc] ];
nλ=Round[Abs[(λ2-λ1)/λinc1]+1];
nang=1;
result=Table[λ=λ1+(j-1)λinc1;
f2[λ,r1+t1,r1,n1,n2,n3,nang],{j,nλ}];
Qs=Table[{λ1+(j-1)λinc1,result[[j,1]]},{j,nλ}];
Qe=Table[{λ1+(j-1)λinc1,result[[j,2]]},{j,nλ}];
Qb=Table[{λ1+(j-1)λinc1,result[[j,3]]},{j,nλ}];
QsE=Table[{1.24/(λ1+(j-1)λinc1),result[[j,1]]},{j,nλ}];
QeE=Table[{1.24/(λ1+(j-1)λinc1),result[[j,2]]},{j,nλ}];
QbE=Table[{1.24/(λ1+(j-1)λinc1),result[[j,3]]},{j,nλ}];
Return[{Qs,Qe,Qb,QsE,QeE,QbE}];
Print["Made it here"],

If[t1≠t2,            (*  find shell thickness dependence  *)
If[t2>t1,tinc1=Abs[tinc],tinc1=-Abs[tinc] ];
nt = Round[Abs[(t2-t1)/tinc1]+1];
nang=1;
result=Table[rr2=r1+(j-1)tinc1;
f2[λ1,rr2,r1,n1,n2,n3,nang],{j,nr}];
Qs=Table[{(j-1)tinc1,result[[j,1]]},{j,nt}];
Qe=Table[{(j-1)tinc1,result[[j,2]]},{j,nt}];
Qb=Table[{(j-1)tinc1,result[[j,3]]},{j,nt}];
Return[{Qs,Qe,Qb}],

If[r1≠r2,            (*  find core radius dependence  *)
If[r2>r1,rinc1=Abs[rinc],rinc1=-Abs[rinc] ];
nr = Round[Abs[(r2-r1)/rinc1]+1];
nang=1;
result=Table[rr1=r1+(j-1)rinc1;
f2[λ1,rr1+t1,rr1,n1,n2,n3,nang],{j,nr}];
Qs=Table[{r1+(j-1)rinc1,result[[j,1]]},{j,nr}];
Qe=Table[{r1+(j-1)rinc1,result[[j,2]]},{j,nr}];
Qb=Table[{r1+(j-1)rinc1,result[[j,3]]},{j,nr}];
Return[{Qs,Qe,Qb}],

nang = Round[Abs[90/θinc]+1];          (*  find angle dependence  *)
{s1,s2}=f2[λ1,r1+t1,r1,n1,n2,n3,nang];
s11norm=(Abs[s1[1]]^2+ Abs[s2[1]]^2 )/2;
Do[θ=θinc(j-1);
s11=(Abs[s1[j]]^2+Abs[s2[j]]^2)/2;
s12=(Abs[s2[j]]^2-Abs[s1[j]]^2)/2;
perp[j]={θ,s11-s12};
perp1[j]={(s11-s12)Cos[π θ/180.0],(s11-s12)Sin[π θ/180.0]};
par[j]={θ,s11+s12};
par1[j]={(s11+s12)Cos[π θ/180.0],(s11+s12)Sin[π θ/180.0]};
pol=-s12/s11;
s33=Re[s2[j] Conjugate[s1[j]]]/s11;
s34=Im[s2[j] Conjugate[s1[j]]]/s11;
s11=s11/s11norm,{j,2 nang-1}] ;
p1=Table[perp[j],{j,2 nang-1}];   (* perp scattering in polar coords *)
pp1=Table[perp1[j],{j,2 nang-1}]; (* perp scattering in x-y coords  *)
p2=Table[par[j],{j,2 nang-1}];      (* par scattering in polar coords *)
pp2=Table[par1[j],{j,2 nang-1}];  (* par scattering in x-y coords  *)
Return[{p1,pp1,p2,pp2}]]]]];

Examples from the literature

In Fig. 4 of Messinger et al., Phys. Rev. B 24 (1981) 649 the far field and near field scattering coefficients are computed for a 44 nm silver sphere in water.

λ1=0.3;                  (*  initial wavelength  *)
λ2=0.7;                  (*  final wavelength  *)
λinc=0.001;                  (*  increment of wavelength  *)
r1=0.022;                  (*  initial radius of sphere in microns  *)
r2=0.022;                  (*  final radius of sphere in microns  *)
rinc = 0.0;              (*  increment of radius  *)
θinc = 1;                   (*  increment of angle  *)     
n2=-11;                   (*  refractive index of silver sphere from J & C *)
n1 = 1.33;                   (* refractive index of medium  *)
{QsAg,QeAg,QbAg,QnAg}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];

ListLinePlot[QnAg,PlotLegend->{"Qn"},LegendPosition-> {.45,.2},LegendSize->0.25,LegendShadow->{.01,-.01},FrameLabel->{"λ (μm)","(Subscript[Q, NF])"},Evaluate[option]]

ListLinePlot[{QsAg, QeAg}, PlotLegend -> {"Qs", "Qe"}, LegendPosition -> {.4, .2}, LegendSize -> 0.25, LegendShadow -> {.01, -.01}, FrameLabel -> {"λ (μm)", "Coefficient for Ag Sphere"}, Evaluate[option]]

In Fig. 5 of Messinger et al., Phys. Rev. B 24 (1981) 649 the far field and near field scattering coefficients are computed for a 200 nm silver sphere in water.

λ1=0.3;                  (*  initial wavelength  *)
λ2=0.7;                  (*  final wavelength  *)
λinc=0.001;                  (*  increment of wavelength  *)
r1=0.1;                  (*  initial radius of sphere in microns  *)
r2=0.1;                  (*  final radius of sphere in microns  *)
rinc = 0.0;              (*  increment of radius  *)
θinc = 1;                   (*  increment of angle  *)     
n2=-11;                   (*  refractive index of silver sphere from J & C *)
n1 = 1.33;                   (* refractive index of medium  *)
{QsAg,QeAg,QbAg,QnAg}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];

ListLinePlot[QnAg,PlotLegend->{"Qn"},LegendPosition-> {.4,.2},LegendSize->0.25,LegendShadow->{.01,-.01},FrameLabel->{"λ (μm)","Near Field Coefficient"},Evaluate[option]]

ListLinePlot[{QsAg, QeAg}, PlotLegend -> {"Qs", "Qe"}, LegendPosition -> {.4, .2}, LegendSize -> 0.25, LegendShadow -> {.01, -.01}, FrameLabel -> {"λ (μm)", "Coefficients"}, Evaluate[option]]

In Fig. 4 of Averitt et al., J. Opt. Soc. Am. B, 16 (1999) 1824-1832 the quasistatic approximation for the extinction cross section is compared to that calculated from Mie theory for nanoshells that are very small compared to the incident wavelength and suspended in water.  The agreement is excellent as it should be.  We show the results of the Mie theory calculation below for that same figure.  The gold shell thickness is fixed at 2 nm while the inner core thickness of (Subscript[Au, 2])S is varied between 2 nm, 8 nm, and 15 nm.  The refractive index vallues for gold are taken from Johnson and Christy, Phys. Rev. B, vol. 6 (1972) 4370-4379.

λ1 = 0.4;            (* initial wavelength in microns *)
λ2 = 1.1;          (* final wavelength in microns  *)
λinc = 0.01;      (* increment of wavelength  *)
t1 = 0.002;       (*  initial thickness of shell  *)
t2 = 0.002;        (*  final thickness of shell  *)
tinc = 0    ;     (*  increment of shell thickness or sphere radius *)
r2 = 0.004;        (*  final radius of inner sphere    *)
rinc = 0.0;        (*  increment of radius  *)
n1 = 1.334;          (*  index of surrounding medium  *)
n2 = -10;              (*  index of gold from Johnson and Christy *)
n3 = 2.3328;           (*  index of inner core  *)
θinc=1;

r1 = 0.002;        (*  initial radius of inner sphere  *)
{Qs1, Qe1, Qb1,QsE, QeE, QbE} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
r1 = 0.008;        (*  initial radius of inner sphere  *)
{Qs2, Qe2, Qb2,QsE, QeE, QbE} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
r1 = 0.015;        (*  initial radius of inner sphere  *)
{Qs3, Qe3, Qb3,QsE, QeE, QbE} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc];  
ListLinePlot[{Qe1,Qe2,Qe3}, LegendPosition -> {.8, -.4}, PlotLegend->{"4 nm","10 nm","17 nm"},LegendSize -> 0.5, LegendShadow -> {.01, -.01}, LegendLabel -> "Total", FrameLabel -> {"λ (nm)", "Extinction Coefficient "}, Evaluate[option]]

In Fig. 2 of Oldenburg et al., Appl. Phys. Lett 75 (1999) 2897-2899 experimental spectra are compared to Mie theory for large nanoshells which exhibit multiple multipole resonances and for which the dipole resonance is shifted well into the IR.  We show the results of the Mie theory calculation below for that same figure.  The gold shell thickness and silica core radius for the four plots are (13 nm, 55 nm), (13 nm, 75 nm), (15 nm, 135 nm), and (13 nm, 210 nm), respectively.   The inner core is silica and the nanoshells are immersed in water.  The refractive index values for gold are taken from the Lorentz-Drude fit by Rakik et al., Appl. Opt. 37 (1998) 5271 but are not modified to include the effects of excess electron scattering from the surfaces of the shell.

λ1 = 0.4;          (* initial wavelength in microns *)
λ2 = 2.5;          (* final wavelength in microns  *)
λinc = 0.01;      (* increment of wavelength  *)
tinc = 0    ;     (*  increment of shell thickness or sphere radius *)
rinc = 0.0;        (*  increment of radius  *)
n1 = 1.334;          (*  index of surrounding medium  *)
n2 = -6;              (*  index of gold from Rakik et al. *)
n3 = 1.5;           (*  index of inner core  *)
θinc=1;

r1 = 0.055;        (*  initial radius of inner sphere  *)
r2 = 0.055;        (*  final radius of inner sphere    *)
t1 = 0.013;       (*  initial thickness of shell  *)
t2 = 0.013;        (*  final thickness of shell  *)
{Qs1, Qe1, Qb1,QsE, QeE, QbE} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
r1 = 0.075;        (*  initial radius of inner sphere  *)
r2 = 0.075;        (*  final radius of inner sphere    *)
t1 = 0.013;       (*  initial thickness of shell  *)
t2 = 0.013;        (*  final thickness of shell  *)
{Qs2, Qe2, Qb2,QsE, QeE, QbE} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
r1 = 0.135;        (*  initial radius of inner sphere  *)
r2 = 0.135;        (*  final radius of inner sphere    *)
t1 = 0.015;       (*  initial thickness of shell  *)
t2 = 0.015;        (*  final thickness of shell  *)
{Qs3, Qe3, Qb3,QsE, QeE, QbE} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
r1 = 0.21;        (*  initial radius of inner sphere  *)
r2 = 0.21;        (*  final radius of inner sphere    *)
t1 = 0.013;       (*  initial thickness of shell  *)
t2 = 0.013;        (*  final thickness of shell  *)
{Qs4, Qe4, Qb4,QsE, QeE, QbE} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
ListLinePlot[{Qe1,Qe2,Qe3,Qe4}, LegendPosition -> {.8, -.4}, PlotLegend->{"55 nm","75 nm","135 nm","210 nm"},LegendSize -> 0.5, LegendShadow -> {.01, -.01}, LegendLabel -> "Core", FrameLabel -> {"λ (nm)", "Extinction Coefficient "}, Evaluate[option]]

In Fig. 1 of Jackson et al., J. Phys. Chem. B 105 (2001) 2743-2746 the extinction spectra for silver and gold nanoshells of the same dimensions are compared.  It is noted that the peak extinction coefficient for silver is only slightly larger than that of gold.

λ1 = 0.4;          (* initial wavelength in microns *)
λ2 = 1.5;          (* final wavelength in microns  *)
λinc = 0.01;      (* increment of wavelength  *)
r1 = 0.05;        (*  initial radius of inner sphere  *)
r2 = 0.05;        (*  final radius of inner sphere    *)
rinc = 0.0;        (*  increment of radius  *)
t1 = 0.01;       (*  initial thickness of shell  *)
t2 = 0.01;        (*  final thickness of shell  *)
tinc = 0    ;     (*  increment of shell thickness or sphere radius *)
n1 = 1.334;          (*  index of surrounding medium - water  *)
n3 = 1.5;           (*  index of inner core - silica  *)
θinc=1;

n2 = -6;              (*  index of gold from Rakik et al. *)
{Qs1, Qe1, Qb1,QsE, QeE, QbE} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
n2 = -7;              (*  index of silver from Rakik et al. *)
{Qs2, Qe2, Qb2,QsE, QeE, QbE} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
ListLinePlot[{Qe1,Qe2}, LegendPosition -> {.3, .2}, PlotLegend->{"Au","Ag"},LegendSize -> 0.3, LegendShadow -> {.01, -.01}, FrameLabel -> {"λ (μm)", "Extinction Coefficient "}, Evaluate[option]]

In Fig. 5 of Nehl et al., Nano Lett. 4 (2004) 2355-2359 the extinction spectra for gold nanoshells on silica and immersed in immersion oil.

λ1 = 0.413;          (* initial wavelength in microns *)
λ2 = 1.24;          (* final wavelength in microns  *)
λinc = 0.01;      (* increment of wavelength  *)
r1 = 0.052;        (*  initial radius of inner sphere  *)
r2 = 0.052;        (*  final radius of inner sphere    *)
rinc = 0.0;        (*  increment of radius  *)
t1 = 0.02;       (*  initial thickness of shell  *)
t2 = 0.02;        (*  final thickness of shell  *)
tinc = 0    ;     (*  increment of shell thickness or sphere radius *)
n1 = 1.5;          (*  index of surrounding medium - immersion oil  *)
n2 = -6;              (*  index of gold from Rakik et al. *)
n3 = 1.5;           (*  index of inner core - silica  *)
θinc=1;

{Qs, Qe, Qb,QsE, QeE, QbE} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
ListLinePlot[QsE, FrameLabel -> {"Energy (eV)", "Scattering Coefficient "}, Evaluate[option]]

Code for generating figures

9.1 Nanoparticles

9.1.2 Mie theory examples

λ1=0.30;                  (*  initial wavelength  *)
λ2=0.7;                  (*  final wavelength  *)
λinc=0.001;                  (*  increment of wavelength  *)
θinc = 1;                   (*  increment of angle  *)   
r1=0.03;                  (*  initial radius of sphere in microns  *)
r2=0.03;                  (*  final radius of sphere in microns  *)
rinc = 0.0;              (*  increment of radius  *)
n1 = 1.0;                   (* refractive index of medium  *)  
n2=-1;                   (*  refractive index of gold sphere  *)
{QsAu,QeAu,QbAu,QnAu}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
n2=-2;                   (*  refractive index of silver sphere  *)
{QsAg,QeAg,QbAg,QnAg}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
n2=-3;                   (*  refractive index of copper sphere  *)
{QsCu,QeCu,QbCu,QnCu}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
n2=-4;                   (*  refractive index of aluminum sphere  *)
{QsAl,QeAl,QbAl,QnAl}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];

Fig2=ListLinePlot[{QeAu,QeAg,QeCu,QeAl},PlotLegend->{"Au","Ag","Cu","Al"},
LegendPosition-> {0.4,0.1},LegendSize->0.35,
LegendShadow->{.01,-.01},FrameLabel->{"λ (μm)","(Subscript[Q, ext])"},
Evaluate[option]]
Fig2BW=ListLinePlot[{QeAu,QeAg,QeCu,QeAl},PlotLegend->{"Au","Ag","Cu","Al"},
LegendPosition-> {0.4,0.1},LegendSize->0.35,
LegendShadow->{.01,-.01},FrameLabel->{"λ (μm)","(Subscript[Q, ext])"},
PlotRange->{All,{0,9}},
Evaluate[optB]];

Export[FileNameJoin[{base,"Fig9.2.JPG"}],Fig2BW,ImageResolution->600];
Clear[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];

λ1= 0.38;          (*  initial wavelength in microns *)
λ2= 0.38;          (*  final wavelength in microns *)
λinc= 0;          (*  increment of wavelength in microns *)
a=0.03;            (*  radius of sphere in microns  *)
nmed = 1;          (*  index of surrounding medium  *)
nsph = -2;              (*  index of silver  *)
x1 = 0.001;             (*  initial x position  *)
x2 = 0.1;        (*  final x position  *)
xinc = 0.001;        (*  increment of x  *)
y1 = 0;             (*  initial y position  *)
y2 = 0;            (*  final y position  *)
yinc = 0;        (*  increment of y  *)
z1 = 0;             (*  initial z position  *)
z2 = 0;            (*  final z position  *)
zinc = 0;        (*  increment of z  *)

{ExMag,EyMag,EzMag,EsMag,EtMag} = MieNearField[λ1,λ2,λinc,a,nmed,nsph,x1,x2,xinc,y1,y2,yinc,z1,z2,zinc]; 
Fig4=ListLinePlot[{ExMag,EzMag,EsMag},
 LegendPosition -> {-0.1, 0.2},PlotLegend->{"|(Subscript[E, x])| (radial)","|(Subscript[E, z])| (tang)","|E|"},
LegendSize -> {0.7,0.3}, LegendShadow -> {.01, -.01},  FrameLabel -> {"x (μm)", "Field amplitudes "},
Evaluate[option]]
Fig4BW=ListLinePlot[{ExMag,EzMag,EsMag},
 LegendPosition -> {-0.1, 0.2},PlotLegend->{"|(Subscript[E, x])| (radial)","|(Subscript[E, z])| (tang)","|E|"},
LegendSize -> {0.7,0.3}, LegendShadow -> {.01, -.01},  FrameLabel -> {"x (μm)", "Field amplitudes "},
Evaluate[optB]];

Export[FileNameJoin[{base,"Fig9.4.JPG"}],Fig4BW,ImageResolution->600];
Clear[λ1,λ2,λinc,a,nmed,nsph,x1,x2,xinc,y1,y2,yinc,z1,z2,zinc];

λ1= 0.3;          (*  initial wavelength in microns *)
λ2= 0.7;          (*  final wavelength in microns *)
λinc= 0.002;          (*  increment of wavelength in microns *)
a=0.03;            (*  radius of sphere in microns  *)
nmed = 1;          (*  index of surrounding medium  *)
x1 = 0.03;             (*  initial x position  *)
x2 = 0.03;        (*  final x position  *)
xinc = 0;        (*  increment of x  *)
y1 = 0;             (*  initial y position  *)
y2 = 0;            (*  final y position  *)
yinc = 0;        (*  increment of y  *)
z1 = 0;             (*  initial z position  *)
z2 = 0;            (*  final z position  *)
zinc = 0;        (*  increment of z  *)

nsph =-1;              (*  index of gold  *)
{ExMagAu,EyMagAu,EzMagAu,EsMagAu,EtMagAu} = MieNearField[λ1,λ2,λinc,a,nmed,nsph,x1,x2,xinc,y1,y2,yinc,z1,z2,zinc]; 
nsph =-2;             (*  index of silver  *)
{ExMagAg,EyMagAg,EzMagAg,EsMagAg,EtMagAg} = MieNearField[λ1,λ2,λinc,a,nmed,nsph,x1,x2,xinc,y1,y2,yinc,z1,z2,zinc]; 
Fig5=ListLinePlot[{EtMagAu,EtMagAg},
 LegendPosition -> {.4, .3},PlotLegend->{"Au","Ag"},
LegendSize -> 0.25,LegendShadow -> {.01, -.01},  FrameLabel -> {"λ (μm)", "Field amplitudes "},
Evaluate[option]]
Fig5BW=ListLinePlot[{EtMagAu,EtMagAg},
 LegendPosition -> {.4, .3},PlotLegend->{"Au","Ag"},
LegendSize -> 0.25,LegendShadow -> {.01, -.01},  FrameLabel -> {"λ (μm)", "Field amplitudes "},
Evaluate[optB]];

Export[FileNameJoin[{base,"Fig9.5.JPG"}],Fig5BW,ImageResolution->600];
Clear[λ1,λ2,λinc,a,nmed,nsph,x1,x2,xinc,y1,y2,yinc,z1,z2,zinc];

λ1=0.3;                  (*  initial wavelength  *)
λ2=0.7;                  (*  final wavelength  *)
λinc=0.002;                  (*  increment of wavelength  *)
r1=0.03;                  (*  initial radius of sphere in microns  *)
r2=0.03;                  (*  final radius of sphere in microns  *)
rinc = 0.0;              (*  increment of radius  *)
θinc = 1;                   (*  increment of angle  *)     
n1 = 1;                   (* refractive index of medium  *)

n2=-1;                   (*  refractive index of gold sphere *)
{QsAu,QeAu,QbAu,QnAu}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];n2=-2;                   (*  refractive index of silver sphere *)
{QsAg,QeAg,QbAg,QnAg}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];

Fig6=ListLinePlot[{QnAu,QnAg},PlotLegend->{"Au","Ag"},LegendPosition-> {.3,.2},LegendSize->{0.3,0.25},
LegendShadow->{.01,-.01},FrameLabel->{"λ (μm)","(Subscript[Q, NF])"},
Evaluate[option]]
Fig6BW=ListLinePlot[{QnAu,QnAg},PlotLegend->{"Au","Ag"},LegendPosition-> {.3,.2},LegendSize->{0.3,0.25},
LegendShadow->{.01,-.01},
FrameLabel->{"λ (μm)","(Subscript[Q, NF])"},
Evaluate[optB]];

Export[FileNameJoin[{base,"Fig9.6.JPG"}],Fig6BW,ImageResolution->600];

λ1=0.3;                  (*  initial wavelength  *)
λ2=0.9;                  (*  final wavelength  *)
λinc=0.001;                  (*  increment of wavelength  *)
rinc = 0.0;              (*  increment of radius  *)
n1 = 1.0;                   (* refractive index of medium  *)
θinc = 1;                   (*  increment of angle  *)     
n2=-2;                   (*  refractive index of silver sphere  *)
r1=0.01;                  (*  initial radius of sphere in microns  *)
r2=0.01;                  (*  final radius of sphere in microns  *)
{QsAg1,QeAg1,QbAg1,QnAg1}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
r1=0.03;                  (*  initial radius of sphere in microns  *)
r2=0.03;                  (*  final radius of sphere in microns  *)
{QsAg2,QeAg2,QbAg2,QnAg2}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
r1=0.06;                  (*  initial radius of sphere in microns  *)
r2=0.06;                  (*  final radius of sphere in microns  *)
{QsAg3,QeAg3,QbAg3,QnAg3}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
r1=0.12;                  (*  initial radius of sphere in microns  *)
r2=0.12;                  (*  final radius of sphere in microns  *)
{QsAg4,QeAg4,QbAg4,QnAg4}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
r1=0.15;                  (*  initial radius of sphere in microns  *)
r2=0.15;                  (*  final radius of sphere in microns  *)
{QsAg5,QeAg5,QbAg5,QnAg5}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];

Fig7=ListLinePlot[{QeAg1,QeAg2,QeAg3,QeAg4,QeAg5},PlotLegend->{"20 nm","60 nm","120 nm","240 nm","300 nm"},LegendPosition-> {0.2,.15},LegendSize->{0.45,0.45},
LegendShadow->{.01,-.01},FrameLabel->{"λ (μm)","Extinction coefficient"},
Evaluate[option]] 
Fig7BW=ListLinePlot[{QeAg1,QeAg2,QeAg3,QeAg4,QeAg5},PlotLegend->{"20 nm","60 nm","120 nm","240 nm","300 nm"},LegendPosition-> {0.2,.15},LegendSize->{0.45,0.45},
LegendShadow->{.01,-.01},FrameLabel->{"λ (μm)","Extinction coefficient"},
Evaluate[optB]];

Export[FileNameJoin[{base,"Fig9.7.JPG"}],Fig7BW,ImageResolution->600];
Clear[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];

Fig8=ListLinePlot[{QnAg1, QnAg2, QnAg3, QnAg4, QnAg5}, PlotLegend -> {"20 nm", "60 nm", "120 nm", "240 nm", "300 nm"}, LegendPosition -> {0.3, 0.1},LegendSize -> {0.45,0.45},
 LegendShadow -> {.01, -.01},FrameLabel -> {"λ (μm)", "Near field coefficient"},
Evaluate[option]] 
Fig8BW=ListLinePlot[{QnAg1, QnAg2, QnAg3, QnAg4, QnAg5}, PlotLegend -> {"20 nm", "60 nm", "120 nm", "240 nm", "300 nm"}, LegendPosition -> {0.3, 0.1},LegendSize -> {0.45,0.45},
 LegendShadow -> {.01, -.01},FrameLabel -> {"λ (μm)", "Near field coefficient"},
Evaluate[optB]];

Export[FileNameJoin[{base,"Fig9.8.JPG"}],Fig8BW,ImageResolution->600];
Clear[λ1,λ2,λinc,a,nmed,nsph,x1,x2,xinc,y1,y2,yinc,z1,z2,zinc];

λ1=0.3;                  (*  initial wavelength  *)
λ2=0.9;                  (*  final wavelength  *)
λinc=0.001;                  (*  increment of wavelength  *)
r1=0.02;                  (*  initial radius of sphere in microns  *)
r2=0.02;                  (*  final radius of sphere in microns  *)
rinc = 0.0;              (*  increment of radius  *)
θinc = 1;                   (*  increment of angle  *)     
n2=-3;                   (*  refractive index of copper sphere  *)
n1 = 1.0;                   (* refractive index of medium  *)
{QsCu1,QeCu1,QbCu1,QnCu1}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
n1 = 1.3;                   (* refractive index of medium  *)
{QsCu2,QeCu2,QbCu2,QnCu2}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
n1 = 1.5;                   (* refractive index of medium  *)
{QsCu3,QeCu3,QbCu3,QnCu3}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
n1 = 1.7;                   (* refractive index of medium  *)
{QsCu4,QeCu4,QbCu4,QnCu4}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
n1 = 1.9;                   (* refractive index of medium  *)
{QsCu5,QeCu5,QbCu5,QnCu5}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
n1 = 2.1;                   (* refractive index of medium  *)
{QsCu6,QeCu6,QbCu6,QnCu6}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];

Fig10=ListLinePlot[{QeCu1,QeCu2,QeCu3,QeCu4,QeCu5,QeCu6},PlotLegend->{"n=1.0","n=1.3","n=1.5","n=1.7","n=1.9","n=2.1"},LegendPosition-> {0.4,0},LegendSize->{0.4,0.5},
LegendShadow->{.01,-.01},FrameLabel->{"λ (μm)","(Subscript[Q, ext])"},
Evaluate[option]] 
Fig10BW=ListLinePlot[{QeCu1,QeCu2,QeCu3,QeCu4,QeCu5,QeCu6},PlotLegend->{"n=1.0","n=1.3","n=1.5","n=1.7","n=1.9","n=2.1"},LegendPosition-> {0.4,0},LegendSize->{0.4,0.5},
LegendShadow->{.01,-.01},FrameLabel->{"λ (μm)","(Subscript[Q, ext])"},
Evaluate[optB]] ;

Export[FileNameJoin[{base,"Fig9.10.JPG"}],Fig10BW,ImageResolution->600];

λ1=0.3;                  (*  initial wavelength  *)
λ2=0.9;                  (*  final wavelength  *)
λinc=0.001;                  (*  increment of wavelength  *)
r1=0.02;                  (*  initial radius of sphere in microns  *)
r2=0.02;                  (*  final radius of sphere in microns  *)
rinc = 0.0;              (*  increment of radius  *)
θinc = 1;                   (*  increment of angle  *)     
n2=-3;                   (*  refractive index of copper sphere  *)
n1 = 1.0;                   (* refractive index of medium  *)
{QsCu1,QeCu1,QbCu1,QnCu1}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
n1 = 1.3;                   (* refractive index of medium  *)
{QsCu2,QeCu2,QbCu2,QnCu2}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
n1 = 1.5;                   (* refractive index of medium  *)
{QsCu3,QeCu3,QbCu3,QnCu3}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
n1 = 1.7;                   (* refractive index of medium  *)
{QsCu4,QeCu4,QbCu4,QnCu4}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
n1 = 1.9;                   (* refractive index of medium  *)
{QsCu5,QeCu5,QbCu5,QnCu5}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
n1 = 2.1;                   (* refractive index of medium  *)
{QsCu6,QeCu6,QbCu6,QnCu6}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];

Fig11=ListLinePlot[{QnCu1,QnCu2,QnCu3,QnCu4,QnCu5,QnCu6},PlotLegend->{"n=1.0","n=1.3","n=1.5","n=1.7","n=1.9","n=2.1"},LegendPosition-> {-0.6,0},LegendSize->{0.3,0.5},
LegendShadow->{.01,-.01},FrameLabel->{"λ (μm)","(Subscript[Q, NF])"},
Evaluate[option]]
Fig11BW=ListLinePlot[{QnCu1,QnCu2,QnCu3,QnCu4,QnCu5,QnCu6},PlotLegend->{"n=1.0","n=1.3","n=1.5","n=1.7","n=1.9","n=2.1"},LegendPosition-> {-0.6,0},LegendSize->{0.3,0.5},
LegendShadow->{.01,-.01},FrameLabel->{"λ (μm)","(Subscript[Q, NF])"},
Evaluate[optB]];

Export[FileNameJoin[{base,"Fig9.11.JPG"}],Fig11BW,ImageResolution->600];
Clear[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];

λ1= 0.623;          (*  initial wavelength in microns *)
λ2= 0.623;          (*  final wavelength in microns *)
λinc= 0;          (*  increment of wavelength in microns *)
a=0.03;            (*  radius of sphere in microns  *)
nsph =-3;              (*  index of copper  *)
x1 = 0;             (*  initial x position  *)
x2 = 0.1;        (*  final x position  *)
xinc = 0.002;        (*  increment of x  *)
y1 = 0;             (*  initial y position  *)
y2 = 0;            (*  final y position  *)
yinc = 0;        (*  increment of y  *)
z1 = 0;             (*  initial z position  *)
z2 = 0;            (*  final z position  *)
zinc = 0;        (*  increment of z  *)

nmed = 1;          (*  index of surrounding medium  *)
{ExMag1,EyMag1,EzMag1,EsMag1,EtMag1} = MieNearField[λ1,λ2,λinc,a,nmed,nsph,x1,x2,xinc,y1,y2,yinc,z1,z2,zinc]; 
nmed = 1.5;          (*  index of surrounding medium  *)
{ExMag2,EyMag2,EzMag2,EsMag2,EtMag2} = MieNearField[λ1,λ2,λinc,a,nmed,nsph,x1,x2,xinc,y1,y2,yinc,z1,z2,zinc]; 
nmed = 2.1;          (*  index of surrounding medium  *)
{ExMag3,EyMag3,EzMag3,EsMag3,EtMag3} = MieNearField[λ1,λ2,λinc,a,nmed,nsph,x1,x2,xinc,y1,y2,yinc,z1,z2,zinc]; 
Fig12=ListLinePlot[{EtMag1,EtMag2,EtMag3},
 LegendPosition -> {.2, .2},
 PlotLegend->{"n=1.0","n=1.5","n=2.1"},
LegendSize -> {0.4,0.3}, LegendShadow -> {.01, -.01},  FrameLabel -> {"x (μm)", "Field amplitudes "},
Evaluate[option]] 
Fig12BW=ListLinePlot[{EtMag1,EtMag2,EtMag3},
 LegendPosition -> {.2, .2},
 PlotLegend->{"n=1.0","n=1.5","n=2.1"},
LegendSize -> {0.4,0.3}, LegendShadow -> {.01, -.01},  FrameLabel -> {"x (μm)", "Field amplitudes "},
Evaluate[optB]];

Export[FileNameJoin[{base,"Fig9.12.JPG"}],Fig12BW,ImageResolution->600];
Clear[λ1,λ2,λinc,a,nmed,nsph,x1,x2,xinc,y1,y2,yinc,z1,z2,zinc];

λ1=0.55;               (*  initial wavelength  *)
λ2=0.55;                  (*  final wavelength  *)
λinc=0;                  (*  increment of wavelength  *)
r1=0.02;                  (*  initial radius of sphere in microns  *)
r2=0.02;                  (*  final radius of sphere in microns  *)
rinc = 0;              (*  increment of radius  *)
n1 = 1.33;                   (* refractive index of medium  *)
n2=0.57+I 2.45;    (*  refractive index of gold sphere  *)
θinc = 1;                   (*  increment of angle  *)     

{p1,pp1,p2,pp2,p3,pp3}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
Fig13a=ListLinePlot[{pp1,pp2,pp3},PlotLegend->{" ⊥"," ∥","⊥+∥"},LegendPosition-> {.44,.16},LegendPosition-> {1-1.01,-0.2+.35},LegendSize->0.3,LegendShadow->{.01,-.01},
FrameLabel->{"\!\(\*
StyleBox[\"x\",\nFontSlant->\"Italic\"]\)","\!\(\*
StyleBox[\"y\",\nFontSlant->\"Italic\"]\)","\!\(\*
StyleBox[\"d\",\nFontSlant->\"Italic\"]\) = 40 nm",""},
LabelStyle->{FontSize -> 10},Evaluate[option]]; 
Fig13aBW=ListLinePlot[{pp1,pp2,pp3},PlotLegend->{" ⊥"," ∥","⊥+∥"},LegendPosition-> {.44,.16},LegendPosition-> {1-1.01,-0.2+.35},LegendSize->0.3,LegendShadow->{.01,-.01},
FrameLabel->{"\!\(\*
StyleBox[\"x\",\nFontSlant->\"Italic\"]\)","\!\(\*
StyleBox[\"y\",\nFontSlant->\"Italic\"]\)","\!\(\*
StyleBox[\"d\",\nFontSlant->\"Italic\"]\) = 40 nm",""},
LabelStyle->{FontSize -> 10},Evaluate[optB]]; 

λ1=0.55;               (*  initial wavelength  *)
λ2=0.55;                  (*  final wavelength  *)
λinc=0;                  (*  increment of wavelength  *)
r1=0.05;                  (*  initial radius of sphere in microns  *)
r2=0.05;                  (*  final radius of sphere in microns  *)
rinc = 0;              (*  increment of radius  *)
n1 = 1.33;                   (* refractive index of medium  *)
n2=0.57+I 2.45;    (*  refractive index of gold sphere  *)
θinc = 1;                   (*  increment of angle  *)     
{p1,pp1,p2,pp2,p3,pp3}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
Fig13b=ListLinePlot[{pp1,pp2,pp3},
FrameLabel->{"\!\(\*
StyleBox[\"x\",\nFontSlant->\"Italic\"]\)","\!\(\*
StyleBox[\"y\",\nFontSlant->\"Italic\"]\)","\!\(\*
StyleBox[\"d\",\nFontSlant->\"Italic\"]\) = 100 nm",""},
LabelStyle->{FontSize -> 10},Evaluate[option]] ;
Fig13bBW=ListLinePlot[{pp1,pp2,pp3},
FrameLabel->{"\!\(\*
StyleBox[\"x\",\nFontSlant->\"Italic\"]\)","\!\(\*
StyleBox[\"y\",\nFontSlant->\"Italic\"]\)","\!\(\*
StyleBox[\"d\",\nFontSlant->\"Italic\"]\) = 100 nm",""},
LabelStyle->{FontSize -> 10},Evaluate[optB]] ;

λ1=0.55;               (*  initial wavelength  *)
λ2=0.55;                  (*  final wavelength  *)
λinc=0;                  (*  increment of wavelength  *)
r1=0.08;                  (*  initial radius of sphere in microns  *)
r2=0.08;                  (*  final radius of sphere in microns  *)
rinc = 0;              (*  increment of radius  *)
n1 = 1.33;                   (* refractive index of medium  *)
n2=0.57+I 2.45;    (*  refractive index of gold sphere  *)
θinc = 1;                   (*  increment of angle  *)     
{p1,pp1,p2,pp2,p3,pp3}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
Fig13c=ListLinePlot[{pp1,pp2,pp3},
FrameLabel->{"\!\(\*
StyleBox[\"x\",\nFontSlant->\"Italic\"]\)","\!\(\*
StyleBox[\"y\",\nFontSlant->\"Italic\"]\)","\!\(\*
StyleBox[\"d\",\nFontSlant->\"Italic\"]\) = 160 nm",""},
LabelStyle->{FontSize -> 10},Evaluate[option]] ;
Fig13cBW=ListLinePlot[{pp1,pp2,pp3},
FrameLabel->{"\!\(\*
StyleBox[\"x\",\nFontSlant->\"Italic\"]\)","\!\(\*
StyleBox[\"y\",\nFontSlant->\"Italic\"]\)","\!\(\*
StyleBox[\"d\",\nFontSlant->\"Italic\"]\) = 160 nm",""},
LabelStyle->{FontSize -> 10},Evaluate[optB]] ;

λ1=0.55;              (*  initial wavelength  *)
λ2=0.55;                  (*  final wavelength  *)
λinc=0;                  (*  increment of wavelength  *)
r1=0.09;                  (*  initial radius of sphere in microns  *)
r2=0.09;                  (*  final radius of sphere in microns  *)
rinc = 0;              (*  increment of radius  *)
n1 = 1.33;                   (* refractive index of medium  *)
n2=0.57+I 2.45;    (*  refractive index of gold sphere  *)
θinc = 1;                   (*  increment of angle  *)     
{p1,pp1,p2,pp2,p3,pp3}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
Fig13d=ListLinePlot[{pp1,pp2,pp3},
FrameLabel->{"\!\(\*
StyleBox[\"x\",\nFontSlant->\"Italic\"]\)","\!\(\*
StyleBox[\"y\",\nFontSlant->\"Italic\"]\)","\!\(\*
StyleBox[\"d\",\nFontSlant->\"Italic\"]\) = 180 nm",""},
LabelStyle->{FontSize -> 10},Evaluate[option]] ;
Clear[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
Fig13dBW=ListLinePlot[{pp1,pp2,pp3},
FrameLabel->{"\!\(\*
StyleBox[\"x\",\nFontSlant->\"Italic\"]\)","\!\(\*
StyleBox[\"y\",\nFontSlant->\"Italic\"]\)","\!\(\*
StyleBox[\"d\",\nFontSlant->\"Italic\"]\) = 180 nm",""},
LabelStyle->{FontSize -> 10},Evaluate[optB]] ;
Clear[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];

Fig13=Framed[GraphicsGrid[{ {Fig13a, Fig13b},{Fig13c, Fig13d}}, ImageSize ->gg2, Frame -> False, Spacings -> {-10,0}],FrameMargins->1]
Fig13BW=Framed[GraphicsGrid[{ {Fig13aBW, Fig13bBW},{Fig13cBW, Fig13dBW}}, ImageSize ->gg2, Frame -> False, Spacings -> {-10,0}],FrameMargins->1];

Export[FileNameJoin[{base,"Fig9.13.JPG"}],Fig13BW,ImageResolution->600];

9.1.3 Quasistatic approximation

λ1 = 0.3;  (* initial wavelength *)
λ2 = 1.0;  (* final wavelength *)
λinc = 0.01;  (* increment of wavelength *)

num = Round[Abs[λ2 - λ1]/λinc] + 1;
Au=Table[{λ=λ1+(i-1)λinc,ϵ=(((Gold[λ]))^2);3Abs[Re[ϵ]/Im[ϵ]]},{i,num}];
Ag=Table[{λ=λ1+(i-1)λinc,ϵ=(((Silver[λ]))^2);3Abs[Re[ϵ]/Im[ϵ]]},{i,num}];
Cu=Table[{λ=λ1+(i-1)λinc,ϵ=(((Copper[λ]))^2);3Abs[Re[ϵ]/Im[ϵ]]},{i,num}];
Al=Table[{λ=λ1+(i-1)λinc,ϵ=(((Aluminum[λ]))^2);3Abs[Re[ϵ]/Im[ϵ]]},{i,num}];
Fig14=ListLinePlot[{Au,Ag,Cu,Al},
 LegendPosition -> {-.6, .2},
 PlotLegend->{"Au","Ag","Cu","Al"},
LegendSize ->{0.25,0.3},LegendShadow -> {.01, -.01},
 FrameLabel -> {"λ (μm)", "3|ϵ'/ϵ''|"},Evaluate[option]] 
Fig14BW=ListLinePlot[{Au,Ag,Cu,Al},
 LegendPosition -> {-.6, .2},
 PlotLegend->{"Au","Ag","Cu","Al"},
LegendSize ->{0.25,0.3},LegendShadow -> {.01, -.01},
 FrameLabel -> {"λ (μm)", "3|ϵ'/ϵ''|"},Evaluate[optB]] ;

Export[FileNameJoin[{base,"Fig9.14.JPG"}],Fig14BW,ImageResolution->600];
Clear[i,λ,λ1,λ2,λinc,num,Au,Ag,Cu,Al];

9.1.4 Ellipsoidal  particles

r1 = 0.0;   (* initial axes ratio b/a *)
r2 = 1.0;  (* final axes ratio b/a *)
rinc = 0.01;
  
(* increment of ratio *)
shape[r_]:=If[r==0,0,If[r==1,1/3,(-(r^2)/((((r^2)-1))^(3/2))) ArcSin[Sqrt[1-1/(r^2)]]+(r^2)/((r^2)-1)]];
num = Round[Abs[r2 - r1]/rinc] + 1;
Ax=Table[{r=r1+(i-1)rinc,shape[r]},{i,num}];
Ay=Table[{r=r1+(i-1)rinc,0.5-0.5shape[r]},{i,num}];

Fig16Prolate=ListLinePlot[{Ax,Ay},
 LegendPosition -> {.1, -.2},
 PlotLegend->{"Ax","Ay, Az"},
LegendSize -> {0.55,0.2},
 LegendShadow -> {.01, -.01}, FrameLabel -> {"Axes ratio (\!\(\*
StyleBox[\"b\",\nFontSlant->\"Italic\"]\)/\!\(\*
StyleBox[\"a\",\nFontSlant->\"Italic\"]\))", "Shape factor"},
Evaluate[option]] ;
Fig16ProlateBW=ListLinePlot[{Ax,Ay},
 LegendPosition -> {.1, -.2},
 PlotLegend->{"Ax","Ay, Az"},
LegendSize -> {0.55,0.2},
 LegendShadow -> {.01, -.01}, FrameLabel -> {"Axes ratio (\!\(\*
StyleBox[\"b\",\nFontSlant->\"Italic\"]\)/\!\(\*
StyleBox[\"a\",\nFontSlant->\"Italic\"]\))", "Shape factor"},
Evaluate[optB]] ;

r1 = 1.0;(* initial axes ratio b/a *)
r2 = 10.0;  (* final axes ratio b/a *)
rinc = 0.01;  (* increment of ratio *)
shape[r_]:=If[r==0,0,If[r==1,1/3,(-(r^2)/((((r^2)-1))^(3/2))) ArcSin[Sqrt[1-1/(r^2)]]+(r^2)/((r^2)-1)]];
num = Round[Abs[r2 - r1]/rinc] + 1;
Ax=Table[{r=r1+(i-1)rinc,shape[r]},{i,num}];
Ay=Table[{r=r1+(i-1)rinc,0.5-0.5shape[r]},{i,num}];

Fig16Oblate=ListLinePlot[{Ax,Ay}, FrameLabel -> {"Axes ratio (\!\(\*
StyleBox[\"b\",\nFontSlant->\"Italic\"]\)/\!\(\*
StyleBox[\"a\",\nFontSlant->\"Italic\"]\))", "Shape factor"},
Evaluate[option]] ;
Fig16OblateBW=ListLinePlot[{Ax,Ay}, FrameLabel -> {"Axes ratio (\!\(\*
StyleBox[\"b\",\nFontSlant->\"Italic\"]\)/\!\(\*
StyleBox[\"a\",\nFontSlant->\"Italic\"]\))", "Shape factor"},
Evaluate[optB]] ;

Fig16=Framed[GraphicsRow[{Fig16Prolate,Fig16Oblate},ImageSize->gr2,Frame->False,Spacings->{-30,0}],FrameMargins->1]
Fig16BW=Framed[GraphicsRow[{Fig16ProlateBW,Fig16OblateBW},ImageSize->gr2,Frame->False,Spacings->{-30,0}],FrameMargins->1];

Export[FileNameJoin[{base,"Fig9.16.JPG"}],Fig16BW,ImageResolution->600];
Clear[r1,r2,rinc,num,Ax,Ay];

aspect1 = 0.0;  (* initial aspect ratio a/b *)
aspect2 =30.0;  (* final aspect ratio a/b *)
ainc = 0.1;  (* increment of aspect ratio *)
ϵd= 1.0;    (* background dielectric constant *)
λ = 0.8;      (* wavelength in microns *)

field[a_,ϵm_,ϵd_]:=
If[a==0,1,If[a==1,Abs[3ϵm/(2ϵd+ϵm)],Abs[ϵm/(ϵd+(ϵm-ϵd)((-a/(((1-(a^2)))^(3/2))) ArcSin[Sqrt[1-(a^2)]]+1/(1-(a^2))))]]];

num = Round[Abs[aspect2 - aspect1]/ainc] + 1;
ϵm=(((Gold[λ]))^2);
Au=Table[{a=aspect1+(i-1)ainc,field[a,ϵm,ϵd]},{i,num}];
ϵm=(((Silver[λ]))^2);
Ag=Table[{a=aspect1+(i-1)ainc,field[a,ϵm,ϵd]},{i,num}];
ϵm=(((Aluminum[λ]))^2);
Al=Table[{a=aspect1+(i-1)ainc,field[a,ϵm,ϵd]},{i,num}];

Fig17=ListLinePlot[{Au,Ag,Al}, LegendPosition -> {.4, .2}, PlotLegend->{"Au","Ag","Al"},LegendSize -> {0.25,0.25}, LegendShadow -> {.01, -.01}, FrameLabel -> {"Aspect ratio (a/b)", "Field enhancement (ζ)"},
Evaluate[option]]
Fig17BW=ListLinePlot[{Au,Ag,Al}, LegendPosition -> {.4, .2}, PlotLegend->{"Au","Ag","Al"},LegendSize -> {0.25,0.25}, LegendShadow -> {.01, -.01}, FrameLabel -> {"Aspect ratio (a/b)", "Field enhancement (ζ)"},
Evaluate[optB]];

Export[FileNameJoin[{base,"Fig9.17.JPG"}],Fig17BW,ImageResolution->600];
Clear[aspect1,aspect2,ainc,num,field,ϵm,ϵd,λ,Au,Ag,Al];

aspect1 = 2.0;  (* initial aspect ratio a/b *)
aspect2 =30.0;  (* final aspect ratio a/b *)
ainc = 0.1;  (* increment of aspect ratio *)

field[a_]:=
1.5(a^2)(1-((-a/(((1-(a^2)))^(3/2))) ArcSin[Sqrt[1-(a^2)]]+1/(1-(a^2))));

num = Round[Abs[aspect2 - aspect1]/ainc] + 1;
γ=Table[{a=aspect1+(i-1)ainc,field[a]},{i,num}];

Fig18=ListLinePlot[γ, FrameLabel -> {"Aspect ratio (a/b)", "χ"},
Evaluate[option]]
Fig18BW=ListLinePlot[γ, FrameLabel -> {"Aspect ratio (a/b)", "χ"},
Evaluate[optB]];

Export[FileNameJoin[{base,"Fig9.18.JPG"}],Fig18BW,ImageResolution->600];
Clear[aspect1,aspect2,ainc,num,γ];

9.3 Nanoshells

λ1=0.3;                  (*  initial wavelength  *)
λ2=1.2;                  (*  final wavelength  *)
λinc=0.001;                  (*  increment of wavelength  *)
r2=0.02;                  (*  final radius of sphere in microns  *)
rinc = 0.0;              (*  increment of radius  *)
θinc = 1;                   (*  increment of angle  *)     
n2=-1;                   (*  refractive index of gold sphere  *)
n1 = 1.0;                   (* refractive index of medium  *)

r1=0.02;                  (*  radius of sphere in microns  *)
{QsAu1,QeAu1,QbAu1,QnAu1}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
r1=0.04;                  (*  radius of sphere in microns  *)
{QsAu2,QeAu2,QbAu2,QnAu2}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
r1=0.06;                  (*  radius of sphere in microns  *)
{QsAu3,QeAu3,QbAu3,QnAu3}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
r1=0.1;                  (*  radius of sphere in microns  *)
{QsAu4,QeAu4,QbAu4,QnAu4}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
r1=0.15;                  (*  radius of sphere in microns  *)
{QsAu5,QeAu5,QbAu5,QnAu5}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
r1=0.2;                  (*  radius of sphere in microns  *)
{QsAu6,QeAu6,QbAu6,QnAu6}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];

Fig20=ListLinePlot[{QeAu1,QeAu2,QeAu3,QeAu4,QeAu5,QeAu6},PlotLegend->{"40 nm","80 nm","120 nm","200 nm","300 nm","400 nm"},
LegendPosition-> {0.8,-.3+.1},LegendLabel -> "Core",
LegendSize->{0.5,0.65},LegendShadow->{.01,-.01},
FrameLabel->{"λ (μm)","(Subscript[Q, ext])"},Evaluate[option]]
Fig20BW=ListLinePlot[{QeAu1,QeAu2,QeAu3,QeAu4,QeAu5,QeAu6},PlotLegend->{"40 nm","80 nm","120 nm","200 nm","300 nm","400 nm"},
LegendPosition-> {0.8,-.3+.1},LegendLabel -> "Core",
LegendSize->{0.5,0.65},LegendShadow->{.01,-.01},
FrameLabel->{"λ (μm)","(Subscript[Q, ext])"},Evaluate[optB]];

Export[FileNameJoin[{base,"Fig9.20.JPG"}],Fig20BW,ImageResolution->600];
Clear[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];

λ1 = 0.5;          (* initial wavelength in microns *)
λ2 = 1.2;          (* final wavelength in microns  *)
λinc = 0.01;      (* increment of wavelength  *)
t2 = 0.005;        (*  final thickness of shell  *)
tinc = 0    ;     (*  increment of shell thickness or sphere radius *)
r1 = 0.06;        (*  initial radius of inner sphere  *)
r2 = 0.06;        (*  final radius of inner sphere    *)
rinc = 0.0;        (*  increment of radius  *)
n1 = 1.33;          (*  index of surrounding medium  *)
n2 = -1;              (*  index of gold  *)
n3 = 1.5;           (*  index of inner core  *)
θinc=1;

t1 = 0.005;       (*  initial thickness of shell  *){Qs1, Qe1, Qb1,QsE, QeE, QbE} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
t1 = 0.007;       (*  initial thickness of shell  *)
{Qs2, Qe2, Qb2,QsE, QeE, QbE} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
t1 = 0.01;       (*  initial thickness of shell  *)
{Qs3, Qe3, Qb3,QsE, QeE, QbE} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
t1 = 0.02;       (*  initial thickness of shell  *)
{Qs4, Qe4, Qb4,QsE, QeE, QbE} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
Fig22=ListLinePlot[{Qe1,Qe2,Qe3,Qe4}, LegendPosition -> {.85, -.4+.4}, PlotLegend->{"5 nm","7 nm","10 nm","20 nm"},
LegendSize -> {0.5,0.5},LegendShadow -> {.01, -.01},
 LegendLabel -> "Shell",FrameLabel -> {"λ (μm)",
 "(Subscript[Q, ext])"},
Evaluate[option]]
Fig22BW=ListLinePlot[{Qe1,Qe2,Qe3,Qe4}, LegendPosition -> {.85, -.4+.4}, PlotLegend->{"5 nm","7 nm","10 nm","20 nm"},
LegendSize -> {0.5,0.5},LegendShadow -> {.01, -.01},
 LegendLabel -> "Shell",FrameLabel -> {"λ (μm)",
 "(Subscript[Q, ext])"},
Evaluate[optB]];

Export[FileNameJoin[{base,"Fig9.22.JPG"}],Fig22BW,ImageResolution->600];
Clear[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc];

λ1 = 0.5;          (* initial wavelength in microns *)
λ2 = 1.2;          (* final wavelength in microns  *)
λinc = 0.01;      (* increment of wavelength  *)
t1 = 0.013;       (*  initial thickness of shell  *)
t2 = 0.013;        (*  final thickness of shell  *)
tinc = 0.02;     (*  increment of shell thickness or sphere radius *)
r2 = 0.055;        (*  final radius of inner sphere    *)
n1 = 1.33;          (*  index of surrounding medium  *)
n2 = -1;              (*  index of gold  *)
n3 = 1.5;           (*  index of inner core  *)
θinc=1;

r1 = 0.055;        (*  initial radius of inner sphere  *)
{Qs1, Qe1, Qb1,QsE, QeE, QbE} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
r1 = 0.075;        (*  initial radius of inner sphere  *)
{Qs2, Qe2, Qb2,QsE, QeE, QbE} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
r1 = 0.095;        (*  initial radius of inner sphere  *)
{Qs3, Qe3, Qb3,QsE, QeE, QbE} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
r1 = 0.115;        (*  initial radius of inner sphere  *)
{Qs4, Qe4, Qb4,QsE, QeE, QbE} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
Fig23=  ListLinePlot[{Qe1,Qe2,Qe3,Qe4},PlotLegend->{"55 nm","75 nm","95 nm","115 nm"}, LegendPosition -> {.85, 0},LegendSize -> {0.5,0.4},
 LegendShadow -> {.01, -.01}, LegendLabel -> "Core radius",
FrameLabel -> {"λ (μm)", "(Subscript[Q, ext])"},
Evaluate[option]]
Fig23BW=  ListLinePlot[{Qe1,Qe2,Qe3,Qe4},PlotLegend->{"55 nm","75 nm","95 nm","115 nm"}, LegendPosition -> {.85, 0},LegendSize -> {0.5,0.4},
 LegendShadow -> {.01, -.01}, LegendLabel -> "Core radius",
FrameLabel -> {"λ (μm)", "(Subscript[Q, ext])"},
Evaluate[optB]];

Export[FileNameJoin[{base,"Fig9.23.JPG"}],Fig23BW,ImageResolution->600];
Clear[λ1,λ2,λinc,t1,t2,tinc,r1,r2,n1,n2,n3,θinc];

λ1 = 0.5;          (* initial wavelength in microns *)
λ2 = 1.1;          (* final wavelength in microns  *)
λinc = 0.01;      (* increment of wavelength  *)
tinc = 0.02;     (*  increment of shell thickness or sphere radius *)
n1 = 1;          (*  index of surrounding medium  *)
n2 = -2;              (*  index of silver  *)
n3 = 1.5;           (*  index of inner core  *)
θinc=0;

r1 = 0.02;        (*  initial radius of inner sphere  *)
r2 = 0.02;        (*  final radius of inner sphere    *)
t1 = 0.002;       (*  initial thickness of shell  *)
t2 = 0.002;        (*  final thickness of shell  *)
{Qs1, Qe1, Qb1,QsE1, QeE1, QbE1} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
r1 = 0.03;        (*  initial radius of inner sphere  *)
r2 = 0.03;        (*  final radius of inner sphere    *)
t1 = 0.003;       (*  initial thickness of shell  *)
t2 = 0.003;        (*  final thickness of shell  *)
{Qs2, Qe2, Qb2,QsE2, QeE2, QbE2} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
r1 = 0.06;        (*  initial radius of inner sphere  *)
r2 = 0.06;        (*  final radius of inner sphere    *)
t1 = 0.006;       (*  initial thickness of shell  *)
t2 = 0.006;        (*  final thickness of shell  *)
{Qs3, Qe3, Qb3,QsE3, QeE3, QbE3} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
r1 = 0.09;        (*  initial radius of inner sphere  *)
r2 = 0.09;        (*  final radius of inner sphere    *)
t1 = 0.009;       (*  initial thickness of shell  *)
t2 = 0.009;        (*  final thickness of shell  *)
{Qs4, Qe4, Qb4,QsE4, QeE4, QbE4} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
r1 = 0.12;        (*  initial radius of inner sphere  *)
r2 = 0.12;        (*  final radius of inner sphere    *)
t1 = 0.012;       (*  initial thickness of shell  *)
t2 = 0.012;        (*  final thickness of shell  *)
{Qs5, Qe5, Qb5,QsE5, QeE5, QbE5} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
Fig24=  ListLinePlot[{Qe1,Qe2,Qe3,Qe4,Qe5},PlotLegend->{"10 nm","30 nm","60 nm","90 nm","120 nm"}, LegendPosition -> {.9, 0},LegendSize -> {0.5,0.5},
 LegendShadow -> {.01, -.01}, LegendLabel -> "Core radius", FrameLabel -> {"λ (μm)", "(Subscript[Q, ext])"},Evaluate[option]]
Fig24BW=  ListLinePlot[{Qe1,Qe2,Qe3,Qe4,Qe5},PlotLegend->{"10 nm","30 nm","60 nm","90 nm","120 nm"}, LegendPosition -> {.9, 0},LegendSize -> {0.5,0.5},
 LegendShadow -> {.01, -.01}, LegendLabel -> "Core radius", FrameLabel -> {"λ (μm)", "(Subscript[Q, ext])"},Evaluate[optB]];

Export[FileNameJoin[{base,"Fig9.24.JPG"}],Fig24BW,ImageResolution->600];
Clear[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc];

λ1 = 0.5;          (* initial wavelength in microns *)
λ2 = 1.1;          (* final wavelength in microns  *)
λinc = 0.01;      (* increment of wavelength  *)
tinc = 0.02;     (*  increment of shell thickness or sphere radius *)
n1 = 1;          (*  index of surrounding medium  *)
n2 = -2;              (*  index of silver  *)
n3 = 1.5;           (*  index of inner core  *)
θinc=0;

r1 = 0.02;        (*  initial radius of inner sphere  *)
r2 = 0.02;        (*  final radius of inner sphere    *)
t1 = 0.002;       (*  initial thickness of shell  *)
t2 = 0.002;        (*  final thickness of shell  *)
{Qs1, Qe1, Qb1,QsE1, QeE1, QbE1} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
r1 = 0.03;        (*  initial radius of inner sphere  *)
r2 = 0.03;        (*  final radius of inner sphere    *)
t1 = 0.003;       (*  initial thickness of shell  *)
t2 = 0.003;        (*  final thickness of shell  *)
{Qs2, Qe2, Qb2,QsE2, QeE2, QbE2} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
r1 = 0.06;        (*  initial radius of inner sphere  *)
r2 = 0.06;        (*  final radius of inner sphere    *)
t1 = 0.006;       (*  initial thickness of shell  *)
t2 = 0.006;        (*  final thickness of shell  *)
{Qs3, Qe3, Qb3,QsE3, QeE3, QbE3} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
r1 = 0.09;        (*  initial radius of inner sphere  *)
r2 = 0.09;        (*  final radius of inner sphere    *)
t1 = 0.009;       (*  initial thickness of shell  *)
t2 = 0.009;        (*  final thickness of shell  *)
{Qs4, Qe4, Qb4,QsE4, QeE4, QbE4} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
r1 = 0.12;        (*  initial radius of inner sphere  *)
r2 = 0.12;        (*  final radius of inner sphere    *)
t1 = 0.012;       (*  initial thickness of shell  *)
t2 = 0.012;        (*  final thickness of shell  *)
{Qs5, Qe5, Qb5,QsE5, QeE5, QbE5} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
Fig25 = ListLinePlot[{Qs1,Qs2,Qs3,Qs4,Qs5},PlotLegend->{"10 nm","30 nm","60 nm","90 nm","120 nm"}, LegendPosition -> {.85, -.1},LegendSize -> {0.5,0.5},
 LegendShadow -> {.01, -.01}, LegendLabel -> "Core radius", FrameLabel -> {"λ (μm)", "(Subscript[Q, sct])"},
 Evaluate[option]]
Fig25BW= ListLinePlot[{Qs1,Qs2,Qs3,Qs4,Qs5},PlotLegend->{"10 nm","30 nm","60 nm","90 nm","120 nm"}, LegendPosition -> {.85, -.1},LegendSize -> {0.5,0.5},
 LegendShadow -> {.01, -.01}, LegendLabel -> "Core radius", FrameLabel -> {"λ (μm)", "(Subscript[Q, sct])"},
 Evaluate[optB]];

Export[FileNameJoin[{base,"Fig9.25.JPG"}],Fig25BW,ImageResolution->600];
Clear[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc];

λ1 = 0.8;          (* initial wavelength in microns *)
λ2 = 0.8;          (* final wavelength in microns  *)
λinc = 0;      (* increment of wavelength  *)
tinc = 0;     (*  increment of shell thickness or sphere radius *)
n1 = 1.33;          (*  index of surrounding medium  *)
n2 = -1;              (*  index of gold shell *)
n3 = 1.5;           (*  index of inner core  *)
θinc=1;

r1 = 0.06;        (*  initial radius of inner sphere  *)
r2 = 0.06;        (*  final radius of inner sphere    *)
t1 = 0.011;       (*  initial thickness of shell  *)
t2 = 0.011;        (*  final thickness of shell  *)
{p1,pp1,p2,pp2} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 

Fig26= ListLinePlot[{p1,p2},PlotLegend->{" ⊥"," ∥"}, LegendPosition -> {.5, -.2},LegendSize ->{0.3,0.2},
 LegendShadow -> {.01, -.01}, LegendLabel -> "Pol.", FrameLabel -> {"Scattering angle (°)", "Irradiance"},
Evaluate[option]]
Fig26BW= ListLinePlot[{p1,p2},PlotLegend->{" ⊥"," ∥"}, LegendPosition -> {.3, -.2},LegendSize ->{0.3,0.2},
 LegendShadow -> {.01, -.01}, LegendLabel -> "Pol.", FrameLabel -> {"Scattering angle (°)", "Irradiance"},
Evaluate[optB]];

Export[FileNameJoin[{base,"Fig9.26.JPG"}],Fig26BW,ImageResolution->600];
Clear[λ1,λ2,λinc,tinc,n1,n2,n3,θinc,r1,r2,t1,t2];

λ1 = 0.8;          (* initial wavelength in microns *)
λ2 = 0.8;          (* final wavelength in microns  *)
λinc = 0;      (* increment of wavelength  *)
tinc = 0;     (*  increment of shell thickness or sphere radius *)
n1 = 1.33;          (*  index of surrounding medium  *)
n2 = -1;              (*  index of gold shell *)
n3 = 1.5;           (*  index of inner core  *)
θinc=1;

r1 = 0.135;        (*  initial radius of inner sphere  *)
r2 = 0.135;        (*  final radius of inner sphere    *)
t1 = 0.015;       (*  initial thickness of shell  *)
t2 = 0.015;        (*  final thickness of shell  *)
{p1,pp1,p2,pp2} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 

Fig27=ListLinePlot[{p1,p2},PlotLegend->{" ⊥"," ∥"}, LegendPosition -> {.4, .35}, LegendSize -> {0.3,0.2}, LegendShadow -> {.01, -.01}, LegendLabel -> "Pol.", FrameLabel -> {"Scattering angle (°)", "Irradiance"},Evaluate[option]]
Fig27BW=ListLinePlot[{p1,p2},PlotLegend->{" ⊥"," ∥"}, LegendPosition -> {.3, .35}, LegendSize -> {0.3,0.2}, LegendShadow -> {.01, -.01}, LegendLabel -> "Pol.", FrameLabel -> {"Scattering angle (°)", "Irradiance"},Evaluate[optB]];

Export[FileNameJoin[{base,"Fig9.27.JPG"}],Fig27BW,ImageResolution->600];
Clear[λ1,λ2,λinc,tinc,n1,n2,n3,θinc,r1,r2,t1,t2];

Figures

9.1 Nanoparticles

9.1.1 Introduction

9.1.2 Mie theory examples

9.1.3 Quasistatic  approximation

9.1.4 Ellipsoidal particles

9.1.5 Localized surface plasmon damping

9.3 Nanoshells

9.4 Other nanoparticle shapes

9.4.1  Nanodisks

9.4.2  Nanorods

9.4.3  Nanotriangles

9.5 Dual nanoparticles

9.5.1 Dual nanodisks

9.5.2 Bow ties

Exercises

1. In the introduction of this chapter, it was asserted that silver nanoparticles embedded in glass impart a yellow color while gold particles impart a red color.  Compute the extinction coefficient for 50 nm spherical silver and gold particles in a medium with a refractive index of 1.5 to see if this statement makes sense.  Note:  the wavelength range of yellow light is ~580-600 nm, and of red light is ~620-780 nm.

λ1=0.30;                  (*  initial wavelength  *)
λ2=0.7;                  (*  final wavelength  *)
λinc=0.001;                  (*  increment of wavelength  *)
θinc = 1;                   (*  increment of angle  *)   
r1=0.03;                  (*  initial radius of sphere in microns  *)
r2=0.03;                  (*  final radius of sphere in microns  *)
rinc = 0.0;              (*  increment of radius  *)
n1 = 1.5;                   (* refractive index of medium  *)  
n2=-1;                   (*  refractive index of gold sphere  *)
{QsAu,QeAu,QbAu,QnAu}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];
n2=-2;                   (*  refractive index of silver sphere  *)
{QsAg,QeAg,QbAg,QnAg}=SphereScattering[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];

Fig2=ListLinePlot[{QeAu,QeAg},PlotLegend->{"Au","Ag"},
LegendPosition-> {0.4,0.1},
LegendSize->0.35,
LegendShadow->{.01,-.01},
FrameLabel->{"λ (μm)","(Subscript[Q, ext])"},
Evaluate[option]] 
Clear[λ1,λ2,λinc,r1,r2,rinc,n1,n2,θinc];

2. For gold nanoshells on a silica core in vacuum, let the outer shell radius be 10% larger than the core radius.  As the core radius is varied from 20 nm to 120 nm what is the maximum extinction coefficient?  How does this compare to silver (see Fig. 9.24)?

λ1 = 0.5;          (* initial wavelength in microns *)
λ2 = 1.1;          (* final wavelength in microns  *)
λinc = 0.01;      (* increment of wavelength  *)
tinc = 0.02;     (*  increment of shell thickness or sphere radius *)
n1 = 1;          (*  index of surrounding medium  *)
n2 = -1;              (*  index of gold  *)
n3 = 1.5;           (*  index of inner core  *)
θinc=0;

r1 = 0.02;        (*  initial radius of inner sphere  *)
r2 = 0.02;        (*  final radius of inner sphere    *)
t1 = 0.002;       (*  initial thickness of shell  *)
t2 = 0.002;        (*  final thickness of shell  *)
{Qs1, Qe1, Qb1,QsE1, QeE1, QbE1} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
r1 = 0.03;        (*  initial radius of inner sphere  *)
r2 = 0.03;        (*  final radius of inner sphere    *)
t1 = 0.003;       (*  initial thickness of shell  *)
t2 = 0.003;        (*  final thickness of shell  *)
{Qs2, Qe2, Qb2,QsE2, QeE2, QbE2} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
r1 = 0.06;        (*  initial radius of inner sphere  *)
r2 = 0.06;        (*  final radius of inner sphere    *)
t1 = 0.006;       (*  initial thickness of shell  *)
t2 = 0.006;        (*  final thickness of shell  *)
{Qs3, Qe3, Qb3,QsE3, QeE3, QbE3} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
r1 = 0.09;        (*  initial radius of inner sphere  *)
r2 = 0.09;        (*  final radius of inner sphere    *)
t1 = 0.009;       (*  initial thickness of shell  *)
t2 = 0.009;        (*  final thickness of shell  *)
{Qs4, Qe4, Qb4,QsE4, QeE4, QbE4} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
r1 = 0.12;        (*  initial radius of inner sphere  *)
r2 = 0.12;        (*  final radius of inner sphere    *)
t1 = 0.012;       (*  initial thickness of shell  *)
t2 = 0.012;        (*  final thickness of shell  *)
{Qs5, Qe5, Qb5,QsE5, QeE5, QbE5} = ShellScattering[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc]; 
Fig24=  ListLinePlot[{Qe1,Qe2,Qe3,Qe4,Qe5},PlotLegend->{"10 nm","30 nm","60 nm","90 nm","120 nm"}, LegendPosition -> {.8+.1, -.3+.3},
 LegendSize -> 0.6,
 LegendShadow -> {.01, -.01},
 LegendLabel -> "Core radius", FrameLabel -> {"λ (μm)", "Extinction coefficient"},Evaluate[option]]
Clear[λ1,λ2,λinc, t1,t2,tinc,r1, r2,rinc, n1, n2, n3,θinc];

References

1. A. D. McFarland and R. P. Van Duyne.  Single silver nanoparticles as real-time optical sensors with zeptomole sensitivity.  Nano Lett. 3 (2003) 1057.
2. M. Moskovits.  Surface-enhanced spectroscopy.  Rev. Mod. Phys. 57 (1985) 783.
3. R. Elghanian, J. J. Storhoff, R. C. Mucic, R. L. Letsinger and C. A. Mirkin.  Selective colorimetric detection of polynucleotides based on the distance-dependent optical properties of gold nanoparticles.  Science 277 (1997) 1078.
4. X. Huang, I. H. El-Sayed, W. Qian, and M. A. El-Sayed.  Cancer cell imaging and photothermal therapy in the near-infrared region by using gold nanorods.  J. Am. Chem Soc. 128 (2006) 2115.
5. M. L. Brongersma, J. W. Hartman, and H. A. Atwater.  Electromagnetic energy transfer and switching in nanoparticle chain arrays below the diffraction limit.  Phys. Rev. B 62 (2000) R16356.
6. D. Ricard, P. Roussignol and C. Flytzanis.  Surface-mediated enhancement of optical phase conjugation in metal colloids.  Opt. Lett. 10 (1985) 511.
7. L. Novotny, R. X. Bian, and X. S. Xie.  Theory of nanometric optical tweezers.  Phys. Rev. Lett. 79 (1997) 645.
8. J. Aizprurua, P. Hanarp, D. S. Sutherland, M. Käll, G. W. Bryant, and F. J. García de Abajo.  Optical properties of gold nanorings.  Phys. Rev. Lett. 90 (2003) 057401.
9. K. -H. Su, Q. -H. Wei, X. Zhang, J. J. Mock, D. R. Smith, and S. Schultz.  Interparticle coupling effects on plasmon resonances of nanogold particles.  Nano Lett. 3 (2003) 1087.
10. H. Xu and M. Käll.  Surface-plasmon-enhanced optical forces in silver nanoaggregates.  Phys. Rev. Lett. 89 (2002) 246802.
11. G. Mie.  Beiträge zur optik trüber medien, speziell kolloidaler metallösungen.  Ann. d. Physik 25 (1908) 377.
12. C. F. Bohren and D. R. Huffman.  Absorption and Scattering of Light by Small Particles (New York:  John Wiley, 1983) chap. 4.
13. B. J. Messinger, K. U. von Raben, R. K. Chang and P. W. Barber.  Local fields at the surface of noble-metal microspheres.  Phys. Rev. B 24 (1981) 649.
14. P. B. Johnson and R. W. Christy.  Optical constants of copper and nickel as a function of temperature.  Phys. Rev. B  11 (1975) 1315.
15. U. Kreibig and M. Vollmer, Optical Properties of Metal Clusters (Berlin:  Springer-Verlag, 1995) p. 125.
16. J. Δ. Jackson.  Classical Electrodynamics, (2^nd) ed. (New York:  John Wiley, 1975), p. 100.
17. P. F. Liao and A. Wokaun.  Lightning rod effect in surface enhanced Raman scattering.  J. Chem Phys. 76 (1982) 751.
18. B. Lamprecht, J. R. Krenn, A. Leitner, and F. R. Aussenegg.  Resonant and off-resonant light-driven plasmons in metal nanoparticles studied by femtosecond-resolution third-harmonic generation.  Phys. Rev. Lett. 83 (1999) 4421.
19. J. Lehmann, M. Merschdorf, W. Pfeiffer, A. Thon, S. Voll, and G. Gerber.  Surface plasmon dynamics in silver nanoparticles studied by femtosecond time-resolved photoemission.  Phys. Rev. Lett. 85 (2000) 2921.
20. T. Zentgraf, A. Christ, J. Kuhl and H. Giessen.  Tailoring the ultrafast dephasing of quasiparticles in metallic photonic crystals.  Phys. Rev. Lett. 93 (2004) 243901.
21. I. D. Mayergoyz, Z. Zhang, and G. Miano.  Analysis of dynamics of excitation and dephasing of plasmon resonance modes in nanoparticles.  Phys. Rev. Lett. 98 (2007) 147401. 
22. U. Freibig and M. Vollmer, ibid, Sect. 2.2.2.
23. M. Perner, P. Bost, U. Lemmer, G. von Plessen, J. Feldmann, U. Becker, M. Mennig, M. Schmitt and H. Schmidt.  Optically induced damping of the surface plasmon resonance in gold colloids.  Phys. Rev. Lett. 78 (1997) 2192.
24. T. Klar, M. Perner, S. Grosse, G. von Plessen, W. Spirkl and J. Feldmann.  Surface-plasmon resonances in single metallic nanoparticles.  Phys. Rev. Lett. 80 (1998) 4249.
25. C. Sönnichsen, T. Franzl, T. Wilk, G. von Plessen, J. Feldmann, O. Wilson and P. Mulvaney.  Drastic reduction of plasmon damping in gold nanorods.  Phys. Rev. Lett. 88 (2002) 077402.
26. P. B. Johnson and R. W. Christy.  Optical constants of the noble metals.  Phys. Rev. B 6 (1972) 4370.
27. V. Pustovit and T. V. Shahbazyan.  Finite-size effects in surface-enhanced Raman scattering in noble-metal nanoparticles:  a semiclassical approach.  J. Opt. Soc. Am. A 23 (2006) 1369.
28. R. M. Cole, J. J. Baumberg, F. J. Garcia de Abajo, S. Mahajan, M. Abdelsalam and P. N. Bartlett.  Understanding plasmons in nanoscale voids.  Nano Lett. 7 (2007) 2094.
29. S. P. Apell, P. M. Echenique and R. H. Ritchie.  Sum rules for surface plasmon frequencies.  Ultramicrosc. 65 (1996) 53.
30. A. E. Neeves and M. H. Birnboim.  Composite structures for the enhancement of nonlinear-optical susceptibility.  J. Opt. Soc. Am. B 6 (1989) 787.
31. S. J. Oldenburg, J. B. Jackson, S. L. Westcott and N. J. Halas.  Infrared extinction properties of gold nanoshells.  Appl. Phys. Lett. 75 (1999) 2897.
32. G. D. Hale, J. B. Jackson, O. E. Shmakova, T. R. Lee and N. J. Halas.  Enhancing the active lifetime of luminescent semiconducting polymers via doping with metal nanoshells.   Appl. Phys. Lett. 78 (2001) 1502.
33. S. J. Oldenburg, S. L. Westcott, R. D. Averitt, and N. J. Halas.  Surface enhanced Raman scattering in the near infrared using metal nanoshell substrates.  J. Chem. Phys. 111 (1999) 4729. 
34. J. B. Jackson, S. L. Westcott, L. R. Hirsch, J. L. West and N. J. Halas.  Controlling the surface enhanced Raman effect via the nanoshell geometry.  Appl. Phys. Lett. 82 (2003) 257.
35. S. R. Sershen, S. L. Westcott, N. J. Halas and J. L. West.  Temperature-sensitive polymer-nanoshell composites for photothermally modulated drug delivery.  J. Biomed. Mater. Res. 51 (2000) 293.
36. L. R. Hirsch, R. J. Stafford, J. A. Bankson, S. R. Sershen, B. Rivera, R. E. Price, J. D. Hazle, N. J. Halas and J. L. West.  Nanoshell-mediated near-infrared thermal therapy of tumors under magnetic resonance guidance.  Proc. Natl. Acad. Sci. U.S.A. 100 (2003) 13549.
37. R. D. Averitt, D. Sarkar and N. J. Halas.  Plasmon resonance shifts of Au-coated (Subscript[Au, 2])S nanoshells:  insight into multicomponent nanoparticle growth.  Phys. Rev. Lett. 78 (1997) 4217.
38. S. J. Oldenburg, R. D. Averitt, S. L. Westcott and N. J. Halas.  Nanoengineering of optical resonances.  Chem. Phys. Lett. 288 (1998) 243.
39. S. Mohapatra, Y. K. Mishra, D. K. Avasthi, D. Kabiraji, J. Ghatak and S. Varma.  Synthesis of gold-silicon core-shell nanoparticles with tunable localized surface plasmon resonance.  Appl. Phys. Lett. 92 (2008) 103105.
40. A. B. R. Mayer, W. Gregner and  R. J. Wannemacher.  Preparation of silver-latex composites.  J. Phys. Chem B 104 (2000) 7278.
41. J. B. Jackson and N. J. Halas.  Silver nanoshells:  variations in morphologies and optical properties.  J. Phys. Chem B 105 (2001) 2743.
42. E. Prodan, C. Radloff, N. J. Halas and P. Nordlander.  A hybridization model for the plasmon response of complex nanostructures.  Science 302 (2003) 419.
43. C. Radloff and N. J. Halas.  Plasmonic properties of concentric nanoshells.  Nano Lett. 4 (2004) 1323.
44. H. Wang, F. Tam, N. K. Grady and N. J. Halas.  Cu nanoshells: effects of interband transitions on the nanoparticle plasmon resonance.  J. Phys. Chem B 109 (2005) 18218.
45. N. Halas.  The optical properties of nanoshells.  Optics and Photonics News (August, 2002) 26.
46. R. D. Averitt, S. L. Westcott and N. J. Halas.  Linear optical properties of gold nanoshells.  J. Opt. Soc. Am. B 16 (1999) 1824.
47. A. L. Aden and M. Kerker.  Scattering of electromagnetic waves from two concentric spheres.  J. Appl. Phys. 22 (1951) 1242.
48. T. Okamoto.  Near-field spectral analysis of metallic beads.  In Near-field Optics and Surface Plasmon Polaritons (Berlin:  Springer, 2001) p. 112 ff.
49. S. L. Westcott, J. B. Jackson, C. Radloff and N. J. Halas.  Relative contributions to the plasmon line shape of metal nanoshells.  Phys. Rev. B 66 (2002) 155431.
50. G. Raschke, S. Brogl, A. S. Susha, A. L. Rogach, T. A. Klar, J. Feldmann, B. Fieres, N. Petkov, T. Bein, A. Nichtl and K. Kurzinger.  Gold nanoshells improve single nanoparticle molecular sensors.  Nano Lett. 4 (2004) 1853.
51. C. L. Nehl, N. K. Grady, G. P. Goodrich, F. Tam, N. J. Halas and J. H. Hafner.  Scattering spectra of single gold nanoshells.  Nano Lett. 4 (2004) 2355.
52. S. J. Oldenburg, G. D. Hale, J. B. Jackson and N. J. Halas.  Light scattering from dipole and quadrupole nanoshell antennas.  Appl. Phys. Lett. 75 (1999) 1063.
53. S. J. Norton and T. Vo-Dinh.  Plasmon resonances of nanoshells of spheroidal shape.  IEEE Trans. Nanotechnol. 6 (2007) 627.
54. S. R. Nicewarner-Peña, R. G. Freeman, B. D. Reiss, L. He, D. J. Peña, I. D. Walton, R. Cromer, C. D. Keating, and M. J. Natan.  Submicrometer metallic barcodes.  Science 294 (2001) 137.
55. C. Sönnichsen, T. Franzl, T. Wilk, G. von Plessen, J. Feldmann, O. Wilson and P. Mulvaney.  Drastic reduction of plasmon damping in gold nanorods.  Phys. Rev. Lett. 88 (2002) 077402.
56. G. Sando, A. D. Berry, P. M. Campbell, A. P. Baronavski and J. C. Owrutsky.  Surface plasmon dynamics of high-aspect-ratio gold nanorods.  Plasmonics 2 (2007) 23.
57. H. J. Huang, C. P. Yu, H. C. Chang, K. P. Chiu, H. M. Chen, R. S. Liu and D. P. Tsai.  Plasmonic optical properties of a single gold nano-rod.  Opt. Exp. 15 (2007) 7132.
58. S.-M. Lee, Y. Jun, S.-N. Cho and J. Cheon.  Single-crystalline star-shaped nanocrystals and their evolution:  programming the geometry of nano-building blocks.  J. Am. Chem. Soc. 124 (2002) 11244.
59. C. L. Nehl, H. Liao and J. H. Hafner.  Optical properties of star-shaped gold nanoparticles.  Nano Lett. 6 (2006) 683.
60. C. Langhammer, Z. Yuan, I. Zoric and B. Kasemo.  Plasmonic properties of supported Pt and Pd nanostructures.  Nano Lett. 6 (2006) 833.
61. P. Hanarp, M. Kall anbd D. S. Sutherland.  Optical properties of short range ordered arrays of gold disks prepared by colloidal lithography.  J. Phys. Chem. B 107 (2003) 5768.
62. C. Langhammer, B. Kasemo and I. Zoric.  Absorption and scattering of light by Pt, Pd, Ag, and Au nanodisks:  absolute cross sections and branching ratios.  J. Chem. Phys. 126 (2007) 194702.
63. C. Langhammer, M. Schwind, B. Kasemo and I. Zoric.  Localized surface plasmon resonances in aluminum nanodisks.  Nano Lett. 8 (2008) 1461.       
64. J. Aizpurua, P. Hanarp, D. S. Sutherland, M. Kall, G. W. Bryant and F. J. Garcia de Abajo.  Optical properties of gold nanorings.  Phys. Rev. Lett. 90 (2003) 057401.
65. C. Aguirre, T. Kaspar, C. Radloff and N. J. Halas.  CTAB mediated reshaping of metallodielectric nanoparticles.  Nano Lett. 3 (2003) 1707.
66. H. Rochholz, N. Bocchio and M. Kreiter.  Tuning resonances on crescent-shaped noble-metal nanoparticles.  New J. Phys. 9 (2007) 53.
67. Y. Sun and Y. Xia.  Shape-controlled synthesis of gold and silver nanoparticles.  Science 298 (2002) 2176.
68. L. J. Sherry, S.-H. Chang, G. C. Schatz, R. P. Van Duyne, B. J. Wiley and Y. Xia.  Localized surface plasmon resonance spectroscopy of single silver nanocubes.  Nano Lett. 5 (2006) 2034.
69. C. L. Haynes and R. P. Van Duyne.  Plasmon-sampled surface-enhanced Raman excitation spectroscopy.  J. Phys. Chem. B 107 (2003) 7426.
70. A. J. Haes and R. P. Van Duyne.  A nanoscale optical biosensor:  sensitivity and selectivity of an approach based on the localized surface plasmon resonance spectroscopy of triangular silver nanoparticles.  J. Am. Chem. Soc. 124 (2002) 10596.
71. G. S. Métraux, Y. C. Cao, R. Jin and C. A. Mirkin.  Triangular nanoframes made of gold and silver.  Nano Lett. 3 (2003) 519.
72. J. Nelayah, M. Kociak, O. Stéphan, F. J. G. de Abajo, M. Tencé, L. Henrard, D. Taverna, I. Pastoriza-Santos, L. M. Liz-Marzán and C. Colliex.  Mapping surface plasmons on a single metallic nanoparticle.  Nat. Phys. 3 (2007) 348.
73. L. J. Sherry, R. Jin, C. A. Mirkin, G. C. Schatz and R. P. Van Duyne.  Localized surface plasmon resonance spectroscopy of single silver triangular nanoprisms.  Nano Lett. 6 (2006) 2060.
74. K. L. Kelly, E. Coronado, L. L. Zhao and G. C. Schatz.  The optical properties of metal nanoparticles:  the influence of size, shape, and dielectric environment.  J. Phys. Chem. B 107 (2003) 668.
75. J. C. Love, B. D. Gates, D. B. Wolfe, K. E. Paul and G. M. Whitesides.  Fabrication and wetting properties of metallic half-shells with submicron diameters.  Nano Lett. 2 (2002) 891.
76. M. Achermann, K. L. Shuford, G. C. Schatz, D. H. Dahanayaka, L. A. Bumm and V. I. Klimov.  Near-field spectroscopy of surface plasmons in flat gold nanoparticles.  Opt. Lett. 32 (2007) 2254.
77. A. K. Sheridan, A. W. Clark, A. Glidle, J. M. Cooper and D. R. S. Cumming.  Multiple plasmon resonances from gold nanostructures.  Appl. Phys. Lett. 90 (2007) 143105.
78. P. Hanarp, M. Käll, and D. S. Sutherland.  Optical properties of short range ordered arrays of nanometer gold disks prepared by colloidal lithography.  J. Phys. Chem. B  107 (2003) 5768.
79. C. Langhammer, Z. Yuan, I. Zoric, and B. Kasemo.  Plasmonic properties of supported Pt and Pd nanostructures.  Nano Lett. 6 (2006) 833.
80. H. Fredriksson, Y. Alaverdyan, A. Dmitriev, C. Langhammer, D. S. Sutherland, M. Zäch and B. Kasemo.  Hole-mask colloidal lithograpy. Adv. Mater. 19 (2007) 4297.
81. T. Inagaki, K. Kagami, and E. T. Arakawa.  Photoacoustic observation of nonradiative decay of surface plasmons in silver.  Phys. Rev. B 24 (1981) 3644.
82. T. Inagaki, K. Kagami, and E. T. Arakawa.  Photoacoustic study of surface plasmons in metals.  Appl. Opt. 21 (1982) 949.
83. K. Sturm, W. Schülke, and J. R. Schmitz.  Plasmon-Fano resonance inside the particle-hole excitation spectrum of simple metals and semiconductors.  Phys. Rev. Lett. 68 (1992) 228.
84. M. van Exter and A. Lagendijk.  Ultrashort surface-plasmon and phonon dynamics.  Phys. Rev. Lett. 60 (1988) 49.
85. S.-S. Chang, C.-W. Shih, C.-D. Chen, W.-C. Lai, and C. R. C. Wang.  The shape transition of gold nanorods.  Langmuir 15 (1999) 701.
86. M. Pelton, M. Liu, K. C. Toussaint, Jr., H. Y. Kim, G. Smith, J. Pesic, P. Guyot-Sionnest and N. F. Scherer.  Plasmon-enhanced optical trapping of individual metal nanorods.  Proc. SPIE 6644 (2007) 66441C.
87. R. Jin, Y. Cao, C. A. Mirkin, K. L. Kelly, G. C. Schatz and J. G. Zheng.  Photoinduced conversion of silver nanospheres to nanoprisms.  Science 294 (2001) 1901.
88. W. Rechberger, A. Hohenau, A. Leitner, J. R. Krenn, B. Lambrecht, and F. R. Aussenegg.  Optical properties of two interacting gold nanoparticles.  Opt. Commun. 220 (2003) 137.
89. P. Jain, W. Huang, and M. A. El - Sayed.  On the universal scaling behavior of the distance decay of plasmon coupling in metal nanoparticles pairs: a plasmon ruler equation.  Nano Lett. 7 (2007) 2080.