Saturday, April 28, 2012

High Current Magnetic Lens Simulation

This is a follow up related to my post from yesterday about electron beam magnetic lenses. I was curious about different configurations and operating parameters, so I ran the simulations again and increased the current by a factor of ten.  Everything else in the simulations are exactly the same.  I'll present the resultant graphs below with some notes.  The associated files can be found here.
Cross section of the field through the coil axis
Cross section of the field just above the plane of the coil
The field plots above look exactly the same as the low current plots.  This is because they are plots of relative intensity, just to show the shape of the field.  If the current in the coil is multiplied by ten, the field will be ten times stronger, but will still have the same shape.

Focussing and diffusion

From the graph above it appears that the field is so strong that the electrons are focussed just before entering the coil but still keep spiralling along at different speeds with each beam approaching the axis twice more at slightly different positions.
x velocity through the simulation
We can see here once again that the electron's x velocity decreases significantly as is transits the lens, most of this is however regained after passing through the lens.


velocity magnitude

Ignoring any errors from the numerical simulation, we can see that the velocity of the electron remains constant throughout its interaction with the lens.

Electrons spiralling along the x axis

In this plot we can see the electron as it spirals along the axis.  I do find it surprising that the electrons don't spiral around the axis, instead they just keep kissing the axis.  I think if this was longer solenoid with a more uniform field there may be different results.

Electron trajectories

The high current simulation shows some interesting results.  The electrons that start furtherest from the axis loose a significant portion of their x velocity and are flung off away from the axis after spiralling around next to it.  I am slightly puzzled that the electrons don't spiral around the axis.  What little information I could find in textbooks on this topic have the electrons spiralling around the axis.  They weren't simulations though, and were more of just a conceptual look at the problem.  I may need a more uniform field to see that effect.  Hmmm.

Friday, April 27, 2012

Electron Beam Magnetic Lens Simulation

Solenoid field intensity

A while back I explained how cathode ray tubes (CRT) work as part of another blog post, and I understood all of it except for one thing, how on earth did the magnetic focussing coil work.  So I made it my goal to simulate the physics and see what came out of it.  GNU Octave turned out to be handy for this as I already knew MATLAB.

To model the problem an electron is fired towards a current loop parallel to the loop's axis and the trajectory is observed.  The loop is in the y-z plane centred at the origin with the electrons fired in the positive x direction.  First things first, this is a basic simulation, it models the electron as a classical particle with a negative charge.  I'm not taking relativity into account, even though the speed of the electrons in the simulation are 20% of the speed of light after being accelerated through 10 kV, so relativistic effects should be starting to become apparent.  The simulation is for a single electron, so there are no repulsive forces between electrons in the beam causing blooming.  Radiation given off by accelerating charges is not taken into account.

I am going to use a modified version of the formula for the magnetic field of a current loop that I calculated in a previous post.  These equations are for a current loop in the x-y plane and will need to be rotated for a loop in the y-z plane.

Magnetic field equations

 These equations along with an ODE solver in Octave are used to calculate the trajectory and plot some graphs and analyse the results.  I won't explain the code it's pretty simple.  You start with the initial position and velocity for an electron, you put that into the solver along with the equations describing the forces acting on the electron and you get an answer.  The comments in the code should go a long way to describing how it works.  Most of it is just for displaying the graphs, so although it looks hard it really isn't.  The files are available here.

The magnetic field generated by the current loop acts upon the electron with the Lorentz force as shown below.

Lorentz force
By combing that with Newtons second law of motion we obtain the formula for acceleration.
 
electron acceleration
This shows that the acceleration of the electron is proportional to the cross product of the velocity and the magnetic field which means it's accelerated perpendicular to its velocity and will tend to travel in arcs.  You can get a conceptual idea of how the coil focuses electrons with a thought experiment.  Firstly the electron comes in parallel to the coil axis and experiences no force from the component of the magnetic field parallel to the axis.  The component of the field perpendicular to the axis accelerates it tangentially and gives it a velocity perpendicular to the axis.  This velocity now causes a radial force moving the electron towards the axis focussing it.

