Guide to modelling permanent magnets

For the free-electron maser design we describe in Prototype free-electron maser design, it is essential to model the magnetic fields in the undulator as well as in the collimating magnets and magnetic lenses. This nearly conducting computational magnetostatics, since the magnetic fields do not have an exact analytical solution for all but the simplest cases. We will give a brief guide for how these simulations are performed.

Introduction to magnetostatics

We will start by giving a basic overview of magnetostatics. Magnetostatics describes the behavior of charges in time-independent magnetic fields. This is only an approximate description of electromagnetism, but works well when the fields vary very gradually (or not at all) in time. The two fundamental equations of magnetostatics are two of the Maxwell equations (Gauss’s law and Ampere’s law) in the time-independent regime, and take the form:

$$\nabla \cdot \mathbf{B}=0,\nabla \times \mathbf{H}=\mathbf{J}$$

Where $\mathbf{J}$ is the total current density, $\mathbf{B}$ is the magnetic field, $\mathbf{H}$ is the demagnetizing field, and the latter two are related by the constitutive relation of magnetostatics:

$$\mathbf{B}={\mu }_0(\mathbf{M}+\mathbf{H})$$

Where $\mathbf{M}$ is the magnetization of the material. Solving these equations for different boundary conditions allows us to find the fields, and therefore the trajectory of particles within these fields. The study of boundary-value problems (BVPs) in magnetostatics, and in particular the specific BVPs that we encounter in modelling free-electron lasers/masers, are discussed within this guide.

General formulation

To start, we are concerned with finding the magnetic field $\mathbf{B}$ outside an arbitrary permanent magnet (we will discuss specific geometries later). It is also possible to find the fields in the interior of the magnet, but this is less useful for our particular applications, so we first consider the external fields only.

Outside a permanent magnet, assuming there are no significant currents present and magnet is surrounded by a non-permeable medium ($\mu \approx {\mu }_0$), then $\nabla \times \mathbf{H}=\mathbf{J}$ simplifies to $\nabla \times \mathbf{H}=0$. Then it is possible to write $\mathbf{H}$ in terms of a potential ${\psi }_m$, where ${\psi }_m$ is the magnetic scalar potential, and satisfies $\mathbf{H}=-\nabla {\psi }_m$. By substituting into Gauss’s law ($\nabla \cdot \mathbf{B}=0$) where $\mathbf{B}={\mu }_0\mathbf{H}$ in the surrounding region, we arrive at:

$$\nabla \cdot \mathbf{B}=\nabla \cdot {\mu }_0(-\nabla {\psi }_m)=-{\mu }_0{\nabla }^2{\psi }_m=0$$

Therefore, after straightforward simplifications, we arrive at Laplace’s equation for the magnetic scalar potential:

$${\nabla }^2{\psi }_m=0$$

This can be solved in the exterior case to arrive at the following analytical solution¹:

$${\psi }_m(\mathbf{r})=\frac{1}{4\pi }{\oint }_S\frac{\mathbf{M}\cdot \hat{\mathbf{n}}}{|\mathbf{r}-{\mathbf{r}}^{\prime }|}d{A}^{\prime }$$

For a uniform permanent magnet, the magnitude of the magnetic field when no external field is applied (that is, $\mathbf{H}=0$) is defined to be equal to its magnetic remanence ${\mathbf{B}}_r$; thus from the constitutive relation $\mathbf{B}={\mu }_0(\mathbf{M}+\mathbf{H})$ for magnetic fields, we have²:

$${\mathbf{B}}_r={\mu }_0\mathbf{M}\leftrightarrow \mathbf{M}=\frac{{\mathbf{B}}_r}{{\mu }_0}\text{(inside magnet)}$$

Therefore, using this definition, the integral can be evaluated explicitly, though an analytical expression can only be found for simple geometries; the integral must be evaluated numerically (for instance, via Monte-Carlo methods or numerical quadrature) for all other cases. Here, we will take the example of a rectangular bar magnet, for which the surface integral can be done analytically. This corresponds with the magnetic geometry and type used in our undulator design, which uses permanent magnets only (as shown below).

Thus, in our case, $\mathbf{M}=(0,M_y,0)$ where $M_y=B_r/{\mu }_0$, and the integral becomes:

$${\psi }_m(\mathbf{r})=\frac{B_r}{4\pi {\mu }_0}{\oint }_S\frac{d{A}^{\prime }}{|\mathbf{r}-{\mathbf{r}}^{\prime }|}=\frac{B_r}{4\pi {\mu }_0}{\oint }_S\frac{d{A}^{\prime }}{\sqrt{(x-x^{\prime }{)}^2+(y-y^{\prime }{)}^2+(z-z^{\prime }{)}^2}}$$

Once the magnetic scalar potential is computed, the magnetic field can then be found from:

$$\mathbf{B}=-{\mu }_0\nabla {\psi }_m$$

Analytical solutions

In general, the integral solution for the magnetic scalar potential can only be solved exactly in the case of very simple and highly-symmetrical geometries, such as the four geometries shown below:

Diagrams originally from supermagnete.de

The below analytical formulas correspond to the above magnets, where $B_r$ is the magnetic remanence as before³, and the $z$ axis points from the south pole of the magnet to its north pole:

