NEAR-FIELD OPTICS

Remove["Global`*"]

Highlights

Abstract

Introduction

Off[General::spell,General::spell1];
Needs["Graphics`ParametricPlot3D`"]
Needs["Graphics`Arrow`"];
Needs["Graphics`PlotField`"];
Needs["Algebra`Trigonometry`"];
Needs["Algebra`ReIm`"];
Needs["LinearAlgebra`CrossProduct`"];
Needs["Calculus`VectorAnalysis`"];

tera  = 10^12  ; micro = 10^-6;
giga  = 10^9   ; nano  = 10^-9; 
mega  = 10^6   ; pico  = 10^-12;
kilo  = 10^3   ; femto = 10^-15; 
milli = 10.^-3 ; atto  = 10^-18;

NSub = {
c   -> 2.99792458 10^8  , e     -> 1.6021773 10^-19,
h   -> 6.6260755 10.^-34, hbar  -> 1.05457 10.^-34 ,
K   -> 1.380658 10^-23  , me    -> 9.109389 10^-31 ,
mu0 -> 4. Pi 10^-7      , eps0  -> 8.85419 10^-12   };

ESi   =0.179 tera; rhoSi   =2330 ; nuSi   =0.28;
ESiO2 =0.06  tera; rhoSiO2 =2200 ; nuSiO2 =0.3;
ESi3N4=0.15  tera; rhoSi3N4=3100 ; nuSi3N4=0.24;
EAu   =0.08  tera; rhoAu   =19300; nuAu   =0.42;
ENi   =0.207 tera; rhoNi   =8900 ; nuNi   =0.41;
EW    =0.345 tera; rhoW    =19300; nuW    =0.28;
EPt   =0.17  tera; rhoPt   =21400; nuPt   =0.39;
EIr   =0.52  tera; rhoIr   =22500; nuIr   =0.3 ;
EHOPG =0.2   tera; rhoHOPG =2250 ; nuHOPG =0.3 ;
EPol  =0.001 tera; rhoPol  =1000 ; nuPol  =0.3 ;

xx     =0.8;
EPtIr  =xx 0.8 EPt+(1-xx) EIr;
rhoPtIr=xx rhoPt  +(1-xx) rhoIr;
nuPtIr =0.3;

opt= {
PlotRange  ->All,
AxesOrigin ->{0,0},
Frame      ->True,
GridLines  ->Automatic,           
DefaultFont->{"Times-Roman",8},
PlotStyle  ->{
Dashing[{.0001}],
Dashing[{.0500}],
Dashing[{.0200}]}};

opt2 = {
PlotRange  ->All,
AxesOrigin ->{0,0}, 
DefaultFont->{"Times-Roman",8},
FrameLabel ->{"x","y"},
Frame      ->True,
GridLines  ->Automatic};

num   = 45;
optBW = {
PlotPoints ->num,
PlotRange  ->All,
Axes       ->True, 
DefaultFont->{"Times-Roman",8},
Axes       ->True,
Frame      ->True,
Mesh       ->False};

optColor = {
PlotPoints ->num,
PlotRange  ->All,
Axes       ->True, 
DefaultFont->{"Times-Roman",8},
Axes       ->True,
Frame      ->True,
Mesh       ->False,
ColorFunction->Hue};

Operating Conditions

lambda = 0.5 micro;
a      = lambda/10;

case = {
k      -> 2 Pi/lambda,
w      -> 2 Pi c/lambda/.NSub,
m0     -> -2 a^3/(3 Pi)};

Clear[atn];
atn[x_,y_]:=If[x==0&&y==0,0,ArcTan[x,y]]

opt1 = optBW;
(*opt1 = optColor;*)

Far-Field Solution

Vector Potential

Ad = 
m0/r^2 (1-I k r)*
Exp[I(k r-w t)] Sin[theta]*
{0,0,1};

Ad0 = m0/r^2 Sin[theta] Sqrt[1+(k r)^2]{0,0,1};
psi = k r+ArcTan[1,-k r]-w t;
Ad  = Ad0 Exp[I psi];

Electric Field

Ed = -1/c D[Ad,t];

Ed0 = -I Simplify[Ed Exp[-I psi]]/.case/.NSub;

Clear[Edx];
Edx[x_,y_,z_] :=
Sin[phi] Ed0[[3]]/.
theta->ArcCos[y/Sqrt[x^2+y^2+z^2]]/.
r    ->Sqrt[x^2+y^2+z^2]/.
phi  ->atn[x,z];

