(* semiinfi.nb *)

Needs["Graphics`ComplexMap`"]

Clear[V,a,b,c,d,e,frc,nec,nef,
      neq,tcf,tcb,mult,tex1]
Off[General::spell1]

		           	(* Input data *)
			
V = 1 (* Potential on upper plate in Volt, 
		 potential on lower plate (= median plane) = 0 *);
		 
a = 3 (* Distance between plates in mm : 
         for 1.1 10^-5 < a < 3 10^5 *); 
                
nec = 30 (* Number of lines used in 3D-image of 
            el. field`s absolute value: absef3D3 *);
            
nef = 11 (* Number of field lines between plates;
            variable needed in pictures: 
            eflandeqpl3, efl3, absefateflu3, absefateflv3, 
            efuatefl3, efvatefl3 *);
            
neq = 7  (* Number of equipotential lines between plates;
            variable needed in pictures: 
            eflandeqpl3, eqpll3, absefateqplu3, 
            absefataqplv3, efuateqpl3, efvateqpl3 *);

tcf = .005 (* Relative thickness of field or equipotential
              lines in color diagrams *);
               
tcb = .011 (* Relative thickness of lines representing 
		   	  boundaries, electrodes, pole shoes, etc. *);
                  
                   (* Global variables *)
                  
mult = V 1000/a; tex1 = "E[V/m]";           
                    
                     (* Functions *) 
                                
cm3[z_,t_] = a/Pi (1 + z + t); (* Conformal map *)
ef3[t_] = -I mult/Conjugate[1 + t]; (* Function of el. field *)
af31[t_] = 1/(1 + t); 
af32[z_] = cm3[z,t[z]]; 
t[z_] = Exp[z];

         (* Variables needed in all pictures *)
         
b = -5.7155 ; d = cm3[-c,t[-c]]; e = cm3[b,t[b]] ; 
frc = FindRoot[e - cm3[-c,-t[-c]],{c,-2,-3}]; c = c/.frc ;
 
Print[SequenceForm[
"Distance between infinitely thin plate and median plane: a = ", a, " mm."]]

(* Section with black & white graphics *)

Print["Electric field on top side of upper plate "]
ppuptop3 = ParametricPlot[
   {Re[cm3[I Pi - x, -t[-x]]],-mult af31[-t[-x]]}, 
   {x, 0.001, c}, 
    PlotRegion -> {{.01,.99},{.01,.99}},
    PlotRange -> {All, {0,5 mult}}, 
    AxesLabel -> {"u[mm]", tex1}]; 
    
(* Calculations for a later picture *)    
ppuptop3inv = ParametricPlot[
   {-Re[cm3[I Pi-x,-t[-x]]],-mult af31[-t[-x]]},{x,0.001,c}, 
    PlotRange -> {{-e, 0},5 mult {-1, 1}}, 
    DisplayFunction -> Identity]; 
    
Print["Electric field on bottom side of upper plate "]
ppupbot3 = ParametricPlot[
   {Re[cm3[I Pi-x,-t[-x]]],-mult af31[-t[-x]]},{x,-b,-0.1}, 
    PlotRange -> {{e, 0}, {0, -5 mult}}, 
    AxesLabel -> {"u[mm]", tex1}]; 
    
Print["Electric field lines & equipotentials "]
eflandeqpl3 = CartesianMap[
    cm3[#1,t[#1]] & ,{-c,b,(b+c)/(nef-1)},{0,Pi,Pi/(neq-1)}, 
    PlotRegion -> {{.01,.99},{.01,.99}},
    PlotRange -> All, 
    AxesLabel -> {"u[mm]", "v[mm]"}, 
    Epilog -> 
      {Thickness[tcb], 
      {Line[{{d, 0}, {e, 0}}], Line[{{e, a}, {0, a}}]}}]; 

Print["Electric field on median plate "]
ppdown3 = ParametricPlot[
    {cm3[x, t[x]], -mult af31[t[x]]},{x, -c, b}, 
    PlotRegion -> {{.01,.99},{.01,.99}},
    AxesLabel -> {"u[mm]", tex1}];

Print["Electric field of upper plate without units"]
Show[GraphicsArray[{ppupbot3,ppuptop3inv}, GraphicsSpacing -> 0 ]/. 
      (Ticks -> Automatic) -> (Ticks -> None) /.
      (AxesLabel -> {"u[mm]",tex1}) -> 
      (AxesLabel -> {None,None}) /. 
      (PlotRange -> {{e,0},{0,-5 mult}}) -> 
      (PlotRange -> {{e,0},5 mult{-1,1}}),
      Epilog -> { Text["Bottom u ->                  Top <- u",{1,.4}]}] ;
      
Print["Electric fieldlines & equipotentials without units"]
Show[eflandeqpl3, Axes -> None] ;

Print["Electric field of median plate without units"] 
Show[ppdown3, 
     PlotRegion -> {{.01,.99},{.01,.99}},
     Ticks -> None, 
     AxesLabel -> {"u","E"}] ;

(* Section with colour graphics *)

Print[ColumnForm[{"Absolute value of electric field",
                  "in 3D seen from above upper plate"}]]
absef3D3 = ParametricPlot3D[Evaluate[
 		     {Re[af32[x + I y]],Im[af32[x + I y]],
 	          Abs[ef3[t[x + I y]]]}], {x,-c,b},{y,0,Pi},
  	      	 PlotRegion -> {{.01,.99},{.01,.99}},
  	      	 PlotRange -> {0, 2 mult},
  		     AxesLabel -> {"u[mm]","v[mm]","|E| [V/m]"}, 
  		     BoxRatios -> {1.5,1,1}, 
  		     PlotPoints -> nec,
		     ViewPoint -> {-1.07, 2.3, 1.9}];

hli =	Table[Hue[N[1 - (k-1)/(neq - 2)]],{k,neq-2}] ;
hlc =  	Join[{GrayLevel[0]}, Reverse[hli], {GrayLevel[0]}] ;
hlt = 	Table[Thickness[tcf],{k,neq}] ;
styleq = Transpose[{hlc,hlt}] ;

Print[neq," Equipotential lines (including plates)"]
eqpl3 = ParametricPlot[ Evaluate[Table[
     {Re[af32[(k-1)/(neq-1)Pi I + x]],
      Im[af32[(k-1)/(neq-1)Pi I + x]]},
		{k,neq}]],{x,-c,b}, 
		PlotRegion -> {{.01,.99},{.01,.99}},
		AxesLabel -> {"u [mm]","v [mm]"},
		PlotStyle -> styleq,
     	Epilog -> {Thickness[tcb],
                 { Line[{{d,0},{e,0}}], Line[{{e,a},{0,a}}] }} ];

Print[ColumnForm[
	  {"Absolute electric fields on equipotential lines as",
      "functions of u (colours match the colours above)"}]]
absefateqplu3 = ParametricPlot[ Evaluate[Table[
	{Re[af32[(k-1)/(neq-1)Pi I + x]],
    Abs[ef3[t[(k-1)/(neq-1)Pi I + x]]]},
		{k,neq}]],{x,-c,b}, 
		PlotRegion -> {{.01,.99},{.01,.99}},
		AxesLabel -> {"u [mm]","|E(u)| [V/m]"},
		PlotStyle -> styleq] ;

Print[ColumnForm[
	  {"Absolute electric fields on equipotential lines as", 
	  "functions of v (colours match the colours above)"}]]
absefateqplv3 = ParametricPlot[ Evaluate[Table[
	{Im[af32[(k-1)/(neq-1)Pi I + x]],
    Abs[ef3[t[(k-1)/(neq-1)Pi I + x]]]},
		{k,neq}]],{x,-c,b}, 
		PlotRegion -> {{.01,.99},{.01,.99}},
		AxesLabel -> {"v [mm]","|E(v)| [V/m]"},
		PlotStyle -> styleq] ;

Print[ColumnForm[
	  {"Electric fields (v-component) on equipotential lines as", 
	  "functions of u (colours match the colours above)"}]]
efvateqpl3 = ParametricPlot[ Evaluate[Table[
	 {Re[af32[(k-1)/(neq-1)Pi I + x]],
	  Im[ef3[t[(k-1)/(neq-1)Pi I + x]]]},
		{k,neq}]],{x,-c,b}, 
		PlotRegion -> {{.01,.99},{.01,.99}},
		AxesLabel -> {"u [mm]","\!\(E\_v\) [V/m]"},
		PlotStyle -> styleq] ;

Print[ColumnForm[
	  {"Electric fields (u-component) on equipotential lines as",
	   "functions of v (colours match the colours above)"}]]
efuateqpl3 = ParametricPlot[ Evaluate[Table[
	 {Im[af32[(k-1)/(neq-1)Pi I + x]],
	  Re[ef3[t[(k-1)/(neq-1)Pi I + x]]]},
		{k,neq}]],{x,-c,b}, 
		PlotRegion -> {{.01,.99},{.01,.99}},
		AxesLabel -> {"v [mm]","\!\(E\_u\) [V/m]"},
		PlotStyle -> styleq] ;

hli =	Table[Hue[N[1 - (k-1)/(nef - 2)]],{k,nef-2}]; 
hlc = 	Join[{GrayLevel[0]}, Reverse[hli], {GrayLevel[0]}];
hlt = 	Table[Thickness[tcf],{k,nef}];
stylfe = Transpose[{hlc,hlt}];
dx  = 	c + b;

Print[nef," Electric field lines"]
efl3 = ParametricPlot[ Evaluate[Table[
	 {Re[af32[-c + (k-1)/(nef-1) dx  + I y]], 
 	  Im[af32[-c + (k-1)/(nef-1) dx  + I y]]},
		{k,nef}]],{y,0,Pi}, 
		PlotRange -> All,
		PlotRegion -> {{.01,.99},{.01,.99}},
		AxesLabel -> {"u [mm]","v [mm]"},
		PlotStyle -> stylfe,
     	Epilog -> {Thickness[tcb],
                 { Line[{{d,0},{e,0}}],        
                   Line[{{e,a},{0,a}}] }} ] ;

Print[ColumnForm[
	  {"Absolute electric fields on field lines as",
	  "functions of u (colours match the colours above)"}]]
absefateflu3 = ParametricPlot[ Evaluate[Table[
 	 {Re[af32[-c + (k-1)/(nef-1) dx + I y]],
	 Abs[ef3[t[-c + (k-1)/(nef-1) dx + I y]]]},
		{k,nef}]],{y,0,Pi}, 
		PlotRegion -> {{.01,.99},{.01,.99}},
		AxesLabel -> {"u [mm]","|E(u)| [V/m]"},
		PlotStyle -> stylfe] ;

Print[ColumnForm[
	  {"Absolute electric fields on field lines as", 
	  "functions of v (colours match the colours above)"}]]
absefateflv3 = ParametricPlot[ Evaluate[Table[
 	 {Im[af32[-c + (k-1)/(nef-1) dx + I y]],
	 Abs[ef3[t[-c + (k-1)/(nef-1) dx + I y]]]},
		{k,nef}]],{y,0,Pi}, 
		PlotRegion -> {{.01,.99},{.01,.99}},
		AxesLabel -> {"v [mm]","|E(v)| [V/m]"},
		PlotStyle -> stylfe] ;

Print[ColumnForm[
	  {"Electric fields (v-component) on field lines as", 
	  "functions of u (colours match the colours above)"}]]
efvatefl3 = ParametricPlot[ Evaluate[Table[
 	{Re[af32[-c + (k-1)/(nef-1) dx + I y]],
	 Im[ef3[t[-c + (k-1)/(nef-1) dx + I y]]]},
		{k,nef}]],{y,0,Pi}, 
		PlotRegion -> {{.01,.99},{.01,.99}},
		AxesLabel -> {"u [mm]","\!\(E\_v\) [V/m]"},
		PlotStyle -> stylfe] ;

Print[ColumnForm[
	  {"Electric fields (u-component) on field lines as", 
	  "functions of v (colours match the colours above)"}]]
efuatefl3 = ParametricPlot[ Evaluate[Table[
 	{Im[af32[-c + (k-1)/(nef-1) dx + I y]],
	 Re[ef3[t[-c + (k-1)/(nef-1) dx + I y]]]},
		{k,nef}]],{y,0,Pi},
		PlotRegion -> {{.01,.99},{.01,.99}},
		AxesLabel -> {"v [mm]","\!\(E\_u\) [V/m]"},
		PlotStyle -> stylfe] ;

Clear[up,vp,frefp,efp,error,imefp,s]

up = 8; vp = 5;
(* above you may insert values of the Cartesian coordinates 
u and v (in mm) of a point of your choice within the area 
shown in Fig. eflandeqpl3 *) 
				
Which[vp < 0,     Print[vp " = vp < 0  is unphysical.",⚠],
True,

frefp = FindRoot[af32[efp] - (up + I vp), 
        {efp,  1.2  Sign[vp - a] + I Pi/2 }];
efp = efp /. frefp; imefp = N[Im[efp]];
error = Abs[N[af32[efp] - (up + I vp)]];

Which[0 > imefp, 	 Print["The calculated point doesn`t fulfill the requirements!",\[FreakedSmiley]], 
     imefp > N[Pi], Print["The calculated point doesn`t fulfill the requirements!",\[FreakedSmiley]],
     error > 10^-6, Print["The error for the calculated point is too big !",☹], 
     True, Print["The calculated point fulfills the requirements.",☺];
     
valef = ef3[t[efp]]//N;

Print[
SequenceForm["At point (up,vp) = (",up ,", ",vp ,") mm : ",
"Potential = ",N[V Im[efp]/Pi], " V;"]];
Print[SequenceForm["Field (\!\(E\_u\),\!\(E\_v\)) = (", 
Re[valef],", ", Im[valef],") V/m "]]; 

caleqpl3 = ParametricPlot[Evaluate[
	  {Re[af32[I Im[efp] + s]], Im[af32[I Im[efp] + s]]}], 
	  {s,-c,b}, 
	  PlotStyle -> Dashing[{.01}], 
	  DisplayFunction -> Identity];
	  
calefl3 = ParametricPlot[Evaluate[
	  {Re[af32[Re[efp] + I s]], Im[af32[Re[efp] + I s]]}],
	  {s,0,Pi}, DisplayFunction -> Identity];
	  
Show[caleqpl3,calefl3, 
      DisplayFunction -> $DisplayFunction,
      PlotRegion -> {{.01,.99},{.01,.99}},
	  PlotRange -> All,
	  AxesLabel -> {"u [mm]","v [mm]"},
 	  Epilog -> { Thickness[tcb],
                { Line[{{d,0},{e,0}}], 
                  Line[{{e,a},{0,a}}] }}]; ]]  

On[General::spell1]