Motion of a Charge in the Field of a Magnetic Dipole

Franz Michel,  Bernhard Schnizer

Institut fuer Theoretische Physik
		Technische Universitaet Graz
		A-8010  Graz, Austria
		Projekt Multimediale Lehre des
		Zentralen Informatikdienstes der TU Graz

		2002-12-13

Begin["dipol`"] ;Clear["\
dipol`*"] ;

Introduction: Equations of Motion

The equations of motion for a point particle (masse m, charge e) in the field of a magnetic dipole 
( poles strength  g, located at   (Overscript[r, →])  =  (0,0, ± p0) =  ± (Subscript[(Overscript[p, →]), 0])  )  are:

		\(TraditionalForm\`\(\(	\)\(m \(\(\(\ \)\(r\&→\)\)\& .. \)\ \  = \ \ \(e\ \((\(\(\ \)\(r\&→\)\)\& . \ \ ×\ \ B\&→\ )\)\  = \ \ \(\(e\)\(\ \)\(g\)\(\ \)\(μ\_0\)\(\ \ \)\(\(\(\ \)\(r\&→\)\)\& . \)\(\ \ \)\(×\)\(\ \ \)\((\ \(r\&→\ \  - \ \ p\&→\_0\)\/\(\(\((x\^2\  + \ \ y\^2\ \  + \ \((z\ \  - \ \ p\_0)\)\^2)\)\^\(3/2\)\)\(\ \)\)\  - \ \ \(r\&→\ \  + \ \ p\&→\_0\)\/\(\(\((x\^2\  + \ \ y\^2\ \  + \ \((z\ \  + \ \ p\_0)\)\^2)\)\^\(3/2\)\)\(\ \)\)\ )\)\(\ \ \ \)\)\)\)\)\)
					
They are solved numerically and this solution is used to plot the trajectory.

Start Button for Animations

start = 
        ButtonBox[
          StyleBox["Run Animation", "Title", 
            FontFamily -> "Helvetica"], 
          ButtonFunction :> 
            FrontEndExecute[{SelectionMove[EvaluationNotebook[], Before, 
                    ButtonCell] SelectionMove[EvaluationNotebook[], 
                    Previous, CellGroup] SelectionAnimate[
                    SelectedNotebook[], AnimationDisplayTime -> zeit]}], 
          ButtonEvaluator -> Automatic, 
          Active -> True];

Plotting field lines

Basic Formulae

Bispherical coordinates are the most convenient ones for plotting the field lines.

<< Calculus`VectorAnalysis`

SetCoordinates[
      Bispherical[θ, η, ϕ]]
CoordinateRanges[]

lc = {CoordinatesToCartesian[{θ, η, ϕ}]  /.  ϕ  ->  
            0 , 
CoordinatesToCartesian[{θ, η, ϕ}]  /.  ϕ  ->  π}

st  =  
        Table[θ  ->   v , {v, π/
              nc,  π (nc  -  1)/nc, π/nc}];

SetOptions[ParametricPlot, 
        DisplayFunction  ->  Identity];

tg  =  
        Table[ParametricPlot[
            Evaluate[{lc[[1, {1, 3}]], lc[[2, {1, 3}]] } /.  
                st[[k]]], {η, -5, 
              5},  
AspectRatio  ->  Automatic,  
            PlotStyle  ->  Hue[(k  -  1)/(nc  -  1)]], {k, 
            Length[st]}];

Show[tg,  DisplayFunction  -> $DisplayFunction,  
        Ticks  ->  None , 
ImageSize  ->  
          850,  Epilog  ->  {PointSize[ .015], Point[{0, 1}], 
            Point[{0, -1}], 
Thickness[ .007], 
            Line[{{0, 1}, {0, -1}}]}];

Parameters

Equations of Motion and Their Solution

Basic Formulae

Position of positive magnetic point pole.  2 po is the length of the dipole.

p0  =  {0,  0,  po};

vr  =  {x[t], y[t], z[t]};        
vv  =  
      D[vr, t];
vb  =  D[vv, t];
v1  =  
        vr  -  
          p0;                                 (* vector from positive pole to charge *) 
v2  =  
        vr  +  
          p0;                                 (* vector from negative pole to charge *) 
b1  =  Sqrt[v1 . v1];
b2  =  Sqrt[v2 . v2];
f1  =  
        v1/b1^3 ;                             (* magnetic field due to positive pole at   p0 *) 
f2  =  
        v2/b2^3 ;                             (* magnetic field due to positive pole at  -  
        p0 *) 
ft  =  α (f1  -  f2);

Equations of motion

eqm  =  Thread[vb  ==  Cross[vv, ft]]

Initial data:

vr0  =  {x0, y0, z0};
vv0  =  {vx0, vy0, vz0};
ad  =  
      Thread[(Join[vr, vv]  /.  t  ->  0)  ==  
          Join[vr0, vv0]]

Parameters

Solving the Equations of Motion

eqn  =  Join[eqm, ad]  //. sp ;

sr  =  {x, y, z};

so  =  
      NDSolve[eqn, sr,  {t, 0, tm},  MaxSteps  ->   10000] // 
        Flatten

Plotting the Field and the Trajectory

Plotting the dipole

fl  =  
        Graphics3D[{PointSize[ .02], Point[p0], Point[-p0], Hue[ .3], 
              Thickness[ .01], Line[{p0, -p0}]}] //. sp;

Plotting parts of trajectory

SetOptions[ParametricPlot3D, 
        DisplayFunction  ->  Identity];

First period

tt  =  100;
t1  =  
      t  /. FindRoot[
          Evaluate[vr[[3]] /. so], {t, tt}]
p1  =  
        ParametricPlot3D[
          Join[Evaluate[vr /. so], {Hue[0]}],  {t, 0, 
            t1},  
  AspectRatio  ->  
            Automatic,  PlotPoints  ->  1000];

Second period

tt  =  200;
t2  =  
      t  /. FindRoot[
          Evaluate[vr[[3]] /. so], {t, tt}]
p2  =  
        ParametricPlot3D[
          Join[Evaluate[vr /. so], {Hue[0.6]}],  {t, t1, 
            t2},  
  AspectRatio  ->  
            Automatic,  PlotPoints  ->  1000];

Third period

tt  =  300;
t3  =  
      t  /. FindRoot[
          Evaluate[vr[[3]] /. so], {t, tt}]
p3  =  
        ParametricPlot3D[
          Join[Evaluate[vr /. so], {Hue[0.1]}],  {t, t2, 
            t3},  
  AspectRatio  ->  
            Automatic,  PlotPoints  ->  1000];

Displaying all three periods together:

tr  =  
        Show[p1, p2, p3, fl,  PlotRange  ->  All,  
          Axes  ->  
            False,  
DisplayFunction  ->  $DisplayFunction , ImageSize  ->  500 ];

Animated Motion of a Charge

Positions

np  =  1;  tl  =  Range[0, t3, 1/np];

tp  =  
        Table[Point[
            Evaluate[(vr /. so) ] /.  t  ->  tl[[k]]], {k, 
            Length[tl]}];

Show[tr, Graphics3D [tp]];

Animation

Table[
        Show[tr, Graphics3D [tp[[k]]]], {k, 
          Length[tp]}];
SelectionMove[EvaluationNotebook[], After, 
        EvaluationCell];
SelectionMove[EvaluationNotebook[], Next, CellGroup];
FrontEndTokenExecute["OpenCloseGroup"];
start /. zeit -> 0.1 // DisplayForm