A real life coil used to focus the beam would have many turns on it, but to start with I am going to simulate it as if all the turns are in the same place.  So instead of 1000 turns carrying 500mA, I can simulate one turn with 500A.

The parameters for the simulation were chosen to be close to reality, but I wasn't too concerned with the electron velocity, or accurate dimensions for the coil as this is just to investigate the lensing effect.  So for the first simulation I selected 500 A for the coil current and a radius of 3cm.  Three particles starting 1mm, 3mm, and 9mm from the coil axis are simulated.

It may be useful to start with a graphical representation of the field.

Cross section of the field through the coil axis
Cross section of the field just above the plane of the coil

Evidence of focussing

In the graph above we can see that the three particles are are forced towards the axis near the coil at x equals 0, and are focussed near the point x equals 0.15 meters.  There is evidence of aberration at the focal point just like an optical lens, but at one point all beams come within approximately 0.2 mm of the coil axis.

x velocity through the simulation
An interesting property of an electron moving in a magnetic field is that it doesn't change speed.  The velocity will vary but the speed will remain constant, this means it can turn but not speed up or down.  Therefore the magnitude of the velocity vector should remain constant.  For a speed change, work must be done on, or done by the electron, and if we have a look at the definition of work it will become obvious why the speed stays constant.


As the force applied to the electron is always perpendicular to the velocity, no work can be done on it or by it.  This means the electron's kinetic energy remains constant.  The plot above shows that after the interaction with the lens, the electron loses velocity in the x direction, meaning it now has velocity in the plane of the coil.

velocity magnitude

One way to verify our simulation is to check the velocity magnitude.  It should not change throughout the simulation.  The graph above does show a slight variance but if you look at the scale on the y axis it is incredibly small, too small to render useful numbers.  This error is well within what I would expect for general inaccuracies in a numerical simulation.  If there were problems like instability in the equation describing the field, or problem in the code I would expect a much larger effect.  This builds confidence in the results.

Electrons curving towards the axis
 As predicted, the electron will move in an arc and move towards the axis.  I did expect them to spiral afterwards though.  Once through the lens they are thrown off in different directions in essentially a straight path.
electrons converging at the focal point
 From the above graphs we get an insight into how the lensing effect works.  It's not exactly as I expected but it's close.

Focus points for different coil currents

One property of the focussing coil in a CRT that's essential is the ability for it to adjust its focus.  As the beam scans along the front of the tube, the length of the path the electrons travel changes.  The graph above shows what happens as the current varies.  A coil current of 500 Amps focusses at around 12 cm, 400 Amps focusses at around 19cm, and 300 Amps focuses at around 33cm.

This was a fun exercise that illustrates a complex problem with a few basic equations.  I plan to do another couple of posts based on this simulation using different parameters and configurations.  It may seem odd to look into how CRT's work since they are all but extinct, but there are plenty of other applications like electron microscopes that use the same technology, and with a bit of tweaking the code could easily model gravitational problems.

Sunday, April 15, 2012

Mathematics of Security Combinations

After my grandfather had a few medical issues, my family had a remote monitoring device installed in case he needed assistance.  For the system to work properly a secure lock box, that the monitoring company had the combination to, was installed at his back door to hold a set of house keys.  If ambulance officers had to gain access to the house the monitoring company could give them the combination.  As the box is no longer in use, I was having a look at it to see how secure it really is.

SUPRA StorAKeY

From the image above you can see that the unit, the SUPRA StorAKeY, consists of 10 numbered buttons, an open button, and a clear button.  The first thing to do when you get the product is set a secret PIN.  To set each number in the PIN, a set screw corresponding to the number is turned.  This has three implications, firstly the PIN can't contain duplicate numbers, secondly the PIN can contain any number of digits from 0 to 10, and thirdly the order of the numbers in the PIN doesn't matter.  So if we were to try and guess the PIN, and knew that it had 3 numbers in it, we are just randomly selecting 3 numbers out of 10, the order doesn't matter.  It's like the lottery.  Imaging that there were 10 numbered balls in a dish and you randomly selected 3 balls from it, we can calculate the probability of getting the right combination as follows.