$$\begin{array}{rl}{\mathbf{B}}_{block}& =\frac{{\mathbf{B}}_r}{\pi }\left[{\tan }^{-1}\left(\frac{LW}{2z\sqrt{4z^2+L^2+W^2}}\right)-{\tan }^{-1}\left(\frac{LW}{2(D+z)\sqrt{4(D+z{)}^2+L^2+W^2}}\right)\right]\\ {\mathbf{B}}_{cylinder}& =\frac{{\mathbf{B}}_r}{2}\left[\frac{D+z}{\sqrt{R^2+(D+z{)}^2}}-\frac{z}{\sqrt{R^2+z^2}}\right]\\ {\mathbf{B}}_{ring}& =\frac{{\mathbf{B}}_r}{2}\left[\frac{D+z}{\sqrt{R_a^2+(D+z{)}^2}}-\frac{z}{\sqrt{R_a^2+z^2}}-\left(\frac{D+z}{\sqrt{R_i^2+(D+z{)}^2}}-\frac{z}{\sqrt{R_i^2+z^2}}\right)\right]\\ {\mathbf{B}}_{sphere}& =\frac{2{\mathbf{B}}_r}{3}\frac{R^3}{(R+z{)}^3}\end{array}$$

There is, however, no general solution and the integral must generally be done numerically or approximated.

Derivation of exact solution for block magnet

Note: The below section uses the coordinate conventions and symbol notation given at the start of this guide. They are different from the conventions used in the analytical formulas for permanent magnets.

Let us assume that the magnets are of rectangular geometry, and have length $L$, width $w$, and depth $d$, where the length is measured along the $y$ axis, the width is measured along the $z$ axis, and the depth is measured along the $x$ axis. Thus, the rectangular cuboid we are integrating over is defined by $x\in [-d/2,d/2],y\in [0,L],z\in [0,w]$. To make this integral more tractable, we can assume that the magnets are “thin” magnets and thus set $d\to 0$, giving us the slightly “simpler” integral:

$${\psi }_m(\mathbf{r})=\frac{B_r}{4\pi {\mu }_0}{\oint }_S\frac{d{A}^{\prime }}{\sqrt{x^2+(y-y^{\prime }{)}^2+(z-z^{\prime }{)}^2}}$$

Where here we need only integrate over one surface (by symmetry, whether we choose this to be the front or back of the magnet doesn’t matter) as opposed to the six integrals that would otherwise need to be done for $d\ne 0$. Thus the surface integral reduces to the following double integral:

$${\psi }_m(\mathbf{r})=\frac{B_r}{4\pi {\mu }_0}{\int }_0^L{\int }_0^w\frac{d{z}^{\prime }d{y}^{\prime }}{\sqrt{x^2+(y-y^{\prime }{)}^2+(z-z^{\prime }{)}^2}}$$

To solve this integral, we first integrate with respect to $z^{\prime }$ with the following substitution:

$$\begin{array}{c}{\int }_0^L{\int }_0^w\frac{d{z}^{\prime }d{y}^{\prime }}{\sqrt{x^2+(y-y^{\prime }{)}^2+(z-z^{\prime }{)}^2}}={\int }_0^L{\int }_z^{z-w}\frac{du}{\sqrt{R^2+u^2}}d{y}^{\prime }\\ u=z-z^{\prime },R^2=x^2+(y-y^{\prime }{)}^2\end{array}$$

Now we may use the following identity:

$$\int \frac{dx}{\sqrt{a^2+x^2}}=\ln \left|x+\sqrt{x^2+a^2}\right|+C=\mathrm{arcsinh}\left(\frac{x}{a}\right)+C,a>0$$

Which allows us to evaluate the integral in $u$ as follows:

$$\begin{array}{rl}{\int }_0^L{\int }_z^{z-w}\frac{du}{\sqrt{R^2+u^2}}d{y}^{\prime }& ={\int }_0^L\left[\mathrm{arcsinh}\left(\frac{z-w}{R}\right)-\mathrm{arcsinh}\left(\frac{z}{R}\right)\right]d{y}^{\prime }\end{array}$$

Now, substituting in the explicit form of $R$, we have:

$$\begin{array}{rl}{\int }_0^L\left[\mathrm{arcsinh}\left(\frac{z-w}{R}\right)-\mathrm{arcsinh}\left(\frac{z}{R}\right)\right]d{y}^{\prime }& ={\int }_0^L\mathrm{arcsinh}\left(\frac{z-w}{\sqrt{x^2+(y-y^{\prime }{)}^2}}\right)d{y}^{\prime }\\ & -{\int }_0^L\mathrm{arcsinh}\left(\frac{z}{\sqrt{x^2+(y-y^{\prime }{)}^2}}\right)d{y}^{\prime }\end{array}$$

With the substitution $v=y-y^{\prime }$ we obtain the two integrals in a slightly simpler form:

$${\int }_y^{y-L}\mathrm{arcsinh}\left(\frac{z}{\sqrt{x^2+v^2}}\right)dv-{\int }_y^{y-L}\mathrm{arcsinh}\left(\frac{z-w}{\sqrt{x^2+v^2}}\right)dv$$

Both of the resulting integrals in the sum are integrals of the form:

$$\int \mathrm{arcsinh}\left(\frac{a}{\sqrt{b^2+x^2}}\right)dx$$

This is an integral that can be evaluated using integration by parts or (more conveniently) by consulting a computer algebra system; the solution is given (up to a constant) by[^7]:

$$K(a,b,x)=a\ln \left(\left|\sqrt{x^2+b^2+a^2}+x\right|\right)-b\arctan \left(\frac{ax}{b\sqrt{x^2+b^2+a^2}}\right)+x\mathrm{arsinh}\left(\frac{a}{\sqrt{x^2+b^2}}\right)$$

Where we will not write out $K$ explicitly in the rest of the solution for simplicity. Thus we have:

$$\begin{array}{rl}{\int }_y^{y-L}\mathrm{arcsinh}\left(\frac{z}{\sqrt{x^2+v^2}}\right)dv& -{\int }_y^{y-L}\mathrm{arcsinh}\left(\frac{z-w}{\sqrt{x^2+v^2}}\right)dv\\ & =K(z,x,v){|}_{v=y}^{v=(y-L)}-K(z-w,x,v){|}_{v=y}^{v=(y-L)}\\ & =\mathcal{F}(\mathbf{r},L)\end{array}$$

Where we define $\mathcal{F}(\mathbf{r},L)=\mathcal{F}(x,y,z,L)$ as:

$$\begin{array}{rl}\mathcal{F}(x,y,z,L)& =K(z,x,v){|}_{v=y}^{v=(y-L)}-K(z-w,x,v){|}_{v=y}^{v=(y-L)}\\ & =K(z,x,y-L)-K(z,x,y)\\ & -K(z-w,x,y-L)+K(z-w,x,y)\end{array}$$

And thus the complete solution can then be written succinctly as:

$${\psi }_m(\mathbf{r})=\frac{B_r\mathcal{F}(\mathbf{r},L)}{4\pi {\mu }_0}$$

While appearing deceptively-simple, this succinct answer expands out to 12 terms(!) in total, making it highly impractical to use for actually compute this solution. A much more numerically-usable solution can be found by not evaluating the integral in $y^{\prime }$, giving us:

$${\psi }_m(\mathbf{r})=\frac{B_r}{4\pi {\mu }_0}(G(\mathbf{r},L,w)-G(\mathbf{r},L,0))$$

Where the function $G=G(\mathbf{r},h,\omega )$ (here $\mathbf{r}=(x,y,z)$) is defined by the following integral:

$$G(\mathbf{r},h,\omega )={\int }_0^h\mathrm{arcsinh}\left(\frac{z-\omega }{\sqrt{x^2+(y-s{)}^2}}\right)ds$$

And thus the magnetic field can be found (via $\mathbf{B}={\mu }_0\mathbf{H}=-{\mu }_0\nabla {\psi }_m$) to be:

$$\mathbf{B}=\frac{B_r}{4\pi }(\nabla G(\mathbf{r},L,0)-\nabla G(\mathbf{r},L,w))$$

The numerical derivatives here can be found by finite differences, but the analytical derivative can also be found by differentiating under the integral sign; that is:

$$\begin{array}{rl}\nabla G(\mathbf{r},h,\omega )& ={\int }_0^h\nabla \left[\left(\mathrm{arcsinh}\left(\frac{z-\omega }{\sqrt{x^2+(y-s{)}^2}}\right)\right)\right]ds\\ & ={\int }_0^h(z-\omega ){\left(1+\frac{(z-\omega {)}^2}{x^2+(y-s{)}^2}\right)}^{-1/2}\nabla \left(\frac{1}{\sqrt{x^2+(y-s{)}^2}}\right)\end{array}$$

Note that $\mathrm{arcsinh}(\varphi )\approx \varphi$ for small $\varphi$, while $\mathrm{arcsinh}(\varphi )\approx \ln (2\varphi )$ for very large $\varphi$. If we consider the magnetic field far away from the magnet, then $1/\sqrt{x^2+(y-s{)}^2}$ becomes very small; hence we can apply the small $\varphi$ approximation, giving us:

$$\begin{array}{rl}G(\mathbf{r},h,\omega )\approx {\int }_0^h\frac{z-\omega }{\sqrt{x^2+(y-s{)}^2}}ds& =(z-\omega )\ln \left[\frac{y+\sqrt{x^2+y^2}}{(y-h)+\sqrt{x^2+(h-y{)}^2}}\right]\end{array}$$

Substituting this result into our analytical expression for ${\psi }_m$ gives us:

$$\begin{array}{rl}{\psi }_m(\mathbf{r})& =\frac{B_r}{4\pi {\mu }_0}(G(\mathbf{r},L,w)-G(\mathbf{r},L,0))\\ & =-\frac{B_r}{4\pi {\mu }_0}\ln \left[\frac{y+\sqrt{x^2+y^2}}{(y-h)+\sqrt{x^2+(h-y{)}^2}}\right]\\ & =\frac{B_r}{4\pi {\mu }_0}\ln \left[\frac{(y-h)+\sqrt{x^2+(h-y{)}^2}}{y+\sqrt{x^2+y^2}}\right]\end{array}$$

However, since we assume $y$ is large, thus $y-h\approx y$ and $(h-y{)}^2\approx y^2$. Thus the above result becomes:

$${\psi }_m(\mathbf{r})\approx \frac{B_r}{4\pi {\mu }_0}\ln \left[\frac{y+\sqrt{x^2+y^2}}{y+\sqrt{x^2+y^2}}\right]=0,\mathbf{r}\gg h,L,w$$