Clear[Edz];
Edz[x_,y_,z_] :=
Cos[phi] Ed0[[3]]/.
theta->ArcCos[y/Sqrt[x^2+y^2+z^2]]/.
r    ->Sqrt[x^2+y^2+z^2]/.
phi  ->atn[x,z];

Clear[Edy];
Edy[x_,y_,z_] := 0;

Electric Vector Field in the $x-z$ Plane

Fig16x1 = 
PlotVectorField[Evaluate[{
Edx[x a,a,z a],
-Edz[x a,a,z a]}],
{x,-2,2,.5},
{z,-2,2,.5},
AxesLabel  ->{"x [a]","z [a]"},
DefaultFont->{"Times-Roman",8},
AxesOrigin->{-2.24,0},
PlotRange  ->All,
Axes       ->True,
Frame      ->True];

Electric Vector Field in the $x-y$ Plane

Fig16x2 = 
PlotVectorField[Evaluate[{
Edx[x a,y a,a],
Edy[x a,y a,a]}],
{x,-2,2,.5},
{y,-2,2,.5},
PlotRange  ->All,
Axes       ->True,
AxesOrigin->{-2.18,0},
AxesLabel  ->{"x [a]","y [a]"},
DefaultFont->{"Times-Roman",8},
Frame      ->True];

Electric Field in the $x-y$ Plane

Fig16x3 = 
DensityPlot[Evaluate[
Abs[Edx[x a,a y,a]]^2],
{x,-2,2},{y,-2,2},
AxesLabel->{"x [a]","y [a]"},
Evaluate[opt1]];

Fig16x4 = 
DensityPlot[Evaluate[
Abs[Edz[x a,y a,a]]^2],
{x,-2,2},{y,-2,2},
AxesLabel->{"x [a]","y [a]"},
Evaluate[opt1]];

Magnetic Field

SetCoordinates[Spherical];
Hd = Curl[Ad];

Hd01 = Simplify[Hd Exp[-I psi]];

Hd0 = {
Hd01[[1]], 
Sqrt[(I ComplexExpand[Hd01][[2,1]])^2+
(ComplexExpand[Hd01][[2,2]])^2],0}/.case/.NSub;

Poynting Vector and Intensity

Sd = c/(4 Pi) Cross[Ed0,Hd0]/.NSub;

Id = Sqrt[Dot[Sd,Sd]];

Fig16x5 = 
SphericalPlot3D[Evaluate[
Id/.r->a 200],
{theta,Pi/100,Pi},{phi,Pi/2,-Pi/2},
Axes        ->True,
Boxed       ->True,
PlotRange   ->{{0,5},{-5,5},{-5,5}},
ViewVertical->{1,0,0},
DefaultFont->{"Times-Roman",8},
AxesLabel   ->{"z","x","y"}];

Clear[Idc];
Idc[x_,y_,z_] :=
Id/.
theta->ArcTan[y,Sqrt[x^2+z^2]]/.
r    ->Sqrt[x^2+z^2+y^2]/.
phi  ->atn[z,x];

Intensity in the $x-z$ Plane

Fig16x6 = 
DensityPlot[Evaluate[
Idc[x a,a,z a]],
{x,-2,2.01},{z,-2,2.01},
AxesLabel->{"x [a]","z [a]"},
Evaluate[opt1]];

Intensity in the $x-y$ Plane

Fig16x7 = 
DensityPlot[Evaluate[
Idc[x a,y a,a]],
{x,-2,2.01},{y,-2,2.01},
AxesLabel->{"x [a]","y [a]"},
Evaluate[opt1]];

Near-Field Solution

sol=
Solve[
{x^2+y^2==a^2 (1-u^2) (1+v^2),
z==a u v},
{u,v}];

uu=sol[[3,1,2]];
vv=sol[[3,2,2]];

Electric Field

clear[Eax];
Eax[x_,y_,z_] =
k z-2/Pi k a u (1+v ArcTan[v]+
1/(3(u^2+v^2))+(x^2-y^2)/(3a^2(u^2+v^2)(1+v^2)^2))/.
{u->uu,v->vv}/.case;

Clear[Eay];
Eay[x_,y_,z_] =
-4 k x y u/(3Pi a (u^2+v^2)(1+v^2)^2)/.
{u->uu,v->vv}/.case;

Clear[Eaz];
Eaz[x_,y_,z_] =
4 k x v/(3 Pi(u^2+v^2)(1+v^2))/.
{u->uu,v->vv}/.case;