Suppose that each time you selected a ball it was part of the correct combination.  For the first ball you have a 3 in 10 chance of a correct selection.  As there are now nine balls left and only two correct ones, the probability of selecting another correct number is 2 in 9.  It then follows that a third correct selection has a probability of 1 in 8.  By multiplying these probabilities we can calculate the overall probability of the correct combination.  (3/10)x(2/9)x(1/8) = 1/120.  So the probability of selecting the right combination is 1 in 120 which implies that there are 120 different combinations.  If we return to the problem at hand, all we need to do is calculate how many combinations there are for each PIN length and add them up to calculate the total number of combinations.

Calculating the number of combinations is easy by using the following well known statistical formula, where n is the size of the set to select from, and k is the number of selections to make.  This is the same as the process we just went through.


We can now easily calculate the number of different combinations.


So for a device that has 10 buttons, all we end up having is 1024 different combinations.  Pretty pathetic really.  If we were to attempt to brute force the PIN, and assume that each combination takes 15 seconds, it would take just over 4 hours to try every combination, but on average you will find the combination half way through the process at around 2 hours.  Therefore if you attempted this you would need 4 hours maximum but it would probably take a lot less time than that.  If you were then to take into account psychology, It would take even less.  I'd wager that most people that use these would use birthdays and years as the PIN, so by trying these first you may get lucky.

I think the thing to take away from this is that these kind of devices are a deterrent and not extremely secure.  If they are in a secluded place where someone can go unnoticed for several hours, they aren't that effective.  If however they were placed in a more prominent location where someone standing there for 4 hours would draw attention, they become more effective.

An interesting side note about combinations is that they can be used to calculate the odds of winning the lotto.  For Gold Lotto in Australia you pick 6 balls from 45 to win, by using the above formula we can calculate that the odds of winning are 1 in 8145060.  PowerBall is a little different, to win you have to select 5 balls from 45 and then select a PowerBall which is a separate draw of 1 ball from 45.  The odds of the first part of this work out to be 1 in 1221759, but when multiplied by the second part of 1 in 45 the odds rise dramatically to 1 in 54979155.  It's a neat little trick, in both draws you have to pick 6 numbers, in both draws you are picking from 45 numbers, so you may think you have a similar chance of winning, but winning PowerBall is significantly less likely.

Wednesday, April 11, 2012

Off Axis Magnetic Field of a Circular Current Loop

Recently I was curious about the magnetic field generated by a circular current loop.  I had calculated the on axis field easily back in my MRI course at university, but we never went any further than that.  There is a good reason for that, it gets ugly, quickly.  Where you can solve the on axis case in about half a page, it take about 10 pages for the off axis case.

You would not believe how hard it was to find a worked solution to this on-line.  I never really did end up finding one, this is pieced together from textbooks, websites, and on-line lecture notes.  So I thought I would put them all together, inject a few of my ideas, and try to present a clear explanation of the solution.





Here is the PDF I used during the presentation and the wxMaxim file I used to calculate the final answer.

OffAxisMagneticField.pdf
OffAxisMagneticField.wxm

It was a good exercise really. I had to learn LaTeX, which I had wanted to do for a long time, and I found out about Maxima, a computer algebra system.  I used wxMaxima to solve the final part of the problem and was quite impressed.  It's layout and red blue interface made me feel like I was using Maple.  The range and quality of free software out there for Linux really amazes me, but I wouldn't mind my own copies of MATLAB, Maple, and Altium Designer.  While we're at it, throw in COMSOL Multiphysics too.


Update 3 August 2014

I've had a question about how you would calculate the circulating current in a thin disc.  My initial train of thought would be to split the disc up into infinitesimal rings calculate the magnetic field and integrate them, like this, but I feel there are some problems with this method.




Let me be clear, I don't actually know the answer, so I'll just list some of my concerns and ideas.  It's important to note that the current distribution in the disc may not be even.  The current circulating near edge has further to travel compared to current nearer the middle.  This may mean that there is less of it to give a constant voltage drop around the disc at all distances from the centre.