Hence, we see that the magnetic scalar potential (and thus magnetic field) of a permanent magnet drops off extremely rapidly to zero at large distances. This is why a fridge magnet is not dangerous, even though the magnetic remanence of its material may be quite high. A different scenario emerges, however, when we consider the magnetic field close to the surface of the magnet, that is, for small $\mathbf{r}$. In this case, we use the large $\varphi$ approximation of $\mathrm{arcsinh}(\varphi )\approx \ln (2\varphi )$, giving us:

$$G(\mathbf{r},h,\omega )\approx {\int }_0^h\ln \left(\frac{2(z-\omega )}{\sqrt{x^2+(y-s{)}^2}}\right)ds\approx {\int }_0^h\ln \left(\frac{2(z-\omega )}{s}\right)ds={\int }_0^h(\ln [2(z-\omega )]-\ln s)ds$$

This integral is straightforward to evaluate, and gives:

$$G(\mathbf{r},h,\omega )=h+h\ln (2(z-\omega ))-h\ln h$$

Therefore, the magnetic scalar potential is given by:

$$\begin{array}{rl}{\psi }_m(\mathbf{r})& =\frac{B_r}{4\pi {\mu }_0}(G(\mathbf{r},L,w)-G(\mathbf{r},L,0))\\ & =\frac{B_{r}h}{4\pi {\mu }_0}\ln \left[\frac{z-w}{z}\right]\\ & =\frac{B_{r}h}{4\pi {\mu }_0}\ln \left(1-\frac{w}{z}\right)\end{array}$$

From which we may find $\mathbf{B}={\mu }_0\mathbf{H}=-{\mu }_0\nabla {\psi }_m$ to be:

$$\mathbf{B}=-\frac{B_{r}hw}{4\pi z^2}\frac{1}{(1-w/z)}$$

It is useful to note that this contains a term that can be expanded as a geometric series:

$$\frac{1}{1-x}=1+x+x^2+x^3+\dots$$

Hence, assuming that $w/z\ll 1$ then one has:

$$\mathbf{B}\approx -\frac{B_{r}hw}{4\pi z^2}$$

Hence the magnetic field generally satisfies an inverse-square law at short distances and can be incredibly strong near the surface of a permanent magnet.

Solving with a CAS

In theory, the surface integral can also be evaluated with the following Mathematica code:

yprime = Subscript[p, y]
zprime = Subscript[p, z]
Integrate[
	 1/Sqrt[x^2 + (y - yprime)^2 + (z - zprime)^2],
	 {yprime, 0, L}, {zprime, 0, w}
 ]

Or with the following SymPy code:

from sympy import *
x, y, z = symbols("x y z", real=True)
Br, w, L = symbols("B_r w L", constant=True, 
		real=True, positive=True)
mu0 = Symbol(r"\mu_0")
yprime, zprime = symbols("y' z'")
prefactor = Br/(2*pi*mu0)
integrand = 1/sqrt(x**2 + (y-yprime)**2 + (z - zprime)**2)
integral = Integral(integrand, (yprime, 0, L), (zprime, 0, w))
psi = prefactor * integral.doit()
Eq(Symbol(r"\psi_m"), psi)

However, both SymPy and Mathematica time out when trying to compute the integral; it is likely that neither can compute it, even if they are given longer times, due to the complexity of the integral.

Far field of a permanent magnet

The dipole model is a very basic approximate model for the magnetic field of a permanent magnets. It models permanent magnets as magnetic dipoles, regardless of their geometry. This assumption is made by assuming that the separation $d$ between the magnets at the top and bottom of the undulator is large compared to the cross-sectional width of the magnets, so the magnets can be thought of as point-like.

Let us select a coordinate system where $y$ be the direction the magnets are pointing and $x$ be the direction of the particle beam. The undulator is composed of pairs of north-north and south-south facing magnets, arrayed along the $x$ axis: let us choose just a single pair to begin our analysis. Each pair is separated by distance $d$, where we let the electron beam be along $y=\frac{d}{2}$ and the two magnets be placed at $(0,0,0)$ and $(0,d,0)$ respectively.

The far-field behavior of a dipole magnet is very different. Indeed, we find that the magnetic field roughly follows an inverse-cube law:

$$|\mathbf{B}|\propto r^{-3}$$

We can derive this in a straightforward way. To start, the magnetic field of a single magnet with constant magnetization $\mathbf{m}$ could be approximated as a perfect dipole at distances “far” from the magnet, which has a magnetic field given by:

$$\mathbf{B}=\frac{{\mu }_0}{4\pi }\left[\frac{3\mathbf{r}(\mathbf{m}\cdot \mathbf{r})}{r^5}-\frac{\mathbf{m}}{r^3}\right]$$

The magnetic moments $\mathbf{m}$ of the pair of magnets are identical (both facing along either $y$ or $-y$, and of the same magnitude). We will consider a pair facing along $y$ to begin our analysis. Previously, we have already found that:

$$\mathbf{m}=\frac{1}{{\mu }_0}{B}_{r}V\hat{\mathbf{y}}$$

Where $B_r$ is the magnetic remanence (an intrinsic material property) of the magnet, and $V$ is the magnet’s volume.

Note: We adjusted the previously-derived formula $\mathbf{m}=\frac{1}{{\mu }_0}{B}_{r}V\hat{\mathbf{x}}$ since the magnetic fields and therefore magnetic moments of the magnets are aligned along $y$.