Clear[Ea];
Ea[x_,y_,z_] = {
Eax[x,y,z],
Eay[x,y,z],
Eaz[x,y,z]};

Electric Field in the $x-y$ Plane

Fig16x8=Table[ 
DensityPlot[Evaluate[
Abs[Eax[x a, y a, 10^jj/1000 a]^2]],
{x,-2,2},{y,-2,2},
AxesLabel->{"x [a]","y [a]"},
PlotLabel->{jj},
Evaluate[opt1],
DisplayFunction->Identity],
{jj,1,3,1}];

Show[GraphicsArray[{
Fig16x8[[1]],Fig16x8[[2]],Fig16x8[[3]]}]];

tbl1 = 
TableForm[{
N[Abs[Eax[0,0,a/100]]^2],
N[Abs[Eax[0,0,a/10]]^2],
N[Abs[Eax[0,0,a/1]]^2]},
TableHeadings-> {
{"a/100","a/10","a"}}]

Fig16x9 = Table[ 
DensityPlot[Evaluate[
Abs[Eay[x a, y a, 10^jj/1000 a]^2]],
{x,-2,2},{y,-2,2},
AxesLabel->{"x [a]","y [a]"},
PlotLabel->{jj},
Evaluate[opt1],
DisplayFunction->Identity],
{jj,1,3,1}];

Show[GraphicsArray[{
Fig16x9[[1]],Fig16x9[[2]],Fig16x9[[3]]}]];

Fig16x10 = Table[ 
DensityPlot[Evaluate[
Abs[Eaz[x a, y a, 10^jj/1000 a]^2]],
{x,-2,2},{y,-2,2},
AxesLabel->{"x [a]","y [a]"},
PlotLabel->{jj},
Evaluate[opt1],
DisplayFunction->Identity],
{jj,1,3,1}];

Show[GraphicsArray[{
Fig16x10[[1]],Fig16x10[[2]],Fig16x10[[3]]}]];

Magnetic Field

Clear[Hax]
Hax[x_,y_,z_] =
-4 x y v/(Pi a^2 (u^2+v^2) (1+v^2)^2)/.
{u->uu,v->vv}/.case;

Clear[Hay];
Hay[x_,y_,z_] =
1-2/Pi (ArcTan[v]+v/(u^2+v^2))+
2(x^2-y^2) v/(Pi a^2 (u^2+v^2) (1+v^2)^2)/.
{u->uu,v->vv}/.case;

Clear[Haz];
Haz[x_,y_,z_] =
-4 a y u/(Pi a^2 (u^2+v^2) (1+v^2))/.
{u->uu,v->vv}/.case;

Clear[Ha];
Ha[x_,y_,z_] = {Hax[x,y,z],Hay[x,y,z],Haz[x,y,z]};

Poynting Vector and Intensity

SetCoordinates[Cartesian];

Clear[Sa];
Sa[x_,y_,z_] := c/(4 Pi)*
CrossProduct[Ea[x,y,z],Ha[x,y,z]]/.NSub

Clear[Ia];
Ia[x_,y_,z_] :=
N[Sqrt[Dot[Sa[x,y,z],Sa[x,y,z]]]]

Fig16x11 = Table[ 
ListDensityPlot[Table[
Ia[x a,y a,10^jj/1000 a],
{x,-2,2,.5},{y,-2,2,.5}],
AxesLabel->{"x [a]","y [a]"},
PlotLabel->{jj},
Evaluate[Take[opt1,-6]],
DisplayFunction->Identity],
{jj,1,3,1}];

Show[GraphicsArray[{
Fig16x11[[1]],Fig16x11[[2]],Fig16x11[[3]]}]];

tbl2 = 
TableForm[{
N[micro Ia[0,0,a/100]],
N[micro Ia[0,0,a/10]],
N[micro Ia[0,0,a/1]]},
TableHeadings-> {
{"a/100","a/10","a"}}]

Discussion of the Models

Electric Field

Fig16x12 =
Plot[Evaluate[{
Abs[Edx[a x,a/1000,a]],
Abs[Eax[a x,a/1000,a]]}],
{x,-2,2},
FrameLabel ->{"x [a]","Ex [a.u]"},
Evaluate[opt]];

Intensity

Fig16x13 =
Plot[Evaluate[{
micro Idc[a x,a/1000,a],
micro Ia[a x,a/1000,a]}],
{x,-2,2.1},
FrameLabel ->{"x [a]","Intensity [a.u]"},
Evaluate[opt]];

Figure 16x11