Another problem with this method is that you're calculating the magnetic field of each ring in isolation from the rest.  To calculate it properly you can consider each ring of current by itself, but you need to leave the conducting material of the disc behind as it will influence the shape of the magnetic field.  There could even be a problem with that as current from different rings may interact, so considering them in isolation may not be valid.

I'm sure someone could answer this problem, but at this point I'd reach for a finite element solver like www.femm.info (free) or www.comsol.com (paid, but beautiful)


Update 3 November 2015

Someone on YouTube asked a few questions about one of the equations in the video and at first I was worried I'd gotten it wrong, but everything's fine.  I'll post the question and then my response under it.  It's a good question and got me thinking about how the parts of the equation interact.

Question

Something seems wrong. On slide 12, for Bz(r,z). The term a^2-z^2-r^2 should probably be a^2-r^2. It is not physical that the field sign inverts for a certain z away from the loop. Moreover, assuming E2(0) = 1 and E1(0) = 0 (?), substituting r=0 will not result in the classical textbook equation of field at the center of a circular loop. Also, a factor π seems missing. (or E2(0) = π?)

Answer

The equation in question is quite complex and involves the interaction of elliptic integrals with other variables and is hard to tell what it should do, so I made a plot of the field for a coil with a radius of 1.55m with 1 Amp flowing through it.  The vertical axis is the z axis and the horizontal one is the radius.  The python code is here.  The arrows in the plot below only indicate direction, not magnitude of the field.

vector field
Vector field of current loop
The question asks about the sign of the field in the z direction, and yes if you're inside the coil the component of the field in the z direction will always be positive (the blue plot at 1m radius). However, as soon as you move outside of the coil (the green plot at 2m radius) the direction of the field starts out positive, becomes negative, and swings back to positive.  So the field can switch direction if the radius of interest is greater than the coil radius.
Plot
Magnitude of the magnetic field at r = 1 and 2 m along the z axis
The questioner almost answered the next part of the question themselves.  Setting r to 0 will make k zero.  The key is to then realise that E1(0) and E2(0) are both equal to pi/2.  From this point the equation reduces to the well know formula for the field on the axis of the coil.

Equations
Reducing the off axis case to an on axis case

To improve confidence in my result I found an independent derivation of this equation in a patent from 1989 for MRI equipment.  The patent number is EP0310212 A2.

corroborating independent equation 


Update 28 November 2015

George asked in the comments about the in plane magnetic field of a current loop, i.e. z=0.  Due to symmetry Br will be equal to zero and there will only be a Bz component.  George supplied an equation that looks like it's on the right track (I haven't verified that but it has the structure that I'd expect).  My method of solving this would be to use the equation I've already calculated and set z to zero.

George also asks about when calculating the solution in Mathmatica it seems to hang.  I'm unfamiliar with Mathmatica, but I know the equation is undefined when r is equal to the radius of the current loop, so that could be an issue.

Also asked was how I calculated the values for the elliptic integrals.  When done numerically I use python with the package scipy.special.  To illustrate this process I've put together a small script that may help others get started.  It shows that in the plane of the coil the field starts out with a positive magnitude near the axis and rises as you approach the coil.  This is familiar behaviour of a solenoid, it has a reasonably uniform field near the centre of the coil.  Outside the coil the magnetic field is in the opposite direction but it tends to zero the further away from the coil you get.  Once again this makes sense.

I don't think I completely answered the question, but I hope it helps in some way.

graph
In plane magnetic field of coil with radius of 1.55 m and current of 1 Amp


Updated 4 December 2015

I've placed a worked solution for the in plane magnetic field of a circular current loop here

http://www.grant-trebbin.com/2015/12/in-plane-magnetic-field-of-current-loop.html

It starts with the Biot Savart law and progresses from there.


Update 13 April 2016

Where I've written E1(k) it should be E1(k^2) or E1(m) where m is k^2.  Same goes for E2.  I got a little mixed up with the notation.  There are a few inconsistencies in notation depending on what source you read.  The mathematics still work though.