Earnshaw’s theorem demonstrates that stable levitation isn’t possible using static magnetic or electric fields. His proof relies on the idea that the potential energy surface at all points in space in any static system has no local minima. In this article, I’ll discuss potential energy and its relationship to Earnshaw’s theorem and verify that the theorem is, in-fact, correct.
Disclaimer: This article was originally published on December 14, 2014 on my original blog. In the intervening time my blog was hacked. This version is a re-upload with some minor updates.
In 1839, a mathematician by the name of Samuel Earnshaw was attempting to verify the existence of the ether. He noted one particular property of electrostatically charged particles suspended in a static (unchanging) electric field. Namely, that no configuration of such fields can yield stable levitation in three dimensional space. In-fact it turned out that the same was true for static magnetic fields and magnetically “charged” particles [1].
What exactly does it mean for something to be stable? For a (electrically or magnetically) charged particle to levitate stably means that a small perturbation applied to that particle in any direction will result in a corrective force in the opposing direction. But … to understand why stability isn’t possible it’s best to look at energy. Any object or particle in space will move to minimize its potential energy. Earnshaw demonstrated that no minimum exists in the potential energy surfaces resulting from static systems.
In this article, I’ll talk about potential energy as it applies to Earnshaw’s Theorem and explain his proof. I’ll specifically be focusing on electric and magnetic potential energy – however Earnshaw’s Theorem holds for any 1/rn type force
The next few sections are broken down as follows. Richard Feynman provides excellent coverage of many basic theoretical elements used in this article in [2].
- What is Energy? The Work-Energy Theorem
- Kinetic v.s. Potential Energy
- Electromagnetic Potential Energy
- Earnshaw’s Theorem (and Proof)
- Electrostatic Potential Energy Surfaces
- Magnetostatic Potential Energy Surfaces
1. What is Energy? The Work-Energy Theorem
So what is energy? Well, the classical physics definition would be in terms of work, W. I’ve defined work in the equation below – this work is defined in the context of moving something from some point A to some point B. This equation is sometimes known as the work-energy theorem. The work energy theorem is applicable only to rigid objects.
Work Energy Theorem
W_{AB} = \int_A^B {\vec{F} \cdot \vec{dr}} \tag{1} W_{AB} = K_{AB} = K_B - K_A \tag{1A}W_{AB} = - U_{AB} = - U_B + U_A \tag{1B}Work is measured in units of newton metres or joules [Nm] = [J]. The value, F is the force in newtons [N] and r is the vector displacement in [m]. The value K is the kinetic energy of the system and the value U is the potential energy of the system. Notice that work is defined positively in terms of kinetic energy but with an added negative sign in terms of potential energy. Potential energy is defined relative to a reference point at infinity. So the force vector, F is pointing towards, in a direction opposite to the displacement vector, r.
The total amount of energy remains constant as long as the forces that cause motion and positional changes are conservative – or rather the following equation holds true for conservative forces.
Total Energy
U_A + K_A = U_B + K_B \tag{2}But what is a conservative force? Or for that matter, a non-conservative force? Conservative forces are those where the path taken by the object from point A to B does not have an effect on the the total work, WAB. For instance, gravity is a conservative force, so is the electric force defined by Coulomb. In a non-conservative force, the path taken by the object is critical and leads to a different total work. For instance, frictional force would be non-conservative – the force would change depending on the path through which the object is moved. A longer path would mean more dissipative friction and more work. We should note that an electric field generated due to a changing magnetic field would also not be conservative!
2. Kinetic v.s. Potential Energy
In a mechanical system, an object has two types of energy: kinetic and potential energy. One particular definition of kinetic energy, K, is provided below and is the result of motion or velocity of a moving mass. This definition has been derived using the Work-Energy Theoerem and Newton’s law. Note that m is the mass, in [Kg], of the moving object/particle traveling at a velocity v measured in [m/s].
Kinetic Energy
W_{AB} = K_B - K_A = \frac{1}{2} m v_{\small{B}}^2 - \frac{1}{2} m v_{\small{A}}^2 \tag{3}Potential energy can be converted to kinetic energy that can then be used to do useful work in our day to day lives. The following are examples of some different types of potential energy:
- Chemical Potential Energy (UC)
- Gravitational Potential Energy (UG)
- Electric Potential Energy (UE)
- Magnetic Potential Energy (UM)
- Elastic (Spring) Potential Energy (US)
When talking about electromagnetic levitation and Earnshaw’s Theorem, I’ll mostly be dealing with electric potential energy, UE and magnetic potential energy, UM. These are defined below in equations 4 and 5 and use our earlier definition of work. The total potential energy, U of the system is the sum of these individual components – the signs of these potential energies must be chosen correctly with regard to the reference selected (any reference can be chosen).
Note that in the equations below, q represents our levitated electric charge, m represents our levitated magnetic dipole’s moment, r represents the vector distance from our reference point, E the electric field, and B our magnetic field. Instead of r, sometimes l is commonly used for the same purpose – this may be the case later in this article when we start a conversation about potential energy and the laplacian.
Electric Potential Energy
U_{E} = - \int{\vec{F}_{\small{E}} \cdot \vec{dr}} = - q \cdot \vec{E} \cdot \vec{r} \tag{4}U_{M} = - \int{\vec{F}_{\small{M}} \cdot \vec{dr}} = - \vec{m} \cdot \vec{B} \tag{5}When dealing with electric and magnetic fields, its useful to remember that both can be generated via voltage (ε) or current (I). We can define a relationship between the voltage and current, the power, as follows: P = ε I. The power is defined in units [J/s]. This means that a relationship can be created between work (W) and power (P) as shown. Likewise, since work, potential energy and kinetic energy are all interchangeable quantities the energy, U and K also has the same relationship to power.
Work and Power
P = \epsilon \; I = \frac{dW}{dt} = \frac{dU}{dt} \tag{6}This understanding of the relationship between power and work can lead to derivations for energy densities of electric and magnetic fields. And these, in turn, can help us understand why stable levitation isn’t possible with purely static configurations of electric charges and magnetic dipoles.
3. Electromagnetic Potential Energy
When talking about electrostatics we have three distinct fields: a material independent field resulting from free charge (E), an intrinsic polarization field resulting from bound charge (P) and a combination of the two combined after taking into account material properties (D). Likewise for magnetostatics we also have three fields: a material independent field also resulting from free charge (H), an intrinsic magnetization field resulting from bound charge (M) and a combination of the two after taking into account material properties (B). We can see these relationships in the following equations:
Electric Field
\vec{D} = \epsilon_0 \cdot \vec{E} + \vec{P} \tag{7}\vec{B} = \mu_0 \cdot (\vec{H} + \vec{M}) \tag{8}The magnetization field, M, for diamagnetic materials with small magnetic susceptibilities, ꞳM grows linearly to the external applied field, H such that: M=ΧM H. Now, looking back at energy again for a moment, we see that potential energy, U, can also be written in terms of potential energy density, u.
Potential Energy
U = \int_{V}{u \cdot dV} \tag{9}We can express this energy as electromagnetic potential energy as shown in the following formulations. VERY IMPORTANT! I believe the equations below only apply to materials that experience linear responses to external magnetic fields.
Potential Energy
U = U_E + U_M = \int_{V}{\frac{1}{2} \left( \vec{B\;} \cdot \vec{H\;} + \vec{D\;} \cdot \vec{E\;} \right) \cdot dV} \tag{10}Potential Energy Density
u = \frac{1}{2} \left( \vec{B\;} \cdot \vec{H\;} + \vec{D\;} \cdot \vec{E\;} \right) \tag{11}So, let’s get to the derivation of the above equation. But … errhhmmm … I’m scared! So so SOO scared – please oh please, call my mommy! … errhmm …. Right! Ignoring that brief moment of panic … let’s do it!
The potential energy density can be divided into two separate equations, one a magnetostatic equation and the other an electrostatic equation. These are shown below. So let’s look at deriving each of these individual equations.
Electric Energy Density
u_{\small{E}} = \frac{1}{2} \vec{D\;} \cdot \vec{E\;} = \frac{1}{2} \epsilon_0 E^2 \tag{11A}Magnetic Energy Density
u_{\small{M}} = \frac{1}{2} \vec{B\;} \cdot \vec{H\;} = \frac{1}{2} \frac{B^2}{\mu_0} \tag{11B}In the electrostatic case: We can start with the definition of work and its relationship to power as we saw in equation 1A.
Let’s consider the case of a parallel plate capacitor as shown in the figure to the right. The capacitor has a charge, Q, that collects on the plates as a result of a voltage source, ε. The current, I = dq/dt, from this voltage source is responsible for the charge accumulation. The plates in our capacitor are separated by distance l.
We can now use the definition of power, P, to derive an expression for the total work, W done charging the capacitor. This work is equivalent to the energy stored in the capacitor, U.
Power
\begin{aligned}
\mathcal{P} =\;&\frac{dW}{dt} = \varepsilon I \\ \;\\
dW =\;&\varepsilon \cdot I \cdot dt \\ \; \\
dW =\;&\frac{q}{C} \cdot \frac{dq}{dt} \\ \; \\
W =\;&\frac{1}{C} \int_{0}^{Q}{q \cdot dq}
\end{aligned}Capacitor Energy
W = U_{\small{E}} = \frac{1}{2} \cdot \frac{Q^2}{C} = \frac{1}{2} C \varepsilon^{2} \tag{11C}From Gauss’ Electric Field Law and Faraday’s Law an expression for the capacitance, C can be derived. This expression defines the capacitance in terms of the geometry of the parallel plate capacitor. The value E represents the electric field.
Gauss’ Electric Field Law
\begin{aligned}
\phi_E =\;&\oint_{\small{S\;}}{\vec{D} \cdot d\vec{A}} = Q \\ \; \\
D \cdot A =\;&Q \\ \; \\
D =\;&\frac{Q}{A} \\ \; \\
E =\;&\frac{Q}{\epsilon_0 A}
\end{aligned}Faraday’s Law
\varepsilon = \oint_{\small{C\;}}{\vec{E}\cdot d\vec{l}} = - \frac{d\phi_{M}}{dt}Voltage
\varepsilon = E \cdot l = \frac{Q l}{\epsilon_0 A} \tag{11D}Capacitance
\because C = \frac{Q}{\varepsilon}\;\;\;\;\therefore C = \frac{\epsilon_0 A}{l} \tag{11E}Now we can substitute our expression for the capacitance and the capacitor energy and voltage to produce an equation for the total energy, UE, and likewise the energy per unit volume, uE = U/V. Note the small u v.s. the big U.
Capacitor Energy
\begin{aligned}
U_{\small{E}} =\;&\frac{1}{2} C \varepsilon^{2} \\ \; \\
u_{\small{E}} =\;&\frac{1}{2} \frac{C \varepsilon^{2}}{V} \\ \; \\
u_{\small{E}} =\;&\frac{1}{2} \frac{\epsilon_0 A}{l} {(E \cdot l)}^2 \frac{1}{V} \\ \; \\
u_{\small{E}} =\;&\frac{1}{2} \frac{\epsilon_0 A}{l} {(E \cdot l)}^2 \frac{1}{A l} \\ \; \\
u_{\small{E}} =\;&\frac{1}{2} \frac{\epsilon_0 A E^2 l^2 }{l A l} \\ \; \\
\end{aligned}
Electric Energy Density
u_{\small{E}} = \frac{1}{2} \epsilon_0 E^2 \tag{11A}The area and length terms cancel out in the capacitor energy equation to produce the electric energy density. The electric energy density can also be represented in vector form as:
u_{\small{E}} = \frac{1}{2} \vec{D} \cdot \vec{E} \tag{11A}This is based on the definition discussed in equation 7 earlier. The energy is a scalar quantity.
The magnetostatic energy density can be derived in a similar manner to the electrostatic energy density. Again we note the definition of power, P, and use that to derive an expression for the total work done to create a magnetic field on a solenoid (inductor) – see the figure above for the configuration of this solenoid.
The magnetic field, B, generated by the solenoid can be obtained via Ampere’s Law as shown below. The magnetic field then is obtained from Ampere’s law and combined with our definition for H = B/u0.
Ampere’s Law
\mathcal{F} = \oint_{C}{\vec{H} \cdot d\vec{l}} = \frac{\partial {\phi}_{\small{E}}}{\partial t} + I \\ \; \\
\oint_{C}{\vec{H} \cdot d\vec{l}} = N\;I \\ \; \\
H\;l = N\;IMagnetic Field
B = \frac{\mu_0\;I\;N}{l} \tag{11F}From Gauss’ Magnetic Field Law we can compute a relationship between the flux, ϕM and the magnetic field, B as shown below. Note that the factor N is the number of turns of the coil in these equations. The idea being that through a given closed loop C (not to be confused with capacitance), the current I cuts through it N times – the same number of times as the loops of wire.
Gauss’ Magnetic Field Law
{\phi}_{\small{M}} = \oint_{S}{\vec{B} \cdot d\vec{A}} = 0 \tag{11G}Magnetic Flux
{\phi}_{\small{M}} = B\;N\;A \tag{11H}From Faraday’s Law we can go one step further and derive an expression for the voltage, ε. And furthermore an expression for the inductance, L.
Faraday’s Law
\begin{aligned}
\varepsilon = &- \oint_{C}{\vec{E} \cdot d\vec{l}} = - \frac{\partial {\phi}_{\small{M}}}{\partial t} \\ \; \\
\varepsilon = &- \frac{\partial {\phi}_{\small{M}}}{\partial t} \\ \; \\
\varepsilon = &- \frac{\partial (B\;N\;A)}{\partial t} = - \frac{\partial (\frac{\mu_0\;I\;N}{l} \;N\;A)}{\partial t} \\ \; \\
\varepsilon = &- \frac{\mu_0 N^2 A}{l} \frac{\partial I}{\partial t}
\end{aligned}Voltage
\varepsilon = - L \frac{\partial I}{\partial t} \tag{11I}Self Inductance
L = \frac{\mu_0 N^2 A}{l} \tag{11J}Finally, from our expression of power and combined with the previous derivations of voltage, inductance and magnetic field, the magnetic energy density can be derived as shown.
Power
\mathcal{P} = \frac{dW}{dt} = \varepsilon I \\ \; \\
dW = L \frac{\partial I}{\partial t} I\; dt = L\;I\;dIInductor Energy
W = \frac{1}{2} L\;I^2 \tag{11K}W = U_{\small{M}} = \frac{1}{2} (\frac{\mu_0 N^2 A}{l})\;(\frac{B\;l}{\mu_0 \; N})^2 \\ \; \\
U_{\small{M}} = \frac{B^2\;l\;A}{2\;\mu_0}Magnetic Energy Density
u_{\small{M}} = \frac{U_{\small{M}}}{V} = \frac{1}{2}\;\frac{B^2}{\mu_0} \tag{11B}And thus we arrive back at our earlier equation of magnetic energy density from equation 11B. The magnetic energy density can also be represented in vector form as shown below:
u_{\small{M}} = \frac{1}{2} \vec{B} \cdot \vec{H}The energy is once again a scalar quantity. We can look at potential energy in terms of energy flow as well. Another theorem, known as Poynting’s Theorem does exactly this by looking at the rate of change of potential energy (per unit volume) du/dt in a system. There is one thing that bothers me about the magnetic energy density … if I were to use the formulation: UM = -m B even after taking into account volume effects, I wouldn’t get the factor of 1/2 present in the energy density equation. Why is that? Both these expressions should be equivalent! I have a feeling this has something to do with the change in energy or the mutual energy of a two dipole system. But who knows … sigh.
4. Earnshaw’s Theorem (and Proof)
Let’s get back to Earnshaw. Earnshaw’s Theorem states that no static configuration of electric charges or magnetic dipoles in three-dimensional space can lead to stable levitation. It turns out that this is because the potential energy surface resulting from a static configuration has a saddle shape. A simple saddle surface is defined below as a visual aid. In this visual, suppose that a particle is located at the equilibrium point in our saddle surface. Any perturbation force applied to the particle will push it out of equilibrium – it will attempt to find the next lowest energy position (blue areas in the figure).
But how do we actually prove Earnshaw’s Theorem? There are a number of methods. Some of these have been shown in [3]. The first step in many of these proofs lies in understanding the shape of the force field and potential energy surfaces resulting from electromagnetic interactions. There a few common tools used for this purpose.
Gradient, Divergence, Curl and Laplacian
Do you remember the first and second derivative tests we used to do back in our early days of secondary school? The goal of these tests was to find the local minima and maximum in a 2D graph or to understand the shape of a graph. Well, we can look at similar mathematical tools for three or more dimensions, i.e. the gradient, divergence, curl and laplacian. I’ve defined these below for the three dimensional case. Note that the values for F are either scalar or vector fields / functions depending on the absence or presence of the arrow above them. And x hat, y hat and z hat are unit vectors in the x, y and z directions. Finally, Fx, Fy and Fz are the x, y and z components of the particular function or field.
Gradient
\nabla F = \frac{\partial F}{\partial x}\hat{x} + \frac{\partial F}{\partial y}\hat{y} + \frac{\partial F}{\partial z}\hat{z} \tag{12}Divergence
\nabla \cdot \vec{F} = \frac{\partial F_x}{\partial x} + \frac{\partial F_y}{\partial y} + \frac{\partial F_z}{\partial z} \tag{13}Curl
\nabla \times \vec{F} = \left| \begin{matrix} \hat{x} & \hat{y} & \hat{z} \\ \frac{\partial}{\partial x} & \frac{\partial}{\partial y} & \frac{\partial}{\partial z} \\ F_x & F_y & F_z \end{matrix}\right| \tag{14}Laplacian
\nabla^2{F} = \nabla \cdot \nabla F = \frac{\partial^2 F}{\partial x^2} + \frac{\partial^2 F}{\partial y^2} + \frac{\partial^2 F}{\partial z^2} \tag{15}When applying these equations to some specific scalar/vector field or function, the following set of conditions help point out interesting features.
Gradient-Free Condition
\nabla F = 0 \tag{12A}Divergence-Free Condition
\nabla \cdot \vec{F} = 0 \tag{13A}Curl-Free Condition
\nabla \times \vec{F} = \vec{0} \tag{14A}Local Minima Condition
\nabla^2{F} > 0 \tag{15A}Potential Energy & Minima/Maxima
Now, let’s suppose we compute the potential energy of our levitating system (either electrostatic or magnetostatic) using equation 11,
u=\frac{1}{2} \left( \vec{B\;} \cdot \vec{H\;} + \vec{D\;} \cdot \vec{E\;} \right) \tag{11}Now, we have a potential energy surface – there’ll be a potential energy value for every point in space where the levitating charge or dipole can move or be placed with respect to some fixed charges/dipoles. Below is a figure to clarify the situation. For now, I assume unit volume for each of the charged/polarized points – this means that U is equivalent to u.
Okay. Well. What will the potential energy surface, U, look like? Let’s start by considering the gradient of our surface. Where are we likely to find a minima? The answer is: the positions where the gradient (first derivative) equals zero or the places which are gradient-free. These positions are known as critical points – see Fermat’s Theorem for details.
Unfortunately, ∇U = 0 when the position in the surface is at a local maxima or where there is neither a maxima or minima, i.e. the surface is flat or a saddle point. If a function or surface is not continuous then the gradient is undefined. In such cases, a minimum may still exist and other strategies must be used for determining critical points.
To find out the specific shape of the pontential energy surface at the critical points, a further step can be taken by calculating the second derivative, or Laplacian. If the condition: ∇2U > 0 is met, then the shape of the surface at the particular point is concave up. So we can definitively say the point is a minima. If ∇2U < 0 the shape is concave down and the point is a maxima.
Unfortunately, if ∇2U = 0 then the results are inconclusive. The point may be an inflection point or it may be an indication of a completely flat graph. Because ∇2U = ∇⋅∇ U is the divergence of the gradient, we may also say the point is divergence-free. This divergence-free case may require higher-order derivative tests. Alternatively examining the values of the laplacian near this inconclusive zone, can often lead to an understanding and degree of confidence about the shape of the graph. Like the gradient, our laplacian can also be undefined as a result of discontinuities.
To prove Earnshaw’s Theorem we can utilize Gauss’ Electric and Magnetic Field Law’s (i.e the Divergence Theorem). These laws stated in differential form below. ρ represents the electric charge density.
Gauss’ Electric Field Law
\nabla \cdot \vec{D} = \rho \tag{16}Gauss’ Magnetic Field Law
\nabla \cdot \vec{B} = 0 \tag{17}Both these equations make statements about the divergence of our electric field, D and our magnetic field B. What would happen if we could change the form of the potential energy formulations represented by U? If we could re-write U in terms of the divergence of the electric and magnetic field, then we can demonstrate that the local minima condition represented by the laplacian ∇2U > 0 doesn’t hold.
The Electrostatic Case
Let’s start by using Gauss’ law (∇⋅D=ρ), Induced EMF (ε=∮E⋅dl) and the Work-Energy theorem (UE=−∫FE⋅dl) to prove that Earnshaw was right about electrostatic levitation being impossible using only static fields. We can see this proof in the equations below. The conclusion in equation 16A clearly demonstrates that the laplacian of the potential energy surface must be less than or equal to zero for any electrostatic configuration of charges. But the laplacian ∇2UE > 0 condition must be met for stable levitation. Thus Earnshaw is proven correct as we see below!
Induced EMF
\varepsilon = - \int_{C}{\vec{E} \cdot \vec{dl}} \\ \; \\
\frac{d\varepsilon}{\vec{dl}} = - \vec{E} = - \frac{\vec{D}}{\epsilon_0} \implies \vec{D} = - \epsilon_0 \; \nabla \varepsilonGauss’ Law
\nabla \cdot \vec{D} = \rho \implies \nabla \cdot (-\epsilon_0 \; \nabla \varepsilon) = \rho \tag{16}\nabla \cdot \nabla \varepsilon = -\frac{\rho}{\epsilon_0}Work Energy Theorem
U_{\small{E}} = - \int_{C}{\vec{F}_E \cdot \vec{dl}} = - \int_{C}{q\;\vec{E} \cdot \vec{dl}} \tag{1}Potential Energy (Unit Charge)
\begin{aligned}
{\bar{U}}_{\small{E}} = &- \int_{C}{\vec{E} \cdot \vec{dl}} \\ \; \\
{\bar{U}}_{\small{E}} = &\varepsilon \implies \nabla {\bar{U}}_{\small{E}} = \nabla \varepsilon \\ \; \\
\nabla \cdot \nabla {\bar{U}}_{\small{E}} = &\nabla \cdot \nabla \varepsilon \\ \; \\
\nabla \cdot \nabla {\bar{U}}_{\small{E}} = & \nabla^2 {\bar{U}}_{\small{E}} = -\frac{\rho}{\epsilon_0} \leq 0
\end{aligned}\implies \nabla^2 {U}_{\small{E}} \leq 0 \tag{16A}The Magnetostatic Case
So stable electrostatic levitation isn’t possible. What about stable magnetostatic levitation? Well to prove Earnshaw right, we will use Gauss’ law for magnetic fields (∇⋅B = 0), induced MMF (F=∮H⋅dl) and the Work-Energy theorem (UM=−∫FM⋅dl). A proof of Earnshaw’s theorem is shown below. Notice that ∇2UM = 0 which implies that once again no local minima exists in the potential energy surface for a static configuration of magnetic dipoles. And Earnshaw’s theorem is proven correct again!
Induced MMF
\begin{aligned}
\mathcal{F} = &\int_{C}{\vec{H} \cdot \vec{dl}} \\ \; \\
\frac{d\mathcal{F}}{\vec{dl}} = &\vec{H} = \frac{\vec{B}}{\mu_0} \implies \vec{B} = \mu_0 \; \nabla \mathcal{F}
\end{aligned}Gauss’ Law
\nabla \cdot \vec{B} = 0 \implies \nabla \cdot (\mu_0 \; \nabla \mathcal{F}) = 0 \tag{17}\nabla \cdot \nabla \mathcal{F} = \nabla^2 \mathcal{F} = 0Work Energy Theorem
U_{\small{M}} = - \int_{C}{\vec{F} \cdot \vec{dl}} = - \vec{m} \cdot \vec{B} \tag{5}Potential Energy (Unit Dipole)
{\bar{U}}_{\small{M}} = - \vec{B} = - \mu_0 \; \nabla \mathcal{F}\\ \; \\
\nabla^2 {\bar{U}}_{\small{M}} = - \mu_0 \; \nabla^2 \mathcal{F} = 0 \implies \nabla^2 {U}_{\small{M}} = 0 \tag{17A}Computing Potential Energy
We’ve talked about potential energy a fair amount now. But how do we actually compute potential energy, and in-fact, the magnetic field variables, B and H. It’s instructive to touch upon how these fields can be calculated using computer programs. In-fact, we’ll use these field calculations in the next section to explore potential energy surfaces in more detail. So let’s have a look.
Computations of this type usually fall under the umbrella of Finite Element Analysis (FEA). Most FEA techniques require the following four operations: Boundary Selection, Meshing, Computation and Visualization. The ideas behind these steps is really quite simple. Electromagnetic fields are a continuous phenomenon everywhere in space. But, usually we are only interested in a specific area of interest. We define this area of interest by making a physical box around it – i.e. our boundary. Now, because fields are continuous and computers are not, we often have to discretize our space within our selected boundary, this process is known as meshing. We fit a mesh with a specific granularity to our space inside the boundary and perform computations at individual elements of that mesh only. Once complete, computational results need to be visualized in an understandable way – in our case, we will use surface and vector plots to do this. Despite our best attempts, computer systems have limitations in processing power. This means that smaller boundary regions and coarser meshes are more desirable for faster computations.
At the microscopic scale, we can compute the electric field, E, and the magnetic field, H, at any point from a fixed magnetic field source using the following formulation. Note that r represents the vector distance between the test charge and each fixed charge. Likewise, r hat is the unit direction vector for this distance. Also note that k = 1/4πε.
Electric Field
\vec{E} = k \cdot \frac{Q}{{||\vec{r}||}^2} \cdot \hat{r} \tag{18}Magnetic Field
\vec{H} = \frac{1}{r^3}\;\left(3\;(\vec{m} \cdot \hat{r}) \hat{r} - \vec{m} \right) \tag{19}With electric fields given a specific fixed charge Q, we use a test charge, q and move it through every point on our mesh and compute the resultant electric field, E at each of the mesh points. In the case where there are many fixed charges, or a charged object, we can use the superposition principle to add each Qi contribution, with i representing the ith charged object. Once the electric field is computed its easy to compute the potential energy for each mesh location using our earlier formulation as shown below.
Electric Potential Energy
U_{E} = - \int{\vec{F\;}_{\small{E}} \cdot \vec{dr\;\;}} = - q \cdot \vec{E\;\;} \cdot \vec{r\;\;} \tag{4}With magnetic fields things are a bit more complicated. This is because, magnetic charge – or rather magnetic monopoles – don’t exist. In order to obtain an expression (as seen earlier) for the magnetic field we need to model magnetic dipoles (which do exist) as two magnetic monopoles. A detailed proof of this is available and provided at [4]. Once we have the expression we can give magnetic fields a similar treatment to that of electric fields. Given a fixed magnetic dipole with magnetization, M, we use a test dipole, m and move it through every point on our mesh. We can then compute the resultant magnetic field, H at each mesh point. Again if many fixed dipoles exist, Mi, superposition can be used to compute the field at a particular mesh point. Once computed, the magnetic field can be used to obtain the potential energy using the formulation we discussed earlier and which has been shown again below.
Magnetic Potential Energy
U_{M} = - \int{\vec{F\;}_{\small{M}} \cdot \vec{dr\;\;}} = - \vec{m} \cdot \vec{B} \tag{5}From our formulation of electrostatic and magnetostatic potential energy its easy to compute the force as shown below.
Electrostatic Force
\vec{F\;}_{\small{E}} = - \nabla U_{E} \tag{20}Magnetic Force
\vec{F\;}_{\small{M}} = - \nabla U_{M} \tag{21}For those interested, I’ve also got the matlab code for a simple FEA approach. At some point I’ll upload this code. But, for now, I’ll suffice to say that the use of the meshgrid function and nested for loops is all that is required to complete the FEA. Now let’s proceed with some visualizations
5. Electrostatic Potential Energy Surfaces
We can study potential energy surfaces visually by modeling simple static systems in a language like Matlab. Let’s try and do this for electrostatic potential energy surfaces and verify, visually, that Earnshaw was right. Suppose that we have no gravitational or magnetic effects and consider the configuration of charges shown in figure 3a. For now, lets deal with a strictly 2D problem, i.e. the system is constrained so that the charges cannot move along the z direction. There are 3 charges in the figure, two being fixed (Q) and one levitated (q) in 2D space. If we were to move this levitated charge over all space, then we could determine the force (and likewise the electric field) on it at all points in space too. Assuming that both Q and q are positive charges. Then we get an electric field as shown in figure 3b.
Now recall that the electric potential energy is given by
U_{\small{E}} = q \cdot \vec{E} \cdot \vec{r}The energy of our free (levitating) charge q at every position in the 2D space can be mapped as shown in figure 3c. The potential energy surface is a saddle surface. We know our charged particle will move to minimize its potential energy. So we can expect it to move along the XY plane towards the darker regions of the graph. Alright, you’re thinking, “Well, that’s that then” … “Earnshaw’s theorem is verified!” Nope. Not at all. Because we can take this step further and introduce a ring of charge. This ring of charge produces a very different conclusion. Below is the potential energy surface of a ring of charge in the same XY plane.
“Haha!”, you say. “Earnshaw’s thereom has been violated!” After all there’s a local minima at the centre of the ring of charge … (by the way, in the above figure there are 10 charges in the ring). In a 2D problem, there is in-fact, stable levitation of a charge. Of course, Earnshaw knew this. However, he was interested in untethered levitation in 3 dimensional space. It turns out that tethering a 3D system such that one degree of freedom is constrained, is one way in which we can violate Earnshaw’s theorem and achieve stable levitation.
Okay Okay … so let’s look at the 3D case. Suppose we have a ring configuration of charges, but, this time in three dimensions (figure 6a). In this case, as in the 2D case, if we take a sample levitating charge and map the forces on that charge at all points in space we can come up with a potential energy surface.
Actually, this would be a volume. So … its a bit hard to visualize. But fear not, with a little bit of work we can find the most interesting parts of the volume and look at slice planes and use those to visualize our potential energy surfaces. So, in figures 5b, 5c and 5d, I’ve taken three planes: The XY-Plane (at Z=0), the XZ-Plane (at Y=0), and the YZ-Plane (at X=0) and mapped their respective potential energies on a surface.
Notice that in the XY-Plane we see exactly the same situation as we saw earlier in the 2D scenario. Yes you can stably levitate a charge – there’s a local minima! YAy! But wait. What about the other two planes? Alright, well actually .. both the XZ and YZ plane potential energy surfaces are similar because of symmetry. You see! Earnshaw was right after all. In 3D space, you do get a saddle surface. So, while you get stability in one plane there is instability in the other planes – no minima!
What?! You’re not convinced? … yeah you’re right. A single ring of charge in 3D is hardly adequate … Alright well, let’s go ahead and look at a solenoid. A cylinder of charge. Here’s the situation – see figure 6a-d:
We see again the stability in the XY-Plane and a saddle surface in the XZ and YZ-Planes. The spiky bits in figure 6d are a result of quantization error and can be ignored – they’re caused by the slightly non-symmetric distribution of charge on that particular plane.
Okay … so what about a disc of charge? Let’s have a look. Figures 7a-d show a uniform distribution of free charge and associated potential energy surface everywhere in space. Notice that there appears to be a ridge at the zero x and y values for all z-values in figures 7c and d. I’m not completely certain why this is. My initial thought is that it has something to do with my approximation of the uniform charge distribution – but this doesn’t quite jive.
Let’s take things one step further. What about a spherical distribution of fixed charge. Surely, any free charge placed in the centre of the sphere will levitate. Figures 8a-d demonstrate the resultant potential energy surfaces for this configuration.
In this final case, we see that inside the charge sphere a free charge would have the same potential energy everywhere. This is an unstable case. The charge has no pre-disposition towards moving in any one given direction (it experiences no force). However, a perturbation force applied to the charge will cause it to undergo some acceleration until it reaches the edge of the sphere and continues to follow the minimum energy path along the edges of the figures shown above.
6. Magnetostatic Potential Energy Surfaces
Magnetostatic systems also result in similar potential energy surfaces to those we saw in electrostatic systems. We see that no stable equilibrium position is possible with a static configuration of magnets also. Consider a system of two fixed dipoles as shown below embedded in 2D space. In this 2D space, we get a B-field as shown.
If we were to now try and map the potential energy surface as we move our levitated dipole through all points in space. We would see the potential energy surface shown on the left below. Notice that the given surface encapsulates the fact that the test (levitated) dipole that we move through space is fixed in orientation. Once again, you think to yourself: “Ahhh ha! Theorem verified!”.
And again you’re wrong. For if we configure the system as a ring of static 2D dipoles such that one half of the ring has dipoles oriented inwards and the other half, dipoles oriented outwards we get a local minima. So Earnshaw’s theorem can have a local minima and hence stable levitation in 2 dimensions.
But alas, even in the 2D case, there’s one important assumption I make. I assume that the orientation (magnetization direction) of the levitated dipole is fixed to point down (as shown in the figure) regardless of where the dipole is positioned. Look at what happens when we rotate the levitated dipole (see below). No more stable equilibrium! … notice that as the dipole is rotated we see changes in potential energy.
Now I should point out that since the natural tendancy of the levitated dipole is to align itself directionally to the ring magnetic field, in reality, we would expect positional changes to re-orient the dipole such that it is always radially aligned. So there is still a possibility that the surface in the figure above is correct. But, I’ll need to investigate further to know for sure.
Well then! Is there any 2D configuration that is stable? We see that tethering one dimension and fixing the orientation of the levitated dipole results in stability. According to suggestions in the literature, in 2D, where one dimension is already tethered, it should be possible to achieve stability – so we shouldn’t need to fix the orientation. So let’s consider another configuration. What I call the pseudo-halbach configuration. In this configuration the orientation of the distinct dipoles in the ring flips back and forth as shown.
As we see, in this case the best possible scenario is that the levitated charge sits on a flat plane with no minima or maxima and undergoes no correcting force. This would still not be the ideal stable case that we’re looking for. But this is the best case scenario, in that, regardless of orientation of the levitating dipole it will experiences the same force and has the same potential energy.
It turns out that even a halbach array doesn’t really do any better than our pseudo array version. This makes sense. As I mentioned earlier, it isn’t an unreasonable assumption to assume that the dipole will re-align itself to maintain an orientation that produces our local minima. So let’s go with this working assumption for now.
Well, with that, let’s get to our 3D potential energy surfaces and confirm that we get the same sorta results. Let’s start by considering our earlier converging-diverging ring configuration whereby half the dipoles in the ring are pointing in and the other half pointing out.
We can also see the rotational potential energy as the levitated dipole is rotated in the XY plane through azimuthal angle, θ and polar angle, ϕ.
With that, I’m going to conclude. Turns out that Earnshaw was right. There isn’t any three dimensional configuration of static magnetic or electric fields that yield stable levitation. At some point in the future I’ll discuss potential energy and diamagnetism which seems to be the exception to the rule.
References
- Earnshaw S. On the Nature of the Molecular Forces which Regulate the Constitution of the Luminiferous Ether. Camb. Phil. Soc. trans [Internet]. 1839. Available from: https://books.google.co.uk/books?id=UnwfSQAACAAJ
- Feynman R. Feynman Lectures in Physics.; 1963.
- Jones W. Earnshaw’s theorem and the stability of matter. European Journal of Physics [Internet]. 1980;1:85. Available from: https://stacks.iop.org/0143-0807/1/i=2/a=004
- The Field of a Small Magnetic Dipole. Physics Insight [internet]. 2010. Available from: https://www.physicsinsights.org/dipole_field_1.html.