We now apply more simplifying assumptions: namely, that since the magnets are “far” from each other and from the electron beam, then the $1/r^5$ term drops off rapidly and can be effectively ignored. In addition, since the magnetic moment of the magnets are on the axis between the magnets (that is, along $y$, the direction), the magnetic field would be primarily concentrated along $y$ as well, so we can effectively ignore $B_x,B_z$. This leaves us with:

$$\begin{array}{rl}\mathbf{B}& =\frac{{\mu }_0}{4\pi }\left[\frac{3\mathbf{r}(\mathbf{m}\cdot \mathbf{r})}{r^5}-\frac{\mathbf{m}}{r^3}\right]\\ & \approx -\frac{{\mu }_0}{4\pi }\frac{\mathbf{m}}{r^3}\\ B_x=B_z\approx 0,& B_y\approx -\frac{{\mu }_0}{4\pi }\frac{\Vert \mathbf{m}\Vert }{(x^2+y^2+z^2{)}^{\frac{3}{2}}}\end{array}$$

Using this approximation, the magnetic fields of the upper magnet and lower magnets in the pair are therefore, respectively:

$$\begin{array}{rl}{\mathbf{B}}_{y,\text{ upper}}& =-\frac{{\mu }_0}{4\pi }\frac{\Vert \mathbf{m}\Vert }{(x^2+z^2{)}^{\frac{3}{2}}}\\ {\mathbf{B}}_{y,\text{ lower}}& =-\frac{{\mu }_0}{4\pi }\frac{\Vert \mathbf{m}\Vert }{{\left(x^2+z^2+{\left(y-d\right)}^2\right)}^{\frac{3}{2}}}\end{array}$$

And by the superposition principle, their combined field is given by:

$$B_y=-\frac{{\mu }_0}{4\pi }\Vert \mathbf{m}\Vert \left(\frac{1}{(x^2+z^2{)}^{\frac{3}{2}}}+\frac{1}{{\left(x^2+z^2+{\left(y-d\right)}^2\right)}^{\frac{3}{2}}}\right)$$

At $y=d/2$, which is the position of the electron beam, we therefore have:

$$\begin{array}{rl}{B}_y(x,d/2,z)& =-\frac{{\mu }_0}{4\pi }\Vert \mathbf{m}\Vert \left(\frac{1}{(x^2+z^2{)}^{\frac{3}{2}}}+\frac{1}{{\left(x^2+z^2+{\left(-\frac{d}{2}\right)}^2\right)}^{\frac{3}{2}}}\right)\\ & =-\frac{{\mu }_0}{4\pi }\Vert \mathbf{m}\Vert \left(\frac{1}{(x^2+z^2{)}^{\frac{3}{2}}}+\frac{1}{{\left(x^2+z^2+\frac{d^2}{4}\right)}^{\frac{3}{2}}}\right)\end{array}$$

We can now use the approximation that $d$ is very large compared to $x,z$, so we have:

$$\frac{1}{{\left(x^2+z^2+\frac{d^2}{4}\right)}^{\frac{3}{2}}}\approx \frac{1}{{\left(\frac{d^2}{4}\right)}^{3/2}}=\frac{8}{d^3}$$

And at $(0,d/2,z)$, that is, right between the two magnets, where the electron beam is located, we find that:

$$B_y=-\frac{{\mu }_0}{4\pi }\Vert \mathbf{m}\Vert \left(\frac{1}{z^3}+\frac{8}{d^3}\right)$$

While this expression becomes infinite (and thus unphysical) at $z=0$, showcasing the limitations of our model, we note that $1/z^3$ rapidly decays by the inverse cube of the distance, so even at a small (but nonzero) distance $z$ we have $1/z^3\approx 0$. Therefore, we may approximate the ambient field as:

$$\begin{array}{rl}{B}_y& =-\frac{{\mu }_0}{4\pi }\Vert \mathbf{m}\Vert \left(\frac{1}{z^3}+\frac{8}{d^3}\right)\\ & \approx -\frac{{\mu }_0}{4\pi }\Vert \mathbf{m}\Vert \frac{8}{d^3}\\ & \approx -\frac{{\mu }_0}{2\pi d^3}\Vert \mathbf{m}\Vert \end{array}$$

Now using $\mathbf{m}=\frac{1}{{\mu }_0}{B}_{r}V\hat{\mathbf{y}}$, we obtain the ambient field strength between the magnets:

$$\begin{array}{rl}{B}_y& =-\frac{{\mu }_0}{2\pi d^3}\left(\frac{1}{{\mu }_0}{B}_{r}V\right)\\ & =-\frac{B_{r}V}{2\pi d^3}\end{array}$$

This yields field strengths of about 1/1000th of a Tesla - not ideal at all, at least for magnets of a small volume. However, recognize that this is only the case for large $d$, meaning that the paired magnets are far apart from each other. If the two magnets are brought close to each other, the dipole approximation breaks down; however, if we make the approximation of small $d$, or at least where at least one dimension of the magnet’s geometry is roughly equal to $d$, then we have $V\approx d^3$, giving us

$$\begin{array}{rl}{B}_y& =-\frac{B_{r}V}{2\pi d^3}\\ & \approx -\frac{B_r{d}^3}{2\pi d^3}\\ & \approx -\frac{B_r}{2\pi }\\ & \approx -\frac{1}{6}{B}_r\end{array}$$

From which we may obtain fairly high field strengths of $0.15\mathrm{T}-0.25T$. Note that if one of the dimensions is greater than $\sqrt[3]{6}\approx 1.8$ times the separation distance $d$, then we have $B_y\approx -B_r$ and we can achieve very high field strengths of up to $1\mathrm{T}$ or even more.

Basic experiment using a simple compass magnet

A rudimentary way to test the correctness of this theoretical result is with a simple toy experiment. Place a small compass needle near a pair of alternating neodynium magnets. Here $\vec{m}$ is the magnetic moment of the test (compass-style) magnet, which has a magnetic remanence $B_{r,c}$ approximately that of Earth’s magnetic field ($25-65\mu \mathrm{G}$, where $1\mu G=1\times 10^{-6}\text{ T}$). Thus its magnetic moment magnitude is:

$$|\vec{m}|=\frac{1}{{\mu }_0}{B}_{r,c}{V}_{compass}$$

Where $V_{compass}$ is the compass needle volume. Meanwhile, ${\mathbf{B}}_0$ is the ambient field generated by the primary magnet. We assume that the test magnet does not disturb the field to a great degree. The magnetic field ${\mathbf{B}}_m$ of each magnet would have the following magnitude, as per our calculations:

$$B_m=\pm \frac{B_{r}V}{2\pi d^3}$$

We would expect the magnetic field to vary sinusoidally along the alternating magnets, but as a first approximation, we can assume the magnetic field to vary linearly along the $x$ axis between them. Thus, the field in between the two magnets is then:

$$B_0=\frac{B_{m^{+}}-B_{m^{-}}}{2}=\frac{B_{r}V}{2\pi d^3}$$

The magnitude of the force on the compass needle due to the ambient field is given by:

$$F_B=\nabla (\vec{m}\cdot \mathbf{B})=\frac{\partial }{\partial d}|\vec{m}|\frac{B_{r}V}{2\pi d^3}=-3|\vec{m}|\frac{B_{r}V}{2\pi d^2}=-3|\vec{m}|B_0/d$$

Assuming small oscillations, we can model this as simple harmonic oscillator, where $F=m_{mag}R{\omega }^2$, where $R$ is the length of the test magnet needle and $m_{mag}$ is its mass. Therefore, we have:

$${\omega }^2=\frac{|F_B|}{m_{mag}R}=\frac{3|\vec{m}|B_0/d}{m_{mag}R}\Rightarrow \omega =\sqrt{\frac{3|\vec{m}|B_0}{m_{mag}Rd}}$$

Note: It should be pretty evident that this can yield at most an order-of-magnitude agreement with theory. An EMF meter and a much more sensitive experiment would be necessary for anything more precise.

Interior-exterior case

We will now consider the more complicated case where we solve for both the interior and exterior magnetic fields. This involves a more complicated boundary-value problem that generally cannot be solved analytically, even in integral form. Instead, numerical methods are necessary for most cases.

To formulate the problem, we must consider the pertinent boundary conditions. First, we impose the boundary condition that the magnetic scalar potential vanish at infinity, that is:

$$\underset{\mathbf{r}\to \infty }{\lim }{\psi }_m(\mathbf{r})=0$$

In practice, this boundary condition is only useful analytically; simulating an infinitely-large domain numerically is impractical (though there are some numerical methods that can do so, such as the boundary element method), so we can simply set a cutoff radius $R$ (which is large compared to the size of the magnet), giving us:

$${\psi }_m(R)\approx 0$$

Next, we consider the boundary conditions at the interface between the interior of the magnet (with permeability $\mu$) and the exterior of the magnet (with permeability $\approx {\mu }_0$). When the relative permeability ${\mu }_r=\mu /{\mu }_0\to \infty$ (which is approximately true for many ferromagnetic materials), then just outside the surface of the magnet⁴, $H_{\parallel }$ (the tangential component) becomes zero and $H_{\perp }$ (the normal component) becomes the only component of $\mathbf{H}$; this is exactly analogous to the perfect conductor in electrostatics. Therefore, at the surface of the magnet, we obtain the following Neumann boundary condition:

$$H_{\parallel }=-\frac{\partial {\psi }_m}{\partial \tau }{|}_{\partial \Sigma }=0$$

Where $\partial \Sigma$ denotes the surface of the magnet, $\tau$ denotes the coordinate tangent to the surface, and where the solution is valid for the entirety of the exterior of the magnet. Meanwhile, we may determine $H_{\perp }$ using the standard boundary conditions for magnetostatics, given by:

$$\mathbf{K}=\mathbf{M}\times \hat{\mathbf{n}},\hat{\mathbf{n}}\times ({\mathbf{H}}^{(2)}-{\mathbf{H}}^{(1)})=\mathbf{K}$$

Where $\mathbf{K}$ is the effective surface current density, ${\mathbf{H}}^{(1)}$ is the field just underneath the surface of the magnet, and ${\mathbf{H}}^{(2)}$ is the field just outside the magnet (assumed to be vacuum or another non-permeable medium). Then by substitution, one has:

$$\begin{array}{c}\hat{\mathbf{n}}\times ({\mathbf{H}}^{(2)}-{\mathbf{H}}^{(1)})=-(\hat{\mathbf{n}}\times \mathbf{M})\\ {\mathbf{H}}^{(1)}-{\mathbf{H}}^{(2)}=\mathbf{M}\end{array}$$

More precisely, the above is only true if ${\mathbf{H}}^{(1)},{\mathbf{H}}^{(2)}$ are nearly normal to the surface (which occurs in the limit of very high relative permeabilities). In the general case, we have:

$${\mathbf{H}}_{\perp }^{(1)}-{\mathbf{H}}_{\perp }^{(2)}=\mathbf{M}$$

Where ${\mathbf{H}}_{\perp }^{(1)},{\mathbf{H}}_{\perp }^{(2)}$ are the normal components of the ${\mathbf{H}}^{(1)},{\mathbf{H}}^{(2)}$ fields, respectively. Now, to determine ${\mathbf{H}}^{(1)}$, we use the constitutive relation $\mathbf{B}={\mu }_0(\mathbf{M}+\mathbf{H})$. Here it is implied that $\mathbf{M},\mathbf{H},\mathbf{B}$ represent the interior fields and magnetizations; that is, $\mathbf{M}={\mathbf{M}}_{interior}$, $\mathbf{H}={\mathbf{H}}^{(1)}$, $\mathbf{B}={\mathbf{B}}_1$. We must now impose a constitutive relation between $\mathbf{B}$ and $\mathbf{H}$ to solve. In general, this relationship is nonlinear for permanent magnets, and thus $\mathbf{B}=f(\mathbf{H})$. However, for many permanent magnets⁵ we can make the linear approximation given by⁶:

$$\mathbf{B}={\mathbf{B}}_r+{\mu }_0{\mu }_{rec}\mathbf{H}$$

Here, ${\mu }_{rec}$ is the recoil permeability of the material and is a material constant, where we can often make the assumption ${\mu }_{rec}\approx 1$ (which is true for a perfect magnet). From these boundary conditions, we can then solve for the fields both inside and outside the magnet.

Magnetic field of an undulator

An undulator, as described in Free-electron maser physics, is a device that uses an oscillating magnetic field to stimulate the emission of synchrotron radiation. In an undulator, dipole magnets are arranged such that the magnetic field lines (shown as arrows in the below image) alternate up and down, which lead charge particles to follow sinusoidal trajectories as they pass between the magnets.

Source: European XFEL

The simplest model of the magnetic field of an undulator has been discussed already at length; it is given by:

$$B_x=B_z=0,B_y=B_0\sin k_{u}z$$

However, this is an approximation and only applies in the on-axis region. In the more general case, the magnetic field of an undulator can be derived from the magnetic scalar potential ${\psi }_m$ obtained by solving the (exterior) magnetostatic equation ${\nabla }^2{\psi }_m=0$, and the solution is given by⁴:

$${\psi }_m(\mathbf{r})=\frac{B_0}{k_u}\sinh (k_{u}y)\sin (k_{u}z)$$

Where ${\psi }_m$ is the magnetic scalar potential. Upon differentiation, the magnetic field is thus:

$$\begin{array}{rl}{B}_x& =0\\ B_y& =-B_0\cosh (k_{u}y)\sin (k_{u}z)\\ B_z& =-B_0\sinh (k_{u}y)\cos (k_{u}z)\end{array}$$

If we take a cross-section of the field in the $xy$ plane at some fixed location at $z=z^{\ast }$ along the optical axis, we have:

$$\begin{array}{rl}{B}_x& =0\\ B_y& =-B_0\cosh (k_{u}y)\sin (k_u{z}^{\ast })\end{array}$$

Alternatively, letting $A=-B_0\sin (k_u{z}^{\ast })$ be a constant, this can be also written as:

$$\begin{array}{rl}{B}_x& =0\\ B_y& =A\cosh (k_{u}y)\end{array}$$

A plot looks something like this (where the vertical direction points from the north pole of one magnet to the south pole of another):

$x$ is the horizontal coordinate, $y$ is the vertical coordinate, in accordance with standard conventions for free-electron laser physics.

Note that this is with an idealized case in which we assume that the magnetic field in the $x$ direction is zero; this assumption is only valid for small $x$. For real magnets, $B_x\ne 0$ and $B_x$ grows larger the further away you are from the magnet. Thus, we can make the ansatz $B_x\sim \sinh (k_{u}x)$, at least for $z=\text{const}$.

For $y\ne \text{const.}$ we can make the ansatz $B_x\sim \sinh x\sinh (y)$, giving us:

$$\begin{array}{rl}{B}_x& =-B_0\sinh (k_{u}x)\sinh (k_{u}y)\\ B_y& =-B_0\cosh (k_{u}y)\sin (k_{u}z)\\ B_z& =-B_0\sinh (k_{u}y)\cos (k_{u}z)\end{array}$$

A cross-section plot of the field for constant $z$ but varying $x$ and $y$ is shown below:

With such a complicated field, finding electron trajectories analytically would be extremely difficult. But it is important to note that this is still not the exact field because it does not model the dipole magnets themselves; it treats them as magnetic point dipoles. While this is acceptable as a mean-field approximation, it is inaccurate at distances very close to the magnets, where the geometry of the magnets does matter. Methods to model realistic undulators can be found in Numerical simulations for magnetostatics.

Magnetic field of a magnetic lens

Magnets are also essential for magnetic focusing: controlling and collimating charged particle beams. The type of magnets used depends on the application: both permanent magnets and electromagnets can be used. Regardless of the type, magnetic lenses are typically one of two categories:

1. Multipole magnetic lenses, which compress, stretch, and/or squeeze particle beams using a combination of multipole magnets
2. Thin solenoid lenses, which use a thin solenoid (basic electromagnet) instead of multipole magnets

Let us first consider the case of multipole magnets. Multipoles are a class of magnets built out of two, four, six, or even eight magnets, arranged in a symmetrical fashion. These are correspondingly referred to as dipoles, quadrupoles, sextupoles, and octupoles, reflecting the number of magnets typically used in each configuration, as shown below:

Diagram of multipole magnets. Note: orange arrows represent magnetic field lines.

A multipole magnet’s field can be assumed to be constant (or zero) along $z$, that is, $B_z=0$. However, $B_x,B_y$ are in general non-constant, meaning that $B_x=B_x(x,y)$ and $B_y=B_y(x,y)$. It is conventional to use complex numbers to represent the magnetic fields of multipole magnets, since they are naturally suited for describing two-dimensional vector-like quantities. Note that the use of complex numbers is purely for mathematical convenience; magnetic fields are of course real-valued. But by making use of this formalism, the magnetic field of a multipole magnet can be written very simply as a power series of an analytic function in the complex plane¹:

$$B_y+i{B}_x=\sum _{n=1}^{\infty }{C}_n(x+iy{)}^{n-1}$$

Where $C_1,C_2,\dots ,C_n$ are constant coefficients that can be real or imaginary. From the above, it follows that:

$$\begin{array}{r}{B}_x(x,y)=\mathrm{Im}\left\{\sum _{n=1}^{\infty }{C}_n(x+iy{)}^{n-1}\right\}\\ B_y(x,y)=\mathrm{Re}\left\{\sum _{n=1}^{\infty }{C}_n(x+iy{)}^{n-1}\right\}\end{array}$$

A dipole magnet (that is, a single magnet with two poles) is represented by the first term ($n=1$) of the series. To describe a quadrupole (four magnets with four poles), we include the second ($n=2$) term. To describe a sextupole magnet, we also include the third ($n=3$) term, hence the name “sextupole” (which comes after the dipole (first) term and the quadrupole (second) term). It is noteworthy that the dipole term corresponds to a constant magnetic field; the quadrupole term corresponds to a magnetic field that varies linearly with distance; and the sextupole term represents a magnetic field that varies quadratically with distance. (Continuing on this progression, the octupole term represents a field that is cubic with respect to the distance.)

Note: For the purposes of our discussion, the precise values of $C_1,C_2,C_3$ don’t matter, since we’re only interested in the qualitative characteristics of the solution, which are not affected by the values of the constants.

In general, dipole magnets are used for beam steering (that is, turning a beam left or right), sextupoles are used for correcting the effects of chromatic aberration⁶ (a form of misalignment caused by energy differences in the beam), and octupoles are used for beam stabilization⁵. However, of the different multipoles, only quadrupole magnets can be used for beam collimation: that is, focusing charged particles and keeping them closely-clustered. This cannot be done with a single quadrupole magnet alone; rather, a pair of quadrupole magnets (also called a quadrupole doublet) is required for this purpose². The two quadrupoles respectively squeeze the charged particles horizontally- and vertically-together, allowing for a formation of a collimated beam, as is evident from the following graphs (from the Feynman lectures):

Original source: Feynman lectures II, Ch. 29.7

A thin solenoid lens, by contrast, is effectively a tightly-wound coil of wire, made such that it has an extremely short length. This has a similar effect as a single loop of wire in that it causes electrons to spiral inwards in helical paths as they pass through the coil. This method is used by deflection yokes on cathode-ray tubes (CRTs) in old television sets, and is notably much simpler; it is not very difficult to simply coil a wire and run a current through it. In addition, while solenoids are generally weaker than natural magnets at low current, it is possible to use a ferromagnetic (e.g. iron) core to greatly increase the magnetic flux density. A basic comparison between the two types is given below; for a more comprehensive treatment, refer to (Hofmann, 2013).

Type of magnetic lens	Advantages	Disadvantages
Quadrupole doublet	Can use permanent magnets (no external power source), widely-used in research	Complex design, not as easy to make, difficult to control
Thin solenoid lens	Simple design, widely-used in commercial equipment, easy to control by adjusting current	Requires constant power to provide current for electromagnet

Numerical simulations

We should mention that in our prior discussion, we generally considered only cases where an analytical solution can be found for the ${\psi }_m$, $\mathbf{H}$, and $\mathbf{B}$ fields. However, this is not possible in most cases, and we can only use numerical simulations. More details on this topic can be found in Numerical simulations for magnetostatics.

Footnotes

1. See the presentation Magnetic Multipoles, Magnet Design by Krafft et. al. at Jefferson Lab ↩ ↩²

2. A sextupole magnet notably cannot be used for beam collimation; it can only twist a beam into a triangular shape. See here for more details. ↩ ↩²

3. Originally from this site (archived version) ↩

4. P. Schmüser et al., Free-Electron Lasers in the Ultraviolet and X-Ray Regime (2nd ed.), Ch. 2 pg. 11 ↩ ↩²

5. See Fermilab’s explanatory article ↩ ↩²

6. If we use two sextupole magnets of opposite polarities (that is, one flipped so all N-S pairs become S-N pairs), the two lenses can act as a combined converging-diverging lens. This means that rather than diverging the beam (as with a single sextupole magnet), the two magnets correct beam aberrations and restore the beam back to a circular shape, hence improving the quality of the beam. ↩ ↩²