ELECTROSTATICS

Remove["Global`*"]

Highlights

Abstract

Introduction

Needs["Graphics`PlotField`"];

Off[General::spell,General::spell1];

Needs["Graphics`Arrow`"];

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}]}};

opt1 = { 
DefaultFont->{"Times-Roman",8},
FrameLabel ->{"x","y"},
Frame      ->True};

Isolated Point Charge

Clear[vpoint];
vpoint[x_,y_] := 1/Sqrt[x^2+y^2];

discGraphics = 
Graphics[
Disk[{0,0},.1],
AxesOrigin ->{0,0},
AspectRatio->Automatic,
PlotLabel  ->"(a)",
opt1];

fieldGraphics = 
PlotGradientField[
-vpoint[x,y],
{x,-1.5,1.5},{y,-1.5,1.5},
ScaleFunction  ->(2+Log[1+#]&),
ScaleFactor    ->0.5,
PlotPoints     ->8,
DisplayFunction->Identity];

discfieldGraphics = 
Show[
{discGraphics,fieldGraphics},
DisplayFunction->Identity];

potentialGraphics = 
ContourPlot[
vpoint[x,y],
{x,-2.5,2.5},{y,-2.5,2.5},
PlotPoints     ->20,
Contours       ->7,
PlotLabel      ->"(b)",
ContourShading->False,
DisplayFunction->Identity,
Evaluate[opt1]];

Fig1 = Show[GraphicsArray[{
discfieldGraphics,
potentialGraphics}]];

Point Charge and Plane

case1 = {d->1};
Qi    = -1;
si    = -2 d/.case1;
Clear[vreal,vimage];
vreal[x_,y_]  :=  1/Sqrt[x^2+y^2];
vimage[x_,y_] := -1/Sqrt[x^2+(y-si)^2];

discGraphics = 
Graphics[{
Disk[{0,0},.1],
Rectangle[{-2,-1.1},{2,-1}]},
AspectRatio->Automatic,
PlotLabel  ->"(a)",
opt1];

fieldGraphics = 
PlotGradientField[
-(vreal[x,y]+vimage[x,y]),
{x,-1.5,1.5},{y,-1,1},
ScaleFunction  ->(2+Log[1+#]&),
ScaleFactor    ->0.5,
PlotPoints     ->8,
DisplayFunction->Identity];

discfieldGraphics = 
Show[
{discGraphics,fieldGraphics},
DisplayFunction->Identity];

potentialGraphics = 
ContourPlot[
vreal[x,y]+vimage[x,y],
{x,-2.5,2.5},{y,-3.5,1.5},
PlotPoints     ->20,
Contours       ->10,
PlotLabel      ->"(b)",
ContourShading->False,
DisplayFunction->Identity,
Evaluate[opt1]];

Fig2 = Show[GraphicsArray[{
discfieldGraphics,
potentialGraphics}]];

Point Charge and Sphere

case2 = {R->1.,d->2.};
Qi    = -R/d/.case2;
si    = d-R^2/d/.case2;
Clear[vreal,vimage];
vreal[x_,y_]  := 1/Sqrt[x^2+y^2];
vimage[x_,y_] := Qi/Sqrt[x^2+(y+si)^2];

discGraphics = 
Graphics[{
Disk[{0,0},.1],
Disk[{0,-d/.case2},1]},
AspectRatio->Automatic,
PlotLabel  ->"(a)",
opt1];

fieldGraphics = 
PlotGradientField[
-(vreal[x,y]+vimage[x,y]),
{x,-2,2},{y,-3.5,1.5},
ScaleFunction  ->(2+Log[1+#]&),
ScaleFactor    ->0.5,
PlotPoints     ->8,
DisplayFunction->Identity];

discfieldGraphics = 
Show[
{discGraphics,fieldGraphics},
DisplayFunction->Identity];

potentialGraphics = 
ContourPlot[
vreal[x,y]+vimage[x,y],
{x,-2.5,2.5},{y,-3.5,1.5},
PlotPoints     ->30,
Contours       ->10,
PlotLabel      ->"(b)",
ContourShading->False,
DisplayFunction->Identity,
Evaluate[opt1]];

Fig3 = Show[GraphicsArray[{
discfieldGraphics,
potentialGraphics}]];

Isolated Sphere

Clear[Csp];
Csp[R_] := 4 Pi eps0 R;

case3 = {R->10 nano};
Print["Csp = ",N[Csp[R]/atto]/.case3/.NSub," [aF]"]

Sphere and Plane

case4 = {R->1.,d->2.};
nn    = 5;

Position of Charges Inside the Sphere

h      = Array[q,nn];
h[[1]] = 0;
h[[2]] = R^2/(2 d);
tbl1a  = Simplify[
Table[h[[n]] =
R^2/(2 d-h[[n-1]]),
{n,3,nn}]];
p11 = Table[h[[n]],{n,1,nn}];
p1  = Take[p11,3]

Magnitude of Charges Inside the Sphere

g      = Array[q,nn];
g[[1]] = 1;
g[[2]] = r;
tbl2a   = Simplify[
Table[g[[n]] = r g[[n-1]]/(1-r g[[n-1]]/g[[n-2]]),
{n,3,nn}]];
c11 = Table[g[[n]],{n,1,nn}]/.r->R/(2 d);
c1 = Take[c11,3]

Position of Charges Outside the Sphere

p22 = Table[2 d-h[[n]],{n,1,nn}];
p2  = Take[p22,3]

Magnitude of Charges Outside the Sphere

c22 = Table[-g[[n]],{n,1,nn}]/.r->R/(2 d);
c2  = Take[c22,3]

Table of Coefficients

Clear[tbl1];
scx = N[c11/.case4,3]; sc4 = Take[scx,4]; 
spx = N[p11/.case4,3]; sp4 = Take[spx,4];
pcx = N[c22/.case4,3]; pc4 = Take[pcx,4];
ppx = N[p22/.case4,3]; pp4 = Take[ppx,4];

tbl1 = 
TableForm[{sc4,sp4,pc4,pp4},
TableHeadings-> {
{"c1","p1","c2","p2"},
{"1","2","3","4"}}]

Potential and Field

Clear[vsp];
vsp[x_,y_] := Sum[
scx[[i]]/Sqrt[x^2+(y+spx[[i]])^2]+
pcx[[i]]/Sqrt[x^2+(y+ppx[[i]])^2],
{i,1,nn}];

discGraphics = 
Graphics[{
Disk[{0,0},1],
Rectangle[{-2,-2.1},{2,-2}]},
AspectRatio->Automatic,
PlotLabel  ->"(a)",
opt1];

fieldGraphics = 
PlotGradientField[
-vsp[x,y],
{x,-1.5,1.5},{y,-2,1.5},
ScaleFunction  ->(2+Log[1+#]&),
ScaleFactor    ->0.5,
PlotPoints     ->8,
DisplayFunction->Identity];

discfieldGraphics = 
Show[
{discGraphics,fieldGraphics},
DisplayFunction->Identity];

potentialGraphics = 
ContourPlot[
vsp[x,y],
{x,-2.5,2.5},{y,-3.5,1.5},
PlotPoints     ->30,
Contours       ->10,
PlotLabel      ->"(b)",
ContourShading->False,
DisplayFunction->Identity,
Evaluate[opt1]];

Fig4 = Show[GraphicsArray[{
discfieldGraphics,
potentialGraphics}]];

Potential Along the Axis of Symmetry

RR = R/.case4;
dd = d/.case4;
Fig5 = Plot[
vsp[0,y],
{y,-dd,-RR},
FrameLabel ->{"y","V(0,y)"},
Evaluate[opt]];

Capacitance

Clear[Csp];
Csp := Function[
{R,s},
{alpha = Log[1+s/R+Sqrt[s^2/R^2+2 s/R]];
cap2   = 1+Sum[Sinh[alpha]/Sinh[n alpha],{n,2,100}];
N[4. Pi eps0 R cap2]//.NSub}[[1]]];

Clear[Cspa];
Cspa[R_,s_] := N[2 Pi eps0 R(2+R/s)]/.NSub;

Example

case5 = {R->10 nano,s->50 nano};
RR = R/.case5;
ss = s/.case5;
Timing[
Print["Csp  = ",Csp[RR,ss]/atto," [aF]"]]
Timing[
Print["Cspa = ",Cspa[RR,ss]/atto," [aF]"]]

Fig6 = Plot[
Cspa[RR,s nano]/Csp[RR,s nano],
{s,RR/nano,ss/nano},
FrameLabel ->{"s [nm]","Cspa/Csp"},
Evaluate[opt]];

Two Spheres

Capacitance: Exact Solution

Clear[R1,R2,s];
Css := Function[
{R1,R2,s},
{temp1 = 4. Pi eps0 R1//.NSub;
 temp2 = 4. Pi eps0 R2//.NSub;
 R1sq  = R1 R1;
 R2sq  = R2 R2;

q  = 1.; da = s + R1 + R2; d1 = da;
qt = 0.; qqt = 0.;

Do[
   {qq  = -R2 q/da;
	dd  = d1 - R2sq/da;
	qt  = qt + q;
	qqt = qqt + qq;
	q   = -R1 * qq /dd;
	da   = d1 - R1sq/dd;},
{i,1,100}];

c11 =  qt temp1;
c12 = qqt temp1;

q  = 1.; da = s + R1 + R2; 
qt = 0.; qqt = 0.;
Do[{
    qq = -R1 q/da;
    dd = d1 - R1sq/da;
    qt = qt + q;
    q  = -R2 qq/dd;
    da  = d1 - R2sq/dd;},
{i,1,100}];

c22 = qt temp2;
c21 = c12;
                     
N[(c11 c22 - c12 c21)/(c11 + c12 + c21 + c22)]}[[1]]]

Example

case6 = {R1->100 nano,R2->10 nano,s->1 nano};
RR1 = R1/.case6;
RR2 = R2/.case6;
ss  = s/.case6;
Print["Css = ",Css[RR1,RR2,ss]/atto," [aF]"]

Matrix Elements

tbl2 = TableForm[{
{N[c11/atto],N[c12/atto]},
{N[c21/atto],N[c22/atto]}},
TableHeadings->{
{"j=1","j=2"},
{"i=1","i=2"}}]

Capacitance: Approximate Solution

Clear[Cssa];
Cssa = Function[
{R1,R2,s},
{R11 = Max[{R1,R2}];
R22 = Min[{R1,R2}];
rb  = 1./(1. + R22/R11);
N[rb Csp[R22,rb s]]}[[1]]];

Example

Timing[
Print["Css  = ",Css[RR1,RR2,ss]/atto," [aF]"]]
Timing[
Print["Cssa = ",Cssa[RR1,RR2,ss]/atto," [aF]"]]

Fig7 = Plot[
Cssa[RR1,RR2,s nano]/Css[RR1,RR2,s nano],
{s,.1,100.},
FrameLabel  ->{"s [nm]","Cssa/Css"},
PlotRange   ->{.9,1.},
Evaluate[opt]];

Electrostatic Force

Clear[F];
F[R_,s_,v_] = -1/2 D[Cspa[R,s],s] v^2;

Fig8 = Plot[
F[100 nano,s nano,1]/nano,
{s,10,100},
FrameLabel  ->{"s [nm]","F [nN]"},
Evaluate[opt]];