Fernow RC Principles of Magnetostatics
Fernow RC Principles of Magnetostatics
Fernow RC Principles of Magnetostatics
RICHARD C. FERNOW
Formerly Brookhaven National Laboratory
Shaftesbury Road, Cambridge CB2 8EA, United Kingdom
One Liberty Plaza, 20th Floor, New York, NY 10006, USA
477 Williamstown Road, Port Melbourne, VIC 3207, Australia
314–321, 3rd Floor, Plot 3, Splendor Forum, Jasola District Centre, New Delhi – 110025, India
103 Penang Road, #05–06/07, Visioncrest Commercial, Singapore 238467
www.cambridge.org
Information on this title: www.cambridge.org/9781009291149
DOI: 10.1017/9781009291156
© Richard C. Fernow 2022
This work is in copyright. It is subject to statutory exceptions and to the provisions
of relevant licensing agreements; with the exception of the Creative Commons version the
link for which is provided below, no reproduction of any part of this work may take
place without the written permission of Cambridge University Press.
An online version of this work is published at doi.org/10.1017/9781009291156 under a
Creative Commons Open Access license CC-BY-NC-ND 4.0 which permits re-use,
distribution and reproduction in any medium for non-commercial purposes providing
appropriate credit to the original work is given. You may not distribute derivative works
without permission. To view a copy of this license, visit
https://creativecommons.org/licenses/by-nc-nd/4.0
All versions of this work may contain content reproduced under license from third parties.
Permission to reproduce this third-party content must be obtained from these third-parties directly.
When citing this work, please include a reference to the DOI 10.1017/9781009291156
First published 2017
Reissued as OA 2022
A catalogue record for this publication is available from the British Library.
ISBN 978-1-009-29114-9 Hardback
ISBN 978-1-009-29116-3 Paperback
Cambridge University Press & Assessment has no responsibility for the persistence
or accuracy of URLs for external or third-party internet websites referred to in this publication
and does not guarantee that any content on such websites is, or will remain,
accurate or appropriate.
Contents
Preface page ix
1 Basic concepts 1
1.1 Current 2
1.2 Magnetic forces 3
1.3 The Biot-Savart law 6
1.4 Divergence of the magnetic field 11
1.5 Circulation of the magnetic field 12
1.6 The Ampère law 14
1.7 Boundary conditions at a current sheet 17
1.8 Inductance 19
1.9 Energy stored in the magnetic field 20
2 Magnetic materials 23
2.1 Magnetization 23
2.2 Magnetic field intensity 26
2.3 Permeability and susceptibility 27
2.4 Types of magnetism 27
2.5 Magnetic circuits 31
2.6 Boundary conditions between regions with different μ 33
2.7 Method of images 34
3 Potential theory 39
3.1 Vector potential 39
3.2 Vector potential in two dimensions 41
3.3 Boundary conditions on A 44
3.4 Vector potential for a localized current distribution 45
3.5 Force on a localized current distribution 49
3.6 Magnetic scalar potential 50
v
vi Contents
Appendices 274
A Symbols and SI units 274
B Vector analysis 275
C Bessel functions 279
D Legendre functions 283
E Complex variable analysis 287
F Complete elliptic integrals 293
Index 297
Preface
ix
x Preface
physics. I also felt it was important to expose a wider audience to a series of very
insightful papers by the late Klaus Halbach of Lawrence Berkeley Laboratory. For
the discussion of numerical methods, I decided to concentrate on a small number of
subjects while including sufficient details to enable readers to begin writing their
own computer codes, if they so desired.
The first three chapters are mostly a survey of basic material. Chapter 1 treats the
theory of magnetic fields from conductors in free space, and Chapter 2 discusses
fields from magnetic materials. Chapter 3 introduces the vector and scalar poten-
tials. It includes general solutions to the Laplace equation and the solution of
boundary value problems.
Chapters 4–6 discuss transverse fields in two dimensions. Chapter 4 looks at
fields from line currents, current sheets, and current blocks. Field quality is
introduced in terms of multipole expansions, and the effects of approximating
ideal current distributions are described. Chapter 5 looks at transverse fields using
complex variable methods. The powerful techniques for computing the fields of
block conductors using contour integration are discussed in detail. Contour inte-
grals are also derived for the fields from image currents and magnetized bodies.
Chapter 6 looks at transverse fields that are determined by the shape of the iron.
The discussion here mainly concerns the iron surfaces used in dipole and quadru-
pole magnets.
Chapters 7–9 discuss some other field configurations. Chapter 7 looks at axial
field arrangements. This includes the fields from current loops, solenoids, and
systems of coils. The solution of the solenoid field using the sheet model is treated
in detail. Chapter 8 considers periodic magnetic channels. First the field from
a helical conductor winding is discussed. Then inverse problems are introduced,
and some of the field configurations used for magnetic wigglers and particle beam-
focusing channels are examined. Chapter 9 begins with a standard treatment of the
properties of permanent magnets. This is followed by a discussion of Halbach’s
model for rare-earth cobalt magnets and his analysis of assemblies of permanent
magnets.
In Chapter 10, we relax the strict conditions of magnetostatics and allow for the
case of slowly varying currents. This leads to brief discussions of some standard
subjects such as Faraday’s law, but also some more engineering-related topics such
as eddy currents and the skin effect. This also seemed an appropriate place to
include a brief discussion of magnetic field measurements using rotating coils.
Chapter 11 discusses numerical methods. No attempt is made here to survey the
thousands of papers devoted to numerical solutions of magnetic field problems.
Instead, three methods for solving the Poisson equation are presented with
a significant amount of detail. This chapter also includes a discussion of the
POISSON code, which is freely available and extremely useful for investigating
Preface xi
2D problems. The chapter ends by returning to the inverse problem and presenting
several examples of using optimization methods.
The appendices collect some important details about mathematical techniques
and special functions used in the book.
The level of the treatment of background magnetostatic topics in the book is
typical of those encountered by undergraduate physics majors. Some of the mate-
rial in Chapters 4–11 will likely be unfamiliar to many readers. However, an
attempt has been made to include sufficient details and references, so interested
physics and engineering majors should be able to follow the discussions.
At a number of places in the text, I have indicated the source of a mathematical
relation in footnotes using the following notation.
Magnetic phenomena have been known since antiquity when a natural ore later
called lodestone was discovered to attract bits of iron. The scientific study of
magnetism dates from around 1600, when William Gilbert summarized experi-
ments on the subject in his treatise De Magnete.[1] However, interest in the
subject greatly increased after 1820, when Hans Christian Øersted reported that
electrical currents could deflect magnetic needles, thereby establishing a con-
nection between the subjects of electricity and magnetism.[2] Almost immedi-
ately, André-Marie Ampère, Jean-Baptiste Biot, and Félix Savart performed a
series of seminal experiments that determined the forces acting between current
loops. Experimental work and theoretical developments continued throughout
the first half of the nineteenth century. A long program of experimental inves-
tigations by Michael Faraday lead him to the conception that the force between
current loops occurred through the action of an intermediary field that existed in
the space around the loops. Faraday’s field concept was developed mathemati-
cally by William Thomson (later Lord Kelvin). This work culminated in a
synthesis of knowledge about electrical and magnetic phenomena by James
Clerk Maxwell in his famous treatise of 1873. Many clarifications of
Maxwell’s ideas and studies of their implications were carried out over the
next twenty years by a small group of followers. Of particular note was the
work of Oliver Heaviside who introduced the use of vector analysis and
reworked the set of equations in Maxwell’s treatise to the four equations we
use today.[3] The resulting Maxwell equations are now accepted as the theore-
tical description underlying electromagnetic phenomena.
Magnetostatics is the study of the fields, forces, and energy associated with
steady currents and magnetic materials. In this chapter, we will review some
basic concepts underlying magnetic effects due to conductor currents in free
space.
1
2 Basic concepts
1.1 Current
Experiments have shown that there exist two kinds of electrical charge q, which are
denoted as positive and negative. A current I exists when there is a net temporal
flow of charge across some arbitrary plane in space.
dq
I¼ : (1.1)
dt
If the current is flowing through a conductor with length L and cross sectional area
A, we can write the current as
ρLA
I¼ ¼ ρv A;
L=v
where ρ is the charge density and v is the velocity of the charges. The current
density J along some direction n is a vector given by
! I
J ¼ n^ ¼ ρ!
v; (1.2)
A
where ^
n is the unit vector perpendicular to A.
If we consider a volume of space V enclosed by a surface S, the conservation of
charge requires that any change in the charge density inside V must be compensated
by a flow of current through the surface or
ð ð
∂ρ !
dV ¼ J ·^n dS:
∂t
Using the Gauss divergence theorem,1 the right-hand side can be written as
ð ð
! !
n dS ¼ r· J dV:
J ·^
Then, since V is arbitrary, we can remove the integrands from the volume integrals
on both sides of the equation and obtain the continuity equation
∂ρ !
þ r· J ¼ 0: (1.3)
∂t
In magnetostatics, we have by definition ∂ρ=∂t ¼ 0, which leads to the relation
!
r· J ¼ 0: (1.4)
1
Readers unfamiliar with vector analysis should review Appendix B.
1.2 Magnetic forces 3
Often we are interested in the current flow in a “central” region far from the ends
of a magnet. If the current and the geometry are uniform along z in this region, we
can simplify the analysis by examining problems in two dimensions. If we consider
a conductor whose thickness is small compared with the distance to the observation
point, we can approximate the conductor as a current sheet.
In addition, we frequently consider line currents or “filaments,” where we
ignore the transverse dimensions of the conductor altogether and use the equiva-
lent current
ð
!
I ¼ J ·^n dS:
The direction of the magnetic field is often represented using Faraday’s concept
of lines of induction.3[4] The lines of induction are defined to be tangent to the
magnetic field at every point in space. It follows that corresponding components of
the lines of induction and the magnetic field are always proportional to each other.
If ds is a small displacement along the line of induction, we have
dx dy dz ds
¼ ¼ ¼ :
Bx By Bz B
In two Cartesian dimensions, the lines can be plotted, for example, by integrating
the equations
2
The vector B is also known as the magnetic induction.
3
Historically, these curves have been referred to as lines of force.
4 Basic concepts
Bx ðx; yÞ
dx ¼ ds
Bðx; yÞ
By ðx; yÞ
dy ¼ ds:
Bðx; yÞ
The magnitude of the magnetic field can be represented by the density of lines in a
given region. The lines of induction do not have to form closed loops.[5, 6] In
particular, the lines become undefined at locations where B = 0.
Now consider two circuits carrying currents Ia and Ib . The force exerted by
circuit a on circuit b is found experimentally to be
! !
! μ0 ! dla R
F ab ¼ Ia Ib dlb ; (1.6)
4π ∯ R3
The double integral of the first term on the right-hand side is then
! ! ! þ þ ! ! þ þ
dl a ðdl b · R Þ ! ðdl b · R Þ ! dR
¼ dl a ¼ dl a :
∯ R3 R3 R2
The last integral vanishes because the scalar integrand is taken over a closed path.
Thus we can express the force as
!! !
! μ0 R ðdl a ·dl b Þ
F ab ¼ Ia Ib : (1.7)
4π ∯ R3
In this form, we see that Newton’s law Fab ¼ Fba is obeyed since R changes
direction for the two cases.
Returning to Equation 1.6, we rewrite the force on circuit b in a form that
explicitly depends on the current in circuit b and on an integration of the elemental
4
We will use SI units exclusively in this book. For more details, see Appendix A.
1.2 Magnetic forces 5
interactions taking place around that circuit. We collect the other factors in
Equation 1.6 into a new vector Ba, which we define as the magnetic field due to
circuit a. Then the force on the circuit can be written as
þ
! ! !
F ab ¼ Ib dl b B a : (1.8)
The force acts at right angles to the direction of Ba. Dropping the subscripts, we see
that the force on a charge q moving with velocity v can be written as
ð
! dq ! ! !
F ¼ dl B ¼ q! v B : (1.9)
dt
where A is the area of the loop. Then the torque acting on the loop can be expressed as
!
τ ¼!
! m B: (1.11)
where we have dropped the subscripts referring to circuit a. The vector R points
from the current element source to the observation (or field) point where the
magnetic field is determined. This relation, known as the Biot-Savart Law, is an
important tool for finding analytic and numerical solutions for the magnetic field
produced by known current distributions. For a surface distribution of current, the
total current in the Biot-Savart law can be generalized to give
ð! !
! μ0 K R
B ¼ dS; (1.13)
4π R3
It is important to keep in mind that the Biot-Savart law and many of the other
mathematical laws that we will subsequently develop ultimately depend on the
validity of the experimental results on magnetic forces.
We consider next several elementary applications of the Biot-Savart law that we
will need to refer to later in this book.
1.3 The Biot-Savart law 7
! μ0 ^
B ¼ ρ ϕ 2 I;
4π
where5
ð∞
dz 1
I¼ ¼ :
0 fρ2 þ z2 g3=2 ρ2
The field is directed azimuthally around the wire and falls off with distance like 1/ρ.
5
GR 2.271.5.
8 Basic concepts
and from the previous example, the field at P due to the current in wire a is
! μ Ia
B a ¼ 0 ^z :
2π ρ
If the current direction in wire b can be either parallel or antiparallel to the current in
wire a, we find that the force per unit length of the wire is
!
dF b μ
¼ 0 Ia Ib ^x : (1.16)
dy 2π ρ
The force between the wires is attractive when the currents are in the same direction
and repulsive when they are antiparallel.
where
ð∞ ð∞
1 2π
I1 ¼ dx dy ¼
∞ ∞ fx2 þ y2 þ z2o g3=2 zo
and
ð∞ ð∞
x
I2 ¼ dx dy ¼ 0:
∞ ∞ fx2 þ y2 þ z2o g3=2
The integral I2 vanishes because the integrand is an odd function and the integration
extends over an even interval. The magnetic field above the sheet is
! μ0
B ¼ Ky ^x :
2
The direction of the field is parallel to the sheet and perpendicular to the current
density. The magnitude of the field is constant and independent of the distance from
the sheet. In the general case, the field above the sheet can be written as
! μ0 !
B ¼ K ^n ; (1.17)
2
10 Basic concepts
where n is the normal to the sheet pointing to the side where B is computed. Note that
the direction of B follows the right-hand rule with respect to the current filaments in
the sheet.
The contributions of the current elements to the field at P lie in a cone surrounding P.
By symmetry, the net field must be in the z direction and
! !
ðdl R Þz ¼ a2 dϕ:
Thus we have ð 2π
μ0 I a2
Bz ¼ dϕ
4π 0 fa2 þ z2o g3=2
(1.18)
μ0 I a 2
¼ :
2fa2 þ z2o g3=2
Note that Bz is proportional to the area of the current loop and falls off at large
distances like z3
o : The field is largest at the center of the loop where the value is
μ0 I
Bz0 ¼ :
2a
1.4 Divergence of the magnetic field 11
where we use primes to indicate the use of source coordinates. The operator r is
defined in terms of field coordinates, while the distance R is a function of both
source and field coordinates. Using the vector identity B.4, we can write
! ! ! ! ! !
r·ðJ 0 R Þ ¼ R ·ðr J 0 Þ J 0 ·ðr R Þ:
!
The first term on the right-hand side vanishes because r J 0 ¼ 0. When the
second term is written in terms of a determinant, we obtain
^x ^y ^z
!
r R ¼ ∂x ∂y ∂z ¼ 0:
x x0 y y0 z z0
This vector relation is one of the fundamental properties of magnetic fields. Its
validity depends on the fact that isolated magnetic charges (monopoles) do not
appear to exist.
If we integrate Equation 1.19 over some volume of space V that is enclosed by a
surface S, we get
ð
!
r· B dV ¼ 0:
Then using the divergence theorem, we find Gauss’s Law for magnetism
ð
!
n dS ¼ 0:
B ·^ (1.20)
12 Basic concepts
using the divergence theorem. Evaluating the last integral on the surface of the
small sphere, we find that
ð
1 ! 1
r · dS ¼ 2 4πR2 ¼ 4π:
R R
We can summarize these results by writing the expression in terms of the Dirac
delta function δ.
2 1
r ¼ 4πδðRÞ: (1.23)
R
Now we can do the first integral in Equation 1.22 using Equation B.3 to give
ð ð !0 ð
!0 0 1 0 0 J 0 1 !0 0
J ·r dV ¼ r · dV r·J dV :
R R R
The first term on the right-hand side can be converted to a surface integral using the
divergence theorem. It vanishes if the surface enclosing the volume in the
!integrals
is sufficiently large. The second integral also vanishes because r0 ·J 0 ¼ 0 for
magnetostatics. Thus we are only left with the second integral in Equation 1.22,
which because of the delta function from Equation 1.23, gives
! !
r B ¼ μ0 J : (1.24)
Thus we have shown that a steady current creates a magnetic field that circulates
around the current. This is a second fundamental vector relation for magnetic
fields.6
6
We will find in Chapter 10 that this relation requires an additional term if the current varies with time.
14 Basic concepts
On the right side, the current density integrated over the surface gives the total
current I. On the left side, we can use Stokes’s theorem from Appendix B to give
ð þ
! ! !!
ðr B Þ· dS ¼ B ·dl ;
where the contour on the right-hand side extends along the perimeter of the surface
S. Thus we have the result
þ
!!
B ·dl ¼ μ0 I: (1.25)
This equation is known as the Ampère law.7 It can be most usefully applied in
highly symmetric cases where, for example, the magnitude of the field is constant
along the integration path.
Again let us consider several elementary examples of using the Ampère law to
derive results that we use later in the book.
which falls off like 1/ρ. Since this result is independent of the radius of the conductor,
it also applies to the field from a current filament, which we previously derived using
the Biot-Savart equation.
When the path of integration is inside the conductor, only part of current is
enclosed by the path and the Ampère law gives
πρ2
Bϕ 2πρ ¼ μ0 I:
πa2
7
According to O. Darrigol,[2] this equation was first given by Maxwell. In that case, we agree with him that it’s
really not appropriate to call it Ampère’s law.
1.6 The Ampère law 15
Figure 1.6 Cylindrical conductor of radius a. The two integration paths are shown
with dotted lines.
μ0 I ρ μ0 J
Bϕ ¼ ¼ ρ; (1.27)
2πa2 2
which increases linearly with ρ.
2π ρ L Bρ ¼ 0:
Since the same argument applies for a cylinder of any radius, we must have Bρ ¼ 0
everywhere for the ideal solenoid.
Now consider the Ampère law applied to the path bdgi. The contributions to the
integral vanish along bd and gi since Bρ ¼ 0. Then we have
Bz ð0ÞL Bz ðr1 ÞL ¼ μ0 n I L;
Bz ðr1 Þ ¼ Bz ð0Þ μ0 n I:
16 Basic concepts
Figure 1.7 Cross-section of an ideal solenoid. The dashed line is the axis of the
solenoid. The dots and crosses refer to the direction of the current.
If we apply the same argument over the path befi, we find that
Bz ðr2 Þ ¼ Bz ð0Þ μ0 n I:
Both of these equations for Bz outside the solenoid have the same right-hand side.
Thus the value of Bz outside the solenoid is constant, independent of radius. However,
if we consider the total flux outside the solenoid, we find
ð 2π ð ∞
ΦB ¼ z r dr d ϕ:
Bout
0 a
The total flux would be infinite if Bz outside the solenoid is any constant other than 0.
This is clearly non-physical, so we must have Bz ¼ 0 everywhere outside the
solenoid.
Since the field vanishes outside the solenoid, applying the Ampère law to the path
bdgi gives
Bz ð0Þ ¼ μ0 nI:
where we write r for the length bc. Thus the field of the ideal solenoid is constant and
along the axis of the solenoid on the inside and it vanishes outside.
1.7 Boundary conditions at a current sheet 17
B1n S B2n S ¼ 0:
Since the surface S is arbitrary, it follows that B1n ¼ B2n or that
! !
ð B2 B1 Þ·^n ¼ 0: (1.29)
B1t L þ B2t L ¼ μ0 K L:
or in general
! ! !
ð B2 B1 Þ ^n ¼ μ0 K: (1.31)
18 Basic concepts
B2t
tan θ2 ¼ :
Bn
1.8 Inductance 19
B1t
tan θ1 ¼
Bn
B2t μ0 K
¼
Bn
μ K
¼ tan θ2 0 :
Bn
We see that in crossing the current sheet, the vector B is refracted in the direction of
the field from the current sheet, i.e., toward the normal for the positive current
density K shown here.
1.8 Inductance
Consider a coil with N turns. We define the flux linkage to be the product of the
magnetic flux going through the coil multiplied by the number of turns. The flux
linkage is proportional to the current flowing through the coils. We define the
coefficient of proportionality to be the self-inductance L of the coil. Thus we have
NΦB
L¼ : (1.32)
I
so the self-inductance is
μ0 N 2 2
L¼ πR : (1.33)
d
This result ignores any end effects present in a real solenoid.
If we now consider two coils, the mutual inductance M is defined as the flux
linkage in the second coil due to the current in the first coil. Thus we have
20 Basic concepts
N2 Φ2;1
M¼ ; (1.34)
I1
where Φ2;1 is the flux in coil 2 due to the current in coil 1. In general, M can be
defined using the Neumann equation [7]
ðð ! !
μ0 dl1 · dl2
M¼ ; (1.35)
4π r12
where r12 is the distance from the current element in the first coil to the current
element in the second coil. Note that this shows that M is a constant times a geometric
factor. The symmetry of this equation between the two coils shows that M of coil 2
due to current in coil 1 is the same as M for coil 1 due to current in coil 2.
μ0 N1 I1 2
Φ2;1 ¼ πR
d1
The force between two coaxial coils can be expressed in terms of the derivative
of their mutual inductance.[8]
∂M
Fz ¼ I 1 I 2 : (1.37)
∂z
dI
IV ¼ I 2 R þ LI :
dt
8
We will reexamine this question more rigorously in Chapter 10.
1.9 Energy stored in the magnetic field 21
The power P = IV generated by the battery is distributed between the power lost to
resistive heating in R and the power in the magnetic field associated with the
inductor L. The energy in the magnetic field is then
ð
Ð Ð dI
WB ¼ P dt ¼ LI dt ¼ LI dI
dt (1.38)
1 2
¼ LI :
2
If we consider the inductor to be a long solenoid, then using Equations 1.28 and 1.33,
1 μ0 N 2 B2 d 2 1 B2
WB ¼ πR2 2 2 ¼ πR2 d
2 d μ0 N 2 μ0
B2
wB ¼ : (1.39)
2 μ0
References
[1] F. Bonaudi, Magnets in particle physics, in S. Turner (ed.), CERN Accelerator School:
Magnetic Measurements and Alignment, Montreux, 1992, p. 2–5.
[2] O. Darrigol, Electrodynamics from Ampère to Einstein, Oxford University Press,
2000, p. 142.
[3] B. Hunt, The Maxwellians, Cornell University Press, 1991, p. 70–71, 245–247.
[4] D. Halliday & R. Resnick, Physics for Students of Science and Engineering, 2nd ed.,
Wiley, 1962, p. 759–760.
[5] K.L. McDonald, Topology of steady current magnetic fields, Am. J. Phys. 22:586,
1954.
[6] M. Lieberherr, The magnetic field lines of a helical coil are not simple loops, Am. J.
Phys. 78:1117, 2010.
22 Basic concepts
[7] G. Harnwell, Principles of Electricity and Magnetism, 2nd ed., McGraw-Hill, 1949,
p. 322.
[8] P. Lorrain & D. Corson, Electromagnetic Fields and Waves, 2nd ed., Freeman, 1970,
p. 360.
[9] D. Tomboulian, Electric and Magnetic Fields, Harcourt, Brace & World, 1965,
p. 267–269.
2
Magnetic materials
We saw in the previous chapter that magnetic fields are produced in a vacuum by
currents in conductors. In this chapter, we will consider magnetic effects asso-
ciated with matter. All materials have spin and orbital motion of charges at the
atomic scale. For most materials, the random orientations of these internal
currents tend to cancel out significant magnetic effects. However, in certain
magnetic materials, such as iron or permanent magnets, these internal currents
do not cancel, and there is a net external effect that also produces or enhances
magnetic fields. In the field theory we have been discussing, the magnetic field
B must be a continuous function of position. Thus, in magnetic materials, the
macroscopic field B must be an average of the rapidly varying local fields
surrounding the atoms in the material.[1]
2.1 Magnetization
When a magnetic material is placed in an external magnetic field, magnetic dipoles
in the material set up internal fields that modify the applied field. We saw in
Equation 1.10 that a current loop has an associated magnetic moment m.
We define the magnetization vector M as the average magnetic moment per unit
volume
! X! m i Ni Ii Ai
M ¼ ¼ ^n i ; (2.1)
i
V V
where Ai is the area of loop i. The volume V must be large enough so the sum is
statistically significant, yet small enough so that we can treat the variation of M as
approximately continuous. For uniform magnetization, all the internal current
loops cancel. However, as shown in Figure 2.1, there is still a net current around
the surface of the material.
23
24 Magnetic materials
It follows that
X
Im A
i i Im
M¼ ¼ ¼ Km ;
AL L
where Im is the Amperian loop current and Km is the magnetization surface current
density. In vector terms,
! !
K m ¼ M ^n : (2.2)
Im I 00 I 0 ∂Mz
¼ ¼ Δy: (2.3)
Δz Δz ∂y
2.1 Magnetization 25
Therefore, using Equation 2.3 we find, in analogy with the Ampère law, that
þ
!!
M ·dl ¼ Im ; (2.4)
where Im is the effective number of internal amp-turns through the loop of integra-
tion. Applying Stokes’s theorem, we can rewrite Equation 2.4 as
ð ð
! ! ! !
ðr M Þ · dS ¼ J m · dS :
where J and K without a subscript refer to the conduction current densities of free
charges.
26 Magnetic materials
only depends on the true conduction current that crosses the path of integration,
despite the presence of any magnetic materials. Applying Stokes’s theorem to
Equation 2.9, we get
ð ð
! ! ! !
ðr H Þ · dS ¼ J · dS :
Since the surface of integration is arbitrary, we find that H satisfies the differential
equation
! !
r H ¼ J ; (2.10)
where J is the conduction current density.
2.4 Types of magnetism 27
Materials where the directions of B and H are parallel are called linear materials.
eB
Δω ¼ ;
2me
where e is the electron charge and me is its mass.[4] The resulting net magnetic
moment Δm is
e2 r2 B
Δm ≃ ; (2.15)
4me
where r is the radius of the atomic orbit. Note that the induced moment is opposite
to the direction of the applied magnetic field. This effect is very small and is
masked by larger effects in paramagnetic and ferromagnetic materials.
Diamagnetic materials have χ < 0 and μr < 1, independent of temperature.
In paramagnetic materials, the magnetic moments from orbital motion and spins
do not cancel, resulting in a small permanent moment. In an external magnetic
field, torques tend to align the magnetic moments with the direction of the field.
In these materials, the induced fields act to increase the magnitude of the applied
field and μr > 1. The degree of alignment is decreased by internal collisions,
vibrations and thermal agitation inside the material. The resulting magnetization
is a spin effect given by the Langevin equation
mH kT
MðHÞ ¼ Nm coth ; (2.16)
kT mH
where N is the number of atoms per unit volume, k is Boltzmann’s constant, and T is
the temperature.[5] Note that the dependence of M on the magnetic intensity H is
temperature dependent and nonlinear in general. However, in the case where
kT 1, the paramagnetic susceptibility is given by Curie’s law
mH
M N m2
χ¼ ¼ : (2.17)
H 3kT
In certain crystalline materials where one of the electron shells is not filled, it is
possible for one or more electrons to have unbalanced electron spins. In these
ferromagnetic materials, it is possible to achieve a very high degree of magnetic
alignment. Below a characteristic temperature known as the Curie temperature,
coupling is possible between neighboring atoms, which can act together in regions
known as domains. In an external magnetic field, the size of favorably oriented
domains can grow. In addition, the magnetization directions in each domain tend to
align with the external field. The dependence of B or M on H is very nonlinear in
ferromagnetic materials. For example, we show a BH curve for a low-carbon steel
alloy in Figure 2.3.[6] The flux density increases extremely rapidly for small values
2.4 Types of magnetism 29
The material exhibits hysteresis because B is not a unique function of H. The value of
B depends on the previous history of how H was varied. The dotted line shows the
increase in B as H is increased, starting from an unmagnetized sample. After the initial
excitation, BðHÞ follows the arrows around the hysteresis loop. This effect arises
because the domain boundaries don’t completely return to their previous locations
when H is reversed. The remanence or remanent field BR is the value of B when H is
returned to 0. The nonzero value for the remanence shows that the material can remain
magnetic when the external driving current is turned off. This leads to the possibility
of making permanent magnets,1 so long as the temperature remains below the Curie
temperature. The coercivity HC is defined2 as the value of H in the negative direction
that is required to get B = 0. The intrinsic coercivity HCi is the reverse field required to
remove the magnetization in a plot of M versus H.[7] The hysteresis loop is symmetric
around the origin. Heat is generated for each cycle around the hysteresis loop.3
Magnetic materials with low values of coercivity are designated as “soft.” In this
case, it is easy for the magnetization to change direction as the external current
changes, so these materials are suitable for ac operation.[8] A number of soft
magnetic materials that have high permeabilities at low values of B are listed in
Table 2.1. Also listed are the values of B corresponding to the maximum perme-
ability, the coercivity, and the saturation value of B. The electrical resistivity ρe of
the material is important for considerations of eddy current4 losses in time-varying
operations. Values for pure iron are also listed for comparison.
1
Permanent magnets are discussed in more detail in Chapter 9.
2
Historically, this quantity is known as the coercive force.
3
The energy loss in the hysteresis loop is discussed in Chapter 10.
4
Eddy currents are discussed in more detail in Chapter 10.
2.5 Magnetic circuits 31
N I ¼ Hi Li þ Hg Lg ; (2.20)
where L is the mean length through the region, and the subscripts i and g refer to the
iron and the gap. On the assumption there is no leakage flux, we have Bi Ai ¼ Bg Ag .
Substituting into Equation 2.20, we get
Bi Bg Bg Ag L i Bg
NI¼ Li þ Lg ¼ þ Lg ;
μ μ0 Ai μ μ 0
so that
Li Lg
N I ¼ B g Ag þ :
μAi μ0 Ag
2.6 Boundary conditions between regions with different μ 33
Since μ μ0 , we can neglect the reluctance in the iron for reasonable values of L and
A and find that the field in the gap is
μ0 N I
Bg ≃ : (2.21)
Lg
However, when dealing with permeable materials, Equation 1.30 for the tangential
component of B is no longer accurate. In the present case, we can use the Ampère
law for H for a path that encompasses both sides of the boundary to find that
where we recall that K is the surface current density. If there is no surface current at
the boundary, then the tangential component of H is continuous.
Consider a linear material with K = 0 and let θ be the angle between the magnetic
field and the normal to the surface, as shown in Figure 2.8. Then we have from
Equations 2.22 and 2.23
μ1 H1 cos θ1 ¼ μ2 H2 cos θ2
H1 sin θ1 ¼ H2 sin θ2 :
If the field is incident from region 1 and μ is much larger in region 2 than region 1,
then θ1 will be smaller than θ2. Thus, the field will make a larger angle with respect
to the normal in the region with larger μ.
As a special case, consider the boundary between vacuum in region 1 and
iron in region 2 when K = 0. Then, using Equation 2.23, we find
Bt2 μ2
¼ :
Bt1 μ1
It is often useful to make the approximation that μ2 in the iron is infinite. Then
the right-hand side of this equation is infinite, which demands that Bt1 ¼ 0.
In this case, the field in the vacuum region would be perpendicular to the iron
surface.
We designate case (b) to be the situation when the observation point is inside the
iron in region 2. In this case, there is no conduction current in the observation
region. We assume the resultant field is due to a virtual current I 00 in region 1.
By symmetry, the original current and the two image currents lie on a line perpen-
dicular to the boundary. We assume that the currents all flow in the same direction
and that the distances of the currents from the boundary are equal and look for
a solution for the magnitude of the currents. Consider a point P along the boundary
at a distance y from the line connecting the currents. In case (a), the boundary
conditions are
ðaÞ d
Ht ¼ ðI I 0 Þ
þ y2 Þ
2πðd2
ðaÞ μ1 y
Bn ¼ ðI þ I 0 Þ;
2πðd2 þ y2 Þ
ðbÞ d
Ht ¼ I 00
2πðd2 þ y2 Þ
ðbÞ μ2 y
Bn ¼ I 00 :
2πðd2 þ y2 Þ
Both cases must give the same solution for any position y along the boundary. Thus,
we obtain the two equations
I I 0 ¼ I 00
μ1 I þ μ1 I 0 ¼ μ2 I 00 :
Solving these two equations for the unknown magnitudes of I 0 and I 00 , we find the
solution [12]
μ2 μ1
I0 ¼ I (2.25)
μ2 þ μ1
and
2μ1
I 00 ¼ I: (2.26)
μ1 þ μ2
In the special case (a) when region 1 is vacuum and region 2 is infinitely permeable
iron, Equation 2.25 reduces to I 0 = I.
36 Magnetic materials
μ0 I μ0 I 0
By ðxÞ ¼ :
2πða xÞ 2πða þ 2d xÞ
ax
EðxÞ ¼ 1 þ :
a þ 2d x
We show the dependence of the enhancement factor on x in Figure 2.11.
The enhancement is greater than 1 in the region to the left of the filament and then
becomes less than 1 in the region between the filament and the iron boundary.
It is also possible to use the method of images to solve problems with more
complicated arrangements of planar surfaces. For example, a line current between
two parallel iron boundaries can be solved using an infinite series of image
2.7 Method of images 37
currents,[13] and a line current near the corner of two perpendicular iron planes can
be solved using three image currents.5
The method of images is also useful when considering a line current near
a circular boundary surface. Consider a line current I0 at radius a in a region with
permeability μ1 near a circular boundary of radius R of a region with permeability
μ2 , as shown in Figure 2.12.
In case (a), the magnetic field at a point P inside the aperture of the magnet (r < a)
can be written as the sum of the fields due to the conduction current I0 at r = a and
the image current I1 at r ¼ r1 inside region 2, where
5
This is discussed in Chapter 5, Section 5.13.
38 Magnetic materials
μ2 μ1
I1 ¼ I0 (2.27)
μ2 þ μ1
R2
r1 ¼ : (2.28)
a
In case (b), the magnetic field at a point Q inside region 2 (r > R) can be written as
the sum of the field from two image currents in region 1, I2 at r = 0 and I3 at r = a,
where
μ2 μ1
I2 ¼ I0 (2.29)
μ2 þ μ1
2μ1
I3 ¼ I0 : (2.30)
μ2 þ μ1
The justifications for these statements come from the solution of the corresponding
boundary value problem, which we will discuss in Chapter 4, Section 4.2.
References
[1] L. Eyges, The Classical Electromagnetic Field, Dover, 1980, p. 142–145.
[2] W. Panofsky & M. Phillips, Classical Electricity and Magnetism, 2nd ed., Addison-
Wesley, 1962, p. 143–144.
[3] J. Roche, B and H, the intensity vectors of magnetism: a new approach to resolving
a century old controversy, Am. J. Phys. 68:438–449, 2000.
[4] D. Tomboulian, Electric and Magnetic Fields, Harcourt, Brace & World, 1965,
p. 213–217.
[5] C. Kittel, Introduction to Solid State Physics, 3rd ed., Wiley, 1968, p. 432.
[6] http://magweb.us/free-bh-curves
[7] C. Rudowicz & H. Sung, Textbook treatments of the hysteresis loop for ferromagnets:
survey of misconceptions and misinterpretations, Am. J. Phys. 71:1080–1083, 2003.
[8] G. Harnwell, Principles of Electricity and Magnetism, 2nd ed., McGraw-Hill, 1949,
p. 403–406.
[9] Magnetically soft materials, American Society for Metals (ASM) Handbook, vol. 2,
1990, p. 761–781.
[10] P. Lorrain & D. Corson, Electromagnetic Fields and Waves, 2nd ed., Freeman, 1970,
p. 405.
[11] P. Hammond, Electric and magnetic images, Inst. Elec. Eng. Monograph 379:306,
1960.
[12] P. Silvester, Modern Electromagnetic Fields, Prentice-Hall, 1968, p. 178–179.
[13] K. Binns & P. Lawrenson, Analysis and Computation of Electric and Magnetic Field
Problems, Pergamon Press, 2nd ed., 1973, p. 38–45.
3
Potential theory
39
40 Potential theory
þ ð ð
!! ! ! ! !
A ·dl ¼ ðr A Þ· dS ¼ B · dS :
Thus the magnetic flux through some surface S is given by the contour integral of
A around the perimeter of S.
þ
!!
ΦB ¼ A ·dl : (3.4)
Any solution for the vector potential is unique, provided that the sources are
confined to a finite region of space.[3]
There is an uncertainty in the defining relation for A in Equation 3.2. Any vector
whose curl vanishes can be added to A without affecting the value of B. For
example, the gradient of a scalar function ψ can be added because
r rψðrÞ ¼ 0:
This freedom can be used to fix the divergence of A. Starting with Equation 3.3, the
divergence of the vector potential can be written as
! μ
r· A ¼ 0 I;
4π
where
ð !
! 1
I ¼ J ðr0 Þr· dV 0 :
j!
r !
0
r j
Primes are used to indicate source coordinates. Note that the r operator is defined
in terms of the observation point (field) coordinates, while the integration and
current density depend on source coordinates. Replacing the r operator with one
defined in terms of source coordinates, we find
ð !
! 0 0 1
I ¼ J ðr Þr · ! !0 dV 0 :
jr r j
We choose to evaluate this integral over a large radius sphere where all the current
sources vanish. Then I ¼ 0 and we conclude that
!
r· A ¼ 0: (3.5)
This result is equivalent to setting the arbitrary scalar function ψ ¼ 0 and is known
as the Coulomb gauge. Equation 3.5 represents a constraint on the components of
the vector A.
Using Equation 1.24, we can relate the vector potential to the conduction current
density J.
! !
r ðr A Þ ¼ μ0 J ; (3.6)
Using Equation 3.5, we can eliminate the gradient term and find that
! !
r2 A ¼ μ0 J : (3.7)
This is the vector Poisson equation for A. The solution is valid both inside and
outside of the conductor. Equation 3.3 is a particular solution of this equation.
When J = 0, Equation 3.7 is called the Laplace equation. Solutions of the Laplace
equation are known as harmonic functions. In Cartesian coordinates, the quantity
!
r2 A has the three scalar components r2 Aα , where α corresponds to x, y, and z. For
example,
r2 Ax ¼ μ0 Jx :
!
In non-Cartesian coordinate systems, the components of r2 A must be found from
Equation B.7.
region of a long magnet, far from its ends. The great power of the two-dimensional
approximation will be most apparent in Chapter 5, where we develop the theory
using complex analysis.
Let us consider an infinitely long current filament along the z axis. Direct
integration using Equation 3.3 to ±∞ diverges. It is possible to perform the
integration over a long, but finite length and then to develop an expression for Az
in a power series in ρ/L.[4] However, we will adopt another approach using the
differential Equation 3.2. The field from the filament was given in Equation 1.15
and so we have
! μI ^
r A ¼ 0 ϕ:
2π ρ
∂Aρ ∂Az μI
¼ 0 :
∂z ∂ρ 2π ρ
The vector potential is constant in the z direction, so the derivative with respect to
z vanishes and we have
μ0 I dρ
dAz ¼ :
2π ρ
The presence of any constant terms in Az is not important since they will be
removed when taking derivatives to find B.
3.2 Vector potential in two dimensions 43
μ0 I
Az ¼ lnðRp Þ lnðRpo Þ lnðRn Þ þ lnðRno Þ ;
2π
where R is the distance from the current element to the observation point P. If both
filaments return at the same reference position, the dependence on R0 drops out and
we find
μ0 I Rp
Az ðx; yÞ ¼ ln : (3.9)
2π Rn
where Kz ¼ dI=ds is the sheet current density. The vector potential for a block
conductor with finite cross-section is
ðð
μ0 0 0 R
Az ðx; yÞ ¼ Jz ðx ; y Þln dx0 dy0 : (3.11)
2π Ro
44 Potential theory
∂2 Az ∂2 Az
r2 Az ¼ þ 2 ¼ 0: (3.13)
∂x2 ∂y
we see that the line integral of A around the perimeter of S must also be conserved.
Therefore the tangential component of A must be conserved on the boundary.
!ð1Þ !ð2Þ
At ¼ At : (3.14)
Azð1Þ ¼ Að2Þ
z ;
3.4 Vector potential for a localized current distribution 45
the tangential component of the curl is along x. Thus Equation 3.15 gives
1 1
ð2Þ
∂y Að2Þ ð1Þ
z ð1Þ ∂y Az ¼ Kz :
μ μ
In this case, we have two constraints that must be satisfied at the boundary.
The first integral vanishes because the current distribution consists of closed
loops. In the second integral, J corresponds to one of the components of A,
46 Potential theory
which we specify with the index i, while r0 is part of the scalar product with r,
which we specify with the index j. Thus the integral has the form
ð
I ¼ r0j Ji0 dV 0 :
The first term vanishes because of Equation 1.4 and in the second term, we have
r0 r0i ¼ ^i. Thus
!
Ji0 ¼ r0 ·ðr0i J 0 Þ
and so we have
ð
!
I ¼ r0 ·ðr0i J 0 Þr0j dV 0 :
This gives
ð ð
!0 0 !
I¼ r0i r0j J dS r0i J 0 ·r0 r0j dV 0 :
The first term vanishes for a surface outside the charge distribution. In the second
term, the gradient in the integrand vanishes except for the j component. Thus
3.4 Vector potential for a localized current distribution 47
ð
I ¼ r0i Jj0 dV 0 :
Substituting this back into Equation 3.16, we find that the vector potential is
ð
!! 1 μ0 ! r !
!
0
Að r Þ≃ 3
r J 0 dV 0 :
2 4π r
We define the magnetic moment of the current distribution as
ð
! ! 0
m ¼½ ! r J ð!
0
r ÞdV 0 : (3.19)
Then the vector potential of the current distribution can be written as [5]
!! μ !m !r
Að r Þ¼ 0 3
: (3.20)
4π r
Thus we find that the elementary form of magnetic matter is a magnetic dipole.
In spherical coordinates, the vector potential for a magnetic dipole is directed in
the azimuthal direction. We can find the magnetic field from the dipole by taking
the curl of A.
1
Br ¼ ∂θ ðAϕ sin θÞ
r sin θ
(3.21)
μ 2m
¼ 0 3 cos θ
4π r
and
1
Bθ ¼ ∂r ðrAϕ Þ
r (3.22)
μ m
¼ 0 3 sin θ:
4π r
48 Potential theory
The remaining component Bϕ ¼ 0. We see that the field of a magnetic dipole falls
off like 1/r3.
We can relate this definition of the magnetic moment, Equation 3.19, with our
discussion in Chapter 1 of the magnetic moment of a planar current loop. We let
! 0 !
J dV →I dl0 :
Then Equation 3.19 gives
þ
! ! 0 !
m ¼½I r dl0 :
!
The magnitude of the quantity ½!
0
r dl0 is the area of a triangular region inside
the current loop. The closed integral then gives the total area A enclosed by the
loop. Thus
!
m ¼ I A ^n ;
Bi ð!
r Þ ≃ Bi ð0Þ þ !
r ·rBi ð0Þ þ :
Then
ð ð h
! ! ! !0 ! 0 ! i
F ≃ B ð0Þ J ð r ÞdV þ J ð!
r Þ ð!
0 0
r ·rÞ B ð0Þ dV 0 :
The first integral vanishes since J consists of closed loops. Using Equation B.2, we
have
0 ! ! ! ! !
rð!
r ·BÞ ¼ !
r ðr B Þ þ B ðr !
r Þ þ ð B ·rÞ !
r þ ð!
0 0 0 0
r ·rÞ B :
The first term on the right-hand side vanishes outside the current distribution, while
the second and third terms vanish because the r operator refers to field coordi-
nates. Thus the force can be written as
ð
! ! 0 !
F ¼ J 0 rð! r · B ÞdV 0 :
The first term on the right side vanishes because r operates on field coordi-
nates, so
50 Potential theory
ð
! ! 0!
F ¼ r J 0 ð!
r · B ÞdV 0 :
Thus
ðh
! ! !i
F ¼ r ½ B ð!
0
r J 0 Þ dV 0 :
Using Equation 3.19, we can write this in terms of the magnetic moment as
! !
F ¼ rðB !
0
m Þ:
Finally, using Equation B.9
! ! ! ! !
m Þ!
F ¼ B ðr·! m ðr· B Þ þ ð!
m ·rÞ B ð B ·rÞ!
0 0 0 0
m:
The first and fourth terms vanish because m is independent of the field coordinates.
Thus dropping the primes, we find that the force on a magnetic dipole in an
inhomogeneous magnetic field is given by
! !
F ¼ ð!
m ·rÞ B : (3.24)
For a continuous magnetization distribution, the force is given by [6]
ð
! ! !
F ¼ ðM ·rÞ B dV: (3.25)
From the curl equation, we know that the magnetic field can be expressed as the
gradient of a scalar function
!
H ¼ rVm ; (3.26)
3.6 Magnetic scalar potential 51
Thus the product of the permeability and the normal derivative of Vm is continuous
across the boundary.
Consider a field point P nearby a current loop, as shown in Figure 3.4. Assign
a normal vector n to the loop according to the right hand rule. Now displace the
loop by the vector du. The cross product of dl with du gives the area of the shaded
area in the figure. In this case, the solid angle at P subtended by the differential area
dS changes by an amount
! ! ! ! !
d S ·^r ð du dl Þ·^r ðdl ^r Þ· du
dΩ ¼ 2 ¼ ¼ :
r r2 r2
The last equation makes use of the fact that the terms in the triple vector product
permute. To get the total change in solid angle, we sum over all the parts of the loop.
þ ! !
! dl r
dΩ ¼ du ·
r3
!
¼ du ·rΩ
Comparing these two equations, we find that the gradient of the solid angle is
given by
52 Potential theory
þ! !
dl r
rΩ ¼ : (3.30)
r3
Expressing B using the Biot-Savart law and using Equation 3.30, we find
þ!
I dl !r
rVm ¼
4π r3
I
¼ rΩ:
4π
The gradients on the two sides of the equation are proportional to each other, so Vm
must depend linearly on Ω.
I
Vm ¼ Ω: (3.31)
4π
This indicates that the magnetic scalar potential is directly related to the solid angle
that a current loop subtends at a field observation point. We can ignore any constant
term since it drops out when calculating the field. Note that Equation 3.31 implies
that the scalar potential is not a single-valued function. If P moves from just in front
of the loop to just behind it, the solid angle changes from 2π to –2π since ^n ·^r
changes sign.
Using the divergence theorem for the first term and dropping the primes gives
"ð ! ð ! #
1 M ·^n r·M
Vm ¼ dS dV : (3.33)
4π R R
We can identify the first term on the right-hand side as coming from fictitious
magnetic charges on the surface of the magnetized body with the surface current
density[8]
!
Km ¼ M ·^n (3.34)
and the second term can be described as due to the volume charge density
!
ρm ¼ r·M : (3.35)
This shows that for the purpose of finding the magnetic field outside a magnetized
body, we can replace the body with equivalent magnetic surface and volume
charges.
∂2 F ∂2 F ∂2 F
þ þ ¼ 0: (3.36)
∂x2 ∂y2 ∂z2
1 d2 X 1 d2 Y 1 d2 Z
þ þ ¼ 0:
X dx2 Y dy2 Z dz2
This equation can only be valid for all values of x, y, or z if each of the three terms is
equal to a constant. Thus we have
1 d2 X
¼ a2
X dx2
1 d2 Y
¼ b2 (3.37)
Y dy2
1 d2 Z
¼ c2 ;
Z dz2
where the three constants have to satisfy the constraint
a2 þ b2 þ c2 ¼ 0: (3.38)
In order to satisfy this equation, the set of constants {a, b, c} must contain both real
and imaginary members. The choice of which constants are real and which are
imaginary depends on the specific conditions for the problem under consideration.
As an example, let us suppose that the constants a and b are imaginary and c is
real. Then we can define a new set of real constants {α, β, γ} such that
a2 ¼ ðiαÞ2 ¼ α2
b2 ¼ ðiβÞ2 ¼ β2
c2 ¼ γ2 :
3.8 General solutions to the Laplace equation 55
In the general case, there will be a set of constants that satisfy Equation 3.37.
αn ¼ fα1 ; α2 ; . . .g
βm ¼ fβ1 ; β2 ; . . .g
1 d 2 Xn
¼ α2n
Xn dx2
1 d2 Ym
¼ β2m
Ym dy2
1 d2 Znm
¼ γ2nm :
Znm dz2
The oscillatory terms could also be written in terms of sines and cosines and the
nonoscillatory terms could be written using hyperbolic sines and cosines. The last
four terms are also a solution of the Laplace equation that allow for continuity of
the potential and the presence of external fields.
Fðρ; ϕ; zÞ ¼ RðρÞΦðϕÞZðzÞ:
d2 Φ
þ n2 Φ ¼ 0 (3.41)
dϕ2
d2 Z
k2 Z ¼ 0; (3.42)
dz2
where k and n are constants that may be real or imaginary.
In order for the azimuthal dependence to be single-valued, the parameter n in
Equation 3.41 must be a real integer. Solutions for the function Φ have the general form
Φn ðϕÞ ¼ Cn cos ðnϕÞ þ Dn sin ðnϕÞ (3.43)
when n ≠ 0 and
Φ0 ðϕÞ ¼ C0 ϕ þ D0 (3.44)
when n = 0.
The equations for R and Z have different solutions, depending on the values for
n and k.
(1) real k ≠ 0
The general form of the z dependence in Equation 3.42 is
(2) imaginary k ≠ 0
If k = i κ where κ is real, then the z solution is oscillatory.
(3) k = 0
The z solution has the form
Z0 ðzÞ ¼ E0 z þ F0 :
The general form of the solution of the Laplace equation in cylindrical coordi-
nates can then be written in the form
X
Fðρ; ϕ; zÞ ¼ C R ðkρÞΦn ðϕÞZk ðzÞ:
k;n kn n
The radial and angular parts of this equation can be separated first in the form
Fðr; θ; ϕÞ ¼ RðrÞYðθ; ϕÞ:
1 1
∂θ ðsin θ ∂θ YÞ þ ∂2ϕ Y þ nðn þ 1Þ Y ¼ 0: (3.50)
sin θ sin 2 θ
58 Potential theory
The solutions for this equation are known as spherical harmonics.[15] The constant
n must be an integer to avoid singularities in Yðθ; ϕÞ at θ = 0 and θ = π. The two
angle coordinates can in turn be separated as
Yðθ; ϕÞ ¼ ΘðθÞΦðϕÞ:
and
d2 Φ
þ m2 Φ ¼ 0; (3.52)
dϕ2
when m ≠ 0 and
Φ0 ðϕÞ ¼ C0 ϕ þ D0
where Pmn and Qn are associated Legendre functions of the first and second
m
The general form of the solution of the Laplace equation in spherical coordinates
can then be written in the form
X
Fðr; θ; ϕÞ ¼ C R ðrÞΘm
n;m nm n n ðθÞΦm ðϕÞ:
8
<0 s=2 ≤ x ≤ w=2
JðxÞ ¼ J0 w=2 ≤ x ≤ w=2
:
0 w=2 ≤ x ≤ s=2:
We replace the parallel side walls of the slot with an infinite set of image conductors,
whose current density is shown in Figure 3.6. The Fourier series representing the
current distribution is
w 2J0 X
∞
1 nkw
JðxÞ ¼ J0 þ sin cosðnkxÞ; (3.53)
s π n¼1 n 2
where k ¼ 2π=s.
Divide the problem space vertically into three regions, as shown in Figure 3.5.
In region 1, there are no currents so the vector potential A1 satisfies the Laplace
equation. To satisfy the boundary conditions, we know that the x dependence has to
correspond with the x dependence of J(x). Thus we have
X
∞
A1 ¼ ðCn enky þ Dn enky Þ cosðnkxÞ;
n¼1
where Cn and Dn are unknown coefficients. The solution A3 for region 3 also has to
satisfy the Laplace equation. Since region 3 extends to infinite values of y, the term
proportional to enky must vanish. Far from the conductor, the field must be uniform
along the x direction, so the potential must contain a term proportional to y. It must
also contain a constant term to guarantee continuity of A. Thus the general form of the
potential in region 3 is
X
∞
A3 ¼ E0 þ E1 y þ Fn enky cosðnkxÞ:
n¼1
Since region 2 contains the conductor, the vector potential A2 has to satisfy the
Poisson equation. The total potential A2 has a general (or homogeneous) part plus
a particular solution to the Poisson equation.
A2 ¼ A2h þ A2p :
3.9 Boundary value problems 61
To match the potential at the boundary with region 3, the general part of the potential
has to include a constant term and a term linear in y. Thus we have
X
∞
A2h ¼ G0 þ G1 y þ ðHn enky þ Mn enky Þ cosðnkxÞ:
n¼1
Since the current density J has a constant term and a periodic term, we look for
a particular potential of the form
X
∞
A2p ¼ G2p y2 þ Ln cosðnkxÞ:
n¼1
Substitute this expression into the Poisson equation, together with Equation 3.53 for
the current density. Since the expansion makes use of orthogonal functions, we can
equate the constant term and each term in the series independently. We find the
coefficients for the particular solution are
μ0 J 0 w
G2p ¼
2s
and
2μ0 J0 nkw
Ln ¼ 3 2 sin :
πn k 2
The three equations for A contain a total of nine unknown coefficients. We find the
values of these coefficients by imposing the boundary conditions.
Case 1: y ¼ 0:
The bottom boundary is an infinite permeability surface, so B ¼ By must be
perpendicular to this surface. Thus
∂A1
Bx ¼ ¼ 0;
∂y
which gives
Cn Dn ¼ 0: (3.54)
Case 2: y = g
The nonperiodic and periodic parts of the equations for the continuity of A and ∂y A
give the four equations
μ0 J 0 w 2
0 ¼ G0 þ G1 g g (3.55)
2s
62 Potential theory
2μ0 J0 nkw
Cn enkg þ Dn enkg ¼ Hn enkg þ Mn enkg þ sin (3.56)
πn3 k2 2
μ0 J 0 w
0 ¼ G1 g (3.57)
s
Cn enkg Dn enkg ¼ Hn enkg Mn enkg : (3.58)
Case 3: y ¼ g þ h
The nonperiodic and periodic parts of the equations for the continuity of A and ∂y A
give the four equations
μ0 J 0 w
G0 þ G1 ðg þ hÞ ðg þ hÞ2 ¼ E0 þ E1 ðg þ hÞ (3.59)
2s
nkðgþhÞ 2μ0 J0 nkw
Hn e nkðgþhÞ
þ Mn e þ 3 2 sin ¼ Fn enkðgþhÞ (3.60)
πn k 2
μ0 J 0 w
G1 ðg þ hÞ ¼ E1 (3.61)
s
Hn enkðgþhÞ Mn enkðgþhÞ ¼ Fn enkðgþhÞ : (3.62)
Case 4: y→∞
At large y, far from the conductor, the field must be uniform, so we have
∂A3
Bx ¼ ¼ E1 :
∂y
From far above, the conductor looks like an infinite current sheet with current density
J0 h, whose strength is reduced by the filling factor w/s and is enhanced by a factor 2 due
to the presence of the bottom permeable surface. Then using Equation 1.17, we find
w
E1 ¼ Bx ¼ ½ μ0 J0 h 2
s (3.63)
μ0 J0 hw
¼ :
s
Equations 3.54–3.63 give ten constraints on the nine unknown coefficients.
However, Equation 3.61 is redundant since it is equivalent to Equations 3.57 and
3.63. Thus we have nine equations in nine unknowns. After solving this system of
equations, the resulting vector potentials are:
X
∞
A1 ¼ α ðenky þ enky ÞcosðnkxÞ; (3.64)
n¼1
3.9 Boundary value problems 63
where
μ0 J0 nkw nkg nkh
α ¼ 3 2 sin e ðe 1Þ; (3.65)
πn k 2
μ0 J 0 w X∞
A2 ¼ ðy gÞ2 þ ðβ1 enky þ β2 enky þ β3 ÞcosðnkxÞ; (3.66)
2s n¼1
where
μ0 J 0 nkw nkðgþhÞ
β1 ¼ 3 2 sin e (3.67)
πn k 2
μ0 J 0 nkw nkðgþhÞ
β2 ¼ 3 2 sin e enkg þ enkg (3.68)
πn k 2
nkw
β3 ¼ 2μ0 J0 sin ; (3.69)
2
μ0 J0 hwð2g þ hÞ μ0 J0 hw X∞
A3 ¼ yþ γenky cosðnkxÞ; (3.70)
2s s n¼1
where
μ0 J 0 nkw nkðgþhÞ
γ ¼ 3 2 sin e þ enkðgþhÞ þ enkg enkg : (3.71)
πn k 2
The series for the vector potential converge rapidly because of the n3 factor in the
denominators of the coefficients.
We find B by taking the curl of the vector potential. The resulting field in the slot is
shown in Figure 3.7. The dotted lines show the location of the conductor.
where P1 is a Legendre polynomial. The potential for the problem including the sphere
must remain finite as r → ∞. Therefore the potential outside the sphere is given by
X
∞
Vmext ¼ V0 þ cn rn1 Pn ðcos θÞ:
n¼0
The potential must also remain finite at r = 0, so the potential inside the sphere has
the form
X
∞
Vmint ¼ V0 þ dn rn Pn ðcos θÞ:
n¼0
3.9 Boundary value problems 65
X
∞
dPn X∞
dPn
cn an2 ¼ dn an1 :
n¼0
dθ n¼0
dθ
Since this must true for any θ, the coefficients must satisfy the relation
cn ¼ dn a2nþ1 : (3.72)
We require that this relation hold for any value of n. For n = 0, we obtain
μ0 c0 a2 ¼ 0:
This shows that c0 ¼ 0 and from Equation 3.72 we find that d0 ¼ 0. For n = 1, we
find that
Using Equation 3.72 for cn, we obtain dn ¼ 0. Using this in Equation 3.72, we find
cn ¼ 0. Thus the only nonvanishing coefficients are c1 and d1 . The solution for the
potential outside the sphere is
μ μ0 a3
Vmext ¼ H0 r cos θ þ H0 cos θ
μ þ 2 μ0 r 2
66 Potential theory
μ μ0
Vmint ¼ H0 r cos θ þ H0 cos θ
μ þ 2 μ0
Since
and
∂ϕ
rϕ·^
n¼ ;
∂n
Equation 3.75 becomes
ð ð ð
∂ϕ
rψ·rϕ dV þ ψr ϕ dV ¼ ψ dS:
2
(3.76)
∂n
Subtracting Equation 3.77 from Equation 3.76, we obtain Green’s second identity
or Green’s theorem [19]
ð ð
∂ϕ ∂ψ
ðψr ϕ ϕr ψÞdV ¼
2 2
ψ ϕ dS: (3.78)
∂n ∂n
where P and Q are continuous functions of x and y and have continuous partial
derivatives.
To apply Green’s theorem to magnetostatics, let us choose the function ψ to be
proportional to the inverse distance between an element of current at r0 and a field
observation point at r.
68 Potential theory
1 1
ψ¼ ¼ 0 : (3.80)
4πR 4πj!
r !r j
In the case where currents are present inside V, we choose ϕ to be one of the
components of the vector potential, e.g., Ax . Then we have
r2 Ax ¼ μ0 Jx
Because of the delta function, we can solve this equation for Ax ðrÞ:
ð ð
! μ0 J x 0 1 1 ∂Ax ∂ 1
Ax ð r Þ ¼ dV þ Ax dS0 : (3.81)
4π R 4π R ∂n ∂n R
The volume integral represents a particular solution of the Poisson equation and only
includes the effects of currents inside V. Any additional currents outside V influence
the value of the surface integral. If we let r be a set of points on the surface S, then this
represents an integral equation for the unknown vector potential.
We can use Green’s theorem to develop integral equation solutions for Dirichlet
and Neumann boundary value problems. We generalize Equation 3.80 used in the
previous derivation by defining the Green’s function [21]
1
Gð!
r ;! ! !0
0
r Þ¼ ! ! 0 þ Lð r ; r Þ;
4πj r r j
where ϕ represents either the magnetic scalar potential or one of the components
of A.
References 69
For Dirichlet boundary value problems, if we choose L such that the Green’s
function
GD ¼ 0 (3.83)
For Neumann boundary value problems for the exterior of the volume V, if we
instead choose L such that
∂GN
¼0 (3.85)
∂n
when r0 is on the surface, we obtain the integral equation
ð !0
! 1 ! ! 0 ∂ϕð r Þ
ϕð r Þ ¼ GN ð r ; r Þ dS0 : (3.86)
4π ∂n
Thus once we have an appropriate Green’s function for a given geometry, the
potential at some field point r can be found by integration over the boundary
surface. The Green’s function is the solution of the Poisson equation for a delta
function source.
References
[1] D. Iencinella & G. Matteucci, An introduction to the vector potential, Eur. J. Phys.
25:249–256, 2004.
[2] J. D. Jackson, Classical Electrodynamics, Wiley, 1962, p. 141.
[3] W. Panofsky & M. Phillips, Classical Electricity and Magnetism, 2nd ed., Addison-
Wesley, 1962, p. 147–151.
[4] P. Lorrain & D. Corson, Electromagnetic Fields and Waves, 2nd ed., Freeman, 1970,
p. 305.
[5] J. D. Jackson, op. cit., p. 145–146.
[6] L. Eyges, The Classical Electromagetic Field, Dover, 1980, p. 166–167.
[7] P. Lorrain & D. Corson, op. cit., p. 92–94.
[8] L. Eyges, op. cit., p. 138–139.
[9] J. D. Jackson, op. cit., p. 47.
[10] J. B. Marion, Classical Electromagnetic Radiation, Academic Press, 1965, p. 58–60.
[11] W. Panofsky & M. Phillips, op. cit., p. 88–89.
[12] W. Smythe, Static and Dynamic Electricity, 2nd ed., McGraw-Hill, 1950, p. 170–182.
[13] Ibid., p. 187–198.
[14] W. Panofsky & M. Phillips, op. cit., p. 81–82.
[15] G. Arfken, Mathematical Methods for Physicists, 3rd ed., Academic Press, 1985,
p. 680–683.
70 Potential theory
71
72 Conductor-dominant transverse fields
! ! 1
B ðr; θÞ ¼ r A ¼ ^r ∂θ Az ^θ ∂r Az : (4.4)
r
An ¼ n Cn Fn
(4.6)
Bn ¼ n Cn En ;
The coefficients An and Bn describe the multipole field content of the transverse field.1
Returning to Equation 4.4, the azimuthal field component is
X
∞
Bθ ðr; θÞ ¼ Cn n rn1 ½En cos nθ þ Fn sin nθ:
n¼1
1
Caveat emptor. The reader should be aware that a number of different definitions are used in the literature to
describe the multipole content of a transverse field.
4.1 General solution to the Laplace equation in two dimensions 73
X
∞
Bθ ðr; θÞ ¼ rn1 ðAn sin nθ þ Bn cos nθÞ: (4.8)
n¼1
On the midplane of the magnet (θ = 0, r = x), Bx and By are given by power series in x
X
∞
μ0 Vm ðr; θÞ ¼ Gn rn ðHn cos nθ þ In sin nθÞ: (4.10)
n¼1
! ∂Vm ^ μ0 ∂Vm
B ðr; θÞ ¼ μ0 rVm ¼ ^r μ0 θ :
∂r r ∂θ
The field components calculated from this potential must equal the same quantities
calculated from the vector potential. We can make Br have the same form as
Equation 4.5 if we demand that
n Gn Hn ¼ n Cn Fn ¼ An
n Gn In ¼ n Cn En ¼ Bn :
Thus we identify
G n ¼ Cn
Hn ¼ Fn
I n ¼ En :
74 Conductor-dominant transverse fields
The resulting magnetic field components are still given by Equations 4.7 and 4.8.
and
ð 2π
I2 ¼ cos nθ cos mθ dθ ¼ π δmn
0
ð 2π (4.14)
¼ sin nθ sin mθ dθ:
0
Likewise, we can multiply Equation 4.8 with sin mθ and find the other set of
multipole components is
ð 2π
1
An ¼ n1 Bθ ðr; θÞ sin nθ dθ: (4.16)
πr 0
2
We will see in the next chapter that this relationship between the vector and scalar potentials follows directly
from the Cauchy-Riemann equations.
3
CRC 497, 502.
4.1 General solution to the Laplace equation in two dimensions 75
A similar analysis using Equation 4.7 for the radial component of the field gives
ð 2π
1
Bn ¼ n1 Br ðr; θÞ sin nθ dθ (4.17)
πr 0
and
ð 2π
1
An ¼ n1 Br ðr; θÞ cos nθ dθ: (4.18)
πr 0
The strength of the multipole fields provides a measure of the field quality.
Limits on the field uniformity are imposed by the application that needs the
magnetic field. The presence of harmonics of the desired field limits the size of
the useful magnet aperture. Sometimes, when examining the field quality of
a magnet, it is more useful to examine the relative magnitude of the multipole
coefficients with respect to the coefficient for the desired multipole. Thus for
a dipole design, for example, one could calculate the dimensionless quantities
Bn r0n1
bn ¼ ;
B1
where r0 is a reference radius, typically ~2/3 of the magnet aperture.
The boundary conditions for the vector potential in polar coordinates at some
radius rb can be determined from the boundary conditions on the magnetic field.
From the condition on the normal component of B, we have
ð1Þ ð2Þ
Br ¼ Br
1 ∂Að1Þ 1 ∂Að2Þ
z
¼ z
:
rb ∂θ rb ∂θ
From this relation, we know that Að2Þ can differ from Að1Þ by at most a constant,
which we can ignore since constants are removed when we take derivatives to
obtain the field. Thus we have
Að2Þ ð1Þ
z ðrb ; θÞ ¼ Az ðrb ; θÞ: (4.19)
1 ∂Að2Þ
z ðrb ; θÞ 1 ∂Að1Þ
z ðrb ; θÞ
þ ¼ K: (4.20)
μ ð2Þ ∂r μ ð1Þ ∂r
76 Conductor-dominant transverse fields
where the distance from the line current at (a, ϕ) to the observation point P located
at (r, θ) is given by
n o1=2
R ¼ ðr cos θ a cos ϕÞ2 þ ðr sin θ a sin ϕÞ2 (4.22)
and r0 is some constant reference radius for the two-dimensional potential. We look
for a harmonic expansion for Az. When r > a we extract a factor of r2 from the
logarithm and find that
a2 a
lnðRÞ ¼ ½ lnðr Þ þ ½ ln 1 þ 2 2 cos ðθ ϕÞ :
2
r r
Thus the vector potential for the line current when r > a is [1]
μI μI X∞
1 an
Az ðr; θÞ ¼ lnðrÞ þ cos ½nðθ ϕÞ (4.23)
2π 2π n¼1 n r
plus a constant term involving r0. The expansion of lnR for the case r < a can be
done in a similar manner by first extracting a factor of a2 from the argument of the
logarithm. This results in the vector potential
μI μI X∞
1 r n
Az ðr; θÞ ¼ lnðaÞ þ cos ½nðθ ϕÞ: (4.24)
2π 2π n¼1 n a
X∞
1 r n X∞
Az1 ¼ μ0 k lnðaÞ þ μ0 k cos nω þ μ0 Cn rn cos nω:
n¼1
n a n¼1
In the region (2) with a < r < R between the line current and the iron boundary, the
vector potential is
X∞
1 an X∞
Az2 ¼ μ0 k lnðrÞ þ μ0 k cos nω þ μ0 Cn rn cos nω:
n¼1
n r n¼1
Lastly, in the region (3) with r > R inside the iron, the vector potential is
X∞
1 an X∞
Az3 ¼ μ k lnðrÞ þ μ k cos nω þ μ Dn rn cos nω;
n¼1
n r n¼1
1
Br ¼ ∂ θ Az
r
1
Hθ ¼ ∂ r Az
μ
μ0 k n n1 μ k n n1
aR þ μ0 Cn Rn1 ¼ a R þ μ Dn Rn1
n n
Cn Rn1 ¼ Dn Rn1 :
This gives two equations in two unknowns, which can be solved to give,
k a n μ μ0
Cn ¼
n R2n μ þ μ0
k n μ μ0
Dn ¼ a :
n μ þ μ0
In the region inside the radius a of the line current, the vector potential can be
expressed as
a 2n
μ0 I μ0 I X
∞
1 r n
Az1 ¼ lnðaÞ þ 1þα cos nω; (4.25)
2π 2π n¼1 n a R
4.2 Harmonic expansion for a line current 79
where
μ μ0
α¼ :
μ þ μ0
The term in brackets shows the enhancement factor due to the iron. For points inside
the iron, the solution is
" #
μI X∞
1 an
Az3 ¼ lnðrÞ þ ð1 αÞ cos nω : (4.26)
2π n¼1
n r
We are now in a position to relate the solution of the previous example with the
method of images for a circular boundary that we discussed in Section 2.7. For the
region interior to the filament, r < a, the vector potential can be written in the form
μ0 I μ IX ∞
1 r n μ IX ∞
1 nα
Az1 ¼ lnðaÞ þ 0 cos nω þ 0 r n cos nω; (4.27)
2π 2π n¼1 n a 2π n¼1 n r1
where
R2
r1 ¼ :
a
Comparing with Equation 4.24, the first two terms give the potential for the true
line current at r = a. Since in magnetostatics constant terms in the potential have no
physical effects, we can arbitrarily add to the potential a constant term
μ0 I
lnðr1 Þ:
2π
Then this term plus the last term in Equation 4.27 give the potential for a line current
at r1 with current α I. Thus the vector potential in region (1) can be written as
Az1 ¼ ALC ðaÞ þ α ALC ðr1 Þ;
where ALC is the vector potential for a line current. For the region inside the iron,
r > R, let us define
2μ0
β¼1α¼ :
μ þ μ0
Using Equation 4.23 for the vector potential of a line current, we find
" !#
μI X∞
1 an
Az3 ¼ ð1 βÞ ln r þ β ln r þ cos nω
2π n¼1
n r
¼ α ALC ð0Þ þ β ALC ðaÞ:
We see that the coefficients α and β are the same as those for the image currents in
Equations 2.27–2.30.
Assume the observation point P is in the x-y plane, as shown in Figure 4.4.
We have
! ¼ dI ^z
K
ds
ρ ¼ r2 þ a2 2 r a cos ðθ ϕÞ
2
R2 ¼ ρ2 þ z2
dS ¼ ds dz;
Figure 4.3 Current sheet with cross-section along the curve C and extending
infinitely along the z direction.
4.3 Field for a current sheet 81
where a ¼ aðsÞ, ϕ ¼ ϕðsÞ and s is the arclength around the sheet. The vector
ρ is the distance between the current element and the field point in the x-y plane.
Thus
ð ð∞
! ðr; θÞ ¼ μ0 KðsÞ ^z ð!ρ þ z ^z Þ
B 3=2
dz ds
4π ∞ fρ2 þ z2 g
ð
μ0
¼ KðsÞ IðρÞ ^z !ρ ds;
4π
where4
ð∞
dz 2
IðρÞ ¼ ¼ : (4.28)
∞ fρ2 þ z2 g3=2 ρ2
Thus we find that the field from the current sheet is given by
ð
! μ0 ^z ! ρ
B ðr; θÞ ¼ KðsÞ ds: (4.29)
2π ρ 2
It will also be useful to have an expression for the vector potential of a current
sheet. Assuming the sheet is composed of parallel line currents and using
Equation 4.21, we have
ð
μ0 ρ
Az ðr; θÞ ¼ KðsÞln ds: (4.30)
2π ρ0
4
GR 2.271.5.
82 Conductor-dominant transverse fields
ds ¼ a dϕ
dI dϕ 1 dI
K¼ ¼ :
dϕ ds a dϕ
!
ρ ¼ ðxo xÞ ^x þ ðyo yÞ ^y
dI
K¼ :
dy
All the current elements on both sheets give positive By because the magnitude of xo
is smaller than a. Thus we can write
μ0
By ¼ K ½ða xo Þ I1 þ ða þ xo Þ I2 ;
2π
where5
ð∞
dy π
I1 ¼ ¼
∞
2
ðxo aÞ þ ðyo yÞ 2 a xo
ð∞
dy π
I2 ¼ ¼ :
∞
2
ðxo þ aÞ þ ðyo yÞ 2 a þ xo
I3 ¼ ½ ln½D∞∞ þ yo I1 :
Using l’Hopital’s rule, we can show that the first term vanishes, and find that
I3 ¼ yo I1
I4 ¼ yo I2 :
5 6
GR 2.172. GR 2.175.1.
84 Conductor-dominant transverse fields
Substituting back in, we find that the horizontal field between the sheets
vanishes.
Bx ¼ 0:
Thus there is a pure dipole field in the region between the current sheets.
dI
¼ I0 cos mϕ;
dϕ
where I0 is the current flowing at the midplane (ϕ = 0). We obtain the vector
potential for the sheet by integrating the weighted distribution of the vector
potential for a line current. Let us first consider the case where the observation
point (r, θ) has r < a. Using Equation 4.24 for the line current and ignoring the
constant term, the vector potential for the multipole sheet is
ð
μI0 X
∞
1 r n 2π
Az ¼ cos mϕ cos ½nðθ ϕÞ dϕ:
2π n¼1 n a 0
Expanding the integrand and using Equations 4.13 and 4.14, the integral has the
value
ð 2π
0 if m ≠ n
cos mϕ cos ½nðθ ϕÞ dϕ ¼ : (4.34)
0 π cos mθ if m ¼ n
μI0 r m
Az ðr; θÞ ¼ cos mθ (4.35)
2m a
and the components of the magnetic field are
μI0 r m1
Br ðr; θÞ ¼ sin mθ
2a a
(4.36)
μI0 r m1
Br ðr; θÞ ¼ cos mθ:
2a a
4.4 Ideal multipole current sheet 85
showing that the cos θ angular distribution also produces a pure vertical field in the
magnet aperture.
For the case when r > a, we use Equation 4.23 for the vector potential of the line
current and obtain
ð ð
μI0 2π μI0 X∞
1 an 2π
Az ¼ cos mϕ ln r dϕþ cos mϕ cos ½nðθ ϕÞ dϕ:
2π 0 2π n¼1 n r 0
The first integral vanishes over a complete circle and the second integral can again
be evaluated using Equation 4.34. Thus the vector potential of the multipole sheet
for the region r > a is
μI0 am
Az ðr; θÞ ¼ cos mθ (4.37)
2m r
and the magnetic field is
μI0 am
Br ðr; θÞ ¼ sin mθ
2 rmþ1
(4.38)
μI0 am
Br ðr; θÞ ¼ cos mθ:
2 rmþ1
To determine the values for C and D, we demand that Br and Hθ are continuous across
the iron boundary.
μ0 I0 a m μI0 a m
þ C m Rm ¼ þ D m Rm
2 R 2 R
Rm Rm
C ¼ D :
μ0 μ
Solving these two equations, we find that the unknown coefficients are
μ0 I0 a m μ μ0
C¼
2 m Rm R μ þ μ 0
μI0 a m μ μ0
D¼ :
2 m Rm R μ þ μ0
Using these values, the vector potential is now known in the three regions. The effect
of the iron can be summarized by defining the iron enhancement factor
μr 1 a 2m
αm ¼ 1 þ ; (4.39)
μr þ 1 R
which agrees with the enhancement factor from Equation 4.25. We can write the
vector potential inside the magnet aperture as
μ0 I0 r m
Az1 ðr; θÞ ¼ αm cos mθ: (4.40)
2m a
μ0 I0 r m1
Br ðr; θÞ ¼ αm sin mθ
2a a
(4.41)
μ I0 r m1
Br ðr; θÞ ¼ 0 αm cos mθ:
2a a
The iron enhancement factor in Equation 4.39 ignores any saturation effects
in the iron. When saturation becomes significant, the enhancement factor for the
dipole field is decreased. In addition, the saturation of the iron does not occur
uniformly. This causes changes in the field from the azimuthal distribution of
the enhanced currents, leading to sextupole and higher multipole errors in the
field in the magnet aperture. Fortunately, there are techniques, such as modify-
ing the iron shape, which can be used to adjust the value of B3 at a fixed
operating current.[5]
4.4 Ideal multipole current sheet 87
Example 4.6: cos mθ sheet in circular iron cavity using the scalar potential
It is instructive to use an alternative method of solving the preceding boundary
value problem. In this case, we will use the scalar potential and take into account
the presence of the current sheet through the addition of another pair of boundary
conditions. Unlike the previous example, the unknown coefficients here take into
account both the field from the sheet and the field from the images in the iron.
We know from the boundary condition across a current sheet, Equation 2.23, that
the angular dependence of the current must match the angular dependence of
the fields on either side of the sheet. Since the fields are given by the derivative
of the potential and the current goes like cos mθ, this implies that we must use
the sine term in Equation 4.2 for the potential. Thus the scalar potentials for the three
regions are
Vm1 ¼ A rm sin mθ
Vm2 ¼ ðB rm þ C rm Þ sin mθ
Vm3 ¼ D rm sin mθ:
A am ¼ B am C am
mðB am þ C am Þ þ m A am ¼ I0 ;
I0
A¼ αm
2mam
I0 am μ μ0
B¼
2mR2m μ þ μ0
I 0 am
C¼
2m
μ0 I0 am
D¼ :
mðμ þ μ0 Þ
Evaluating the field components inside the magnet aperture, we again obtain
Equation 4.41.
88 Conductor-dominant transverse fields
μI X∞
1 r m
Az ðr; θÞ ¼ cos ½mðθ ϕÞ:
2π m¼1 m a
μI X∞ m1
r
Bθ ðr; θÞ ¼ cos ½mðθ ϕÞ : (4.42)
2π m¼1 am
μI X∞
1
Bn ¼ π cos mϕ δmn
2π m¼1 am
2
(4.43)
μI
¼ cos nϕ :
2πan
For a current sheet, we can generalize this by integrating over the current
distribution.
ð
μ 2π dI
Bn ¼ cos nϕ dϕ: (4.44)
2πan 0 dϕ
We can find the skew multipoles in a similar manner using Equation 4.16.
ð
μ 2π dI
An ¼ sin nϕ dϕ: (4.45)
2πan 0 dϕ
pairs in the multipole magnet. Thus N = 1 refers to a dipole magnet because it has one
pair of poles, N = 2 to a quadrupole, etc. The current in the magnet goes in opposite
directions on the adjacent sides of a pole, as shown in Figure 4.6. This corresponds to
the way these coils are usually wound, with the cable bending around the pole and
returning in the opposite direction. As we can see from Equations 4.44 and 4.45, the
multipole coefficients are weighted sums of the current distribution. The Bn coeffi-
cients are referred to as the normal multipoles. They are largest when the magnitude
of the current reaches a maximum on the midplane. The An coefficients are referred
to as the skew multipoles, which are largest when the current changes sign at the
midplane. The field in a skew multipole of order N has the same pattern as the normal
multipole rotated by π/2N. For example, a normal dipole has a vertical field, while
a skew dipole has a horizontal field.
Consider the current distribution for an ideal multipole of order N given by
dI
¼ I0 cos Nϕ:
dϕ
dI
¼ I0 cos 2ϕ:
dϕ
Thus
ð ϕ2
μ I0
B2 ¼ 0 2 cos2 ð2ϕÞ dϕ
2π a ϕ1
μ0 I 0
¼ :
2 a2
Let S( ) = sin( ). Then Bn can be written as the sum of eight sine terms.
μ0 I
Bn ¼ ½ Sðnϕ2 Þ Sðnϕ1 Þ Sðnπ nϕ1 Þ þ Sðnπ nϕ2 Þ
2πnan
Sðnπ þ nϕ2 Þ þ Sðnπ þ nϕ1 Þ Sðnϕ1 Þ þ Sðnϕ2 Þ :
If n is even, the terms in the brackets cancel, so Bn = 0. For odd values of n, we get
2μ0 I
Bn ¼ ½sin nϕ2 sin nϕ1 : (4.46)
πnan
A similar calculation shows that An = 0 for both odd and even values of n. Thus the
allowed harmonics for this current configuration are just the Bn , for odd values of n.
In addition to the desired harmonic B1 that represents the dipole, the approxima-
tion of the ideal multipole distribution above introduces other harmonics, which
represent errors to the desired field. The dominant allowed error here is the
sextupole term B3 . The angle ϕ1 is typically set as close to 0 as possible in order
to maximize the dipole strength. The angle ϕ2 can then be used to remove the
sextupole component from the field. Since
2μ0 I
B3 ≃ sin 3ϕ2 ;
3πa3
we can eliminate the sextupole term by setting ϕ2 ¼ π=3. Removing other allowed
higher harmonics from the field requires additional degrees of freedom. A number
92 Conductor-dominant transverse fields
where
π π
β¼2 ϕ ¼ ϕ:
2N N
To get a net contribution to a normal multipole of order n
dI
Bn ∝ ðϕÞ cos nϕ
dϕ
from the current on both sides of the pole, the cosine function must also change sign.
h π i
cos n ϕ ¼ cos nϕ
N
nπ nπ
cos cos nϕ þ sin sin nϕ ¼ cos nϕ:
N N
A similar argument shows that the allowed skew multipoles also satisfy
Equation 4.47.
The symmetry of the current distribution is directly related to the allowed
harmonics. Let KðϕÞ ¼ dI=dϕ, for example, and assume the current distribution
is up-down symmetric, so that
KðϕÞ ¼ KðϕÞ:
Thus the fact that the dipole approximation had An = 0 follows from the up-down
symmetry of the current distribution that we used. The consequences of some other
symmetries for current distributions are listed in Table 4.1. These symmetries are
inevitably violated to some extent in building an actual magnet and this leads to the
presence of “nonallowed” multipoles in the fields. Random errors in the construc-
tion of the magnet can introduce values for any multipole.[9] If these nondesired
multipoles exceed their tolerances, they must be removed by modifications in the
manufacturing process or by introducing correction coils.
Assume an observation point P is located at (r, θ) in the x-y plane and the current
element is at (a, ϕ), as shown in Figure 4.8. Let the vector ρ be the distance between
the current element and the field point in the x-y plane and let σ be the current
density in the conductor. Then we have
! ¼ σ ^z
J
! !
R ¼ ρ þ z ^z
dV ¼ dS dz:
where IðρÞ is given by Equation 4.28. We find that the field from the current block is
given by
ð
! μ0 ^z ! ρ
B ðr; θÞ ¼ σ dS: (4.48)
2π ρ 2
The vector potential for a current block can be found by integrating the vector
potential for the line current, Equation 4.21, over the area of the block.
ð
μ0 ρ
Az ðr; θÞ ¼ σ ln dS: (4.49)
2π ro
4.7 Field for a block conductor 95
μ0 J
Bx ¼ ½r2 sin θ2 þ r1 sin θ1 ¼ 0
2
μJ μJ
By ¼ 0 ½r1 cos θ1 r2 cos θ2 ¼ 0 2c;
2 2
where 2c is the separation between the centers of the two circles. Thus the field in the
overlap region is a pure dipole. The strength of the field is proportional to the separation
between the circles. In the region where the two coils overlap, the net current is zero.
Thus the conductor in the overlap region can be removed without affecting the field there.
Next we examine the coil thickness t as a function of the azimuthal angle ϕ, as
shown in Figure 4.10. The coil thickness at some angle ϕ is the distance between the
points P1 and P2 . Point P1 is determined by the intersection of circle 1
ðx þ cÞ2 þ y2 ¼ a2
96 Conductor-dominant transverse fields
We find that
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
c þ a2 m2 c2 þ m2 a2
x1 ¼ :
1 þ m2
Point P2 is determined by the intersection of circle 2
ðx cÞ2 þ y2 ¼ a2
with the straight line. We find that the expression for x2 is the same as the one for x1,
except that the first term in the numerator is +c instead of –c. Thus we have
2c
Δx ¼ x2 x1 ¼
1 þ tan2 ϕ
Δy ¼ y2 y1 ¼ Δx tan ϕ:
Thus the overlapping circular conductors represent another form of a cosine current
distribution. Quadrupole fields can be designed in the same manner using over-
lapping elliptical conductors.[8]
4.7 Field for a block conductor 97
Example 4.9: on-axis field of an annular sector with constant current density
Consider the annular sector conductor shown in Figure 4.11. For a field point on the
axis we have
!
ρ ¼ a cos ϕ ^x a sin ϕ ^y :
! μσ
B ð0; 0Þ ¼ 0 ða2 a1 Þ½ðsin ϕ2 sin ϕ1 Þ^y þ ðcos ϕ2 cos ϕ1 Þ^x :
2π
(4.50)
We see that the field for the case of constant current density is directly proportional to
the radial thickness of the conductor.
! μ0 σ ½ðxo xÞ ^y ðyo yÞ ^x
dB ¼ dS: (4.51)
2π ðxo xÞ2 þ ðyo yÞ2
Looking at By , we have
ðb ðd
μ0 σ ðxo xÞ
By ¼ dx dy:
2π a c ðxo xÞ2 þ ðyo yÞ2
7 8
GR 2.175.1. GR 2.733.1.
4.7 Field for a block conductor 99
μ0 σ
By ¼ ½I2 ðc; d; bÞ I2 ðc; d; aÞ
4π
μσ yo c
¼ 0 ðyo cÞ ln½ðxo bÞ2 þ ðyo cÞ2 þ 2jxo bj tan1
4π jxo bj
2 2 1 yo d
ðyo dÞ ln½ðxo bÞ þ ðyo dÞ 2jxo bj tan
jxo bj
2 2 1 yo c
ðyo cÞ ln½ðxo aÞ þ ðyo cÞ 2jxo aj tan
jxo aj
2 2 1 yo d
þ ðyo dÞ ln½ðxo aÞ þ ðyo dÞ þ 2jxo aj tan :
jxo aj
From the symmetry of Equation 4.51, we see that the result for Bx is the negative of
the result for By with the substitutions
x↔y
xo ↔yo
ða; bÞ↔ðc; dÞ:
Thus Bx is given by
μ0 σ 2 2 1 xo a
Bx ¼ ðxo aÞ ln½ðxo aÞ þ ðyo dÞ þ 2jyo dj tan
4π jyo dj
xo b
ðxo bÞ ln½ðxo bÞ2 þ ðyo dÞ2 2jyo dj tan1
jyo dj
2 2 1 xo a
ðxo aÞ ln½ðxo aÞ þ ðyo cÞ 2jyo cj tan
jyo cj
2 2 1 xo b
þ ðxo bÞ ln½ðxo bÞ þ ðyo cÞ þ 2jyo cj tan :
jyo cj
Using these expressions for the field, we show in Figure 4.13 a scan of the vertical
field component along the x axis for a square conductor centered at the origin.
The calculation using these equations fails when the observation point is located at
one of the four corners of the rectangle.
Figure 4.13 Vertical field from a square conductor with sides of length 2 centered
at the origin.
The skew multipoles can be found in a similar manner and are given by
ðð
μ sin nϕ
An ¼ Jða; ϕÞ n1 da dϕ: (4.53)
2π a
Case 1: r < a1
The vector potential, Equation 4.35, for the ideal multipole sheet is
μI0 r m
Azsh ðr; θÞ ¼ cos mθ:
2m a
When m = 2, we have
μJ0 r2 a2
Az ¼ cos 2θ ln : (4.55)
4 a1
Case 2: r > a2
When the observation point is always larger than the radius of any of the multipole
sheets, we use Equation 4.37 for the vector potential of the ideal multipole sheets.
μI0 am
Azsh ðr; θÞ ¼ cos mθ:
2m r
Figure 4.14 Radial scan of the vertical component of the magnetic field of an ideal
dipole (m =1, a1 = 0.2 m, a2 = 0.3 m, J0 = 100 A/mm2).
The components of the magnetic field in each of these regions is easily com-
puted from these vector potentials. Figure 4.14 shows a radial scan of the
vertical component of the magnetic field along the midplane for an ideal dipole.
The current block extends from 0.2 m to 0.3 m in this example. The field is
constant and purely vertical inside the aperture, as expected. The field changes
direction inside the current block. The strength of the field slowly falls off
outside the block.
9
We will consider H inside a magnetized body in Chapter 9.
4.9 Field from a magnetized body 103
The gradient operator only acts on unprimed coordinates, so the second and fourth
terms on the right-hand side vanish.
In the first term, we have using Equation B.6
!
ρ 1 1
!
r 2 ¼ 2 r ρ þr 2 ! ρ: (4.61)
ρ ρ ρ
The first term here vanishes because ρ is radial and so its curl vanishes. In the
second term
1 1 1
r 2 ¼ ^x ∂x 2 þ ^y ∂y 2
ρ ρ ρ
2 !
¼ 4 ρ:
ρ
Thus the cross-product of the last two factors in Equation 4.61 is 0, and the second
term also vanishes.
10
GR 2.271.5.
104 Conductor-dominant transverse fields
Thus only the third term in Equation 4.60 survives and we have
ð !
! 1 !0 ρ
Hm ¼ ð M ·rÞ 2 dS0 : (4.62)
2π ρ
Writing out the vector ρ in terms of its x and y components, taking the derivatives,
combining terms, and dropping the primes, we find that the two-dimensional field
from an iron element is given by
ð "! ! #
! 1 M 2ðM · ! ρÞ !
Hm ¼ ρ dS: (4.63)
2π ρ2 ρ4
4.10 Superconductors
High-field magnets are usually energized using superconducting cables. In the
superconducting state, the resistivity vanishes, so a large current can flow through
the magnet coils without losing power due to Joule heating. In a superconducting
material, attractive forces between pairs of electrons are transmitted through
vibrations in the material lattice.[11] The operating conditions for superconductiv-
ity lie below the surface of a three-dimensional space of temperature, magnetic
field, and current density. The limiting values on each of the three axes are called
the critical values. When the superconductor is not in the superconducting state, it
is said to be in the normal state. There are several classes of superconducting
materials. Type I materials exhibit the Meissner effect, where any external magnetic
field is excluded from the interior of the superconductor. A magnetization is
generated in the superconductor that just cancels the external field. This remains
true as the external field is increased until it reaches the critical field Hc. Type I
superconductors are perfect diamagnetic materials.
Some alloys of intermetallic compounds form what are called type II super-
conductors. Two important examples are NbTi and Nb3 Sn. The critical current
density for these materials at 4 K is shown as a function of the magnetic flux density
in Figure 4.15.[12] NbTi is a useful material for fields up to ~9 T at 4 K, whereas
Nb3 Sn can be used up to ~22 T. These materials have two critical magnetic fields,
Hc1 and Hc2 , which can be much larger than Hc1 . For H < Hc1, the material
4.11 End fields 105
10 000
JE [A / mm2]
1000
Nb3Sn
NbTi
100
5 10 15 20 25
B [T]
completely excludes the external flux and it behaves like a type I superconductor.
For Hc1 < H < Hc2, flux begins to penetrate into the material and the magnitude of
the magnetization begins to fall, but electrically it still has zero resistivity.
The magnetic flux enters the superconductor in the form of discrete, quantized
flux lines known as fluxoids. The field still vanishes in the material surrounding the
fluxoids. The fluxoids can move because of Lorentz forces, creating heat. To stop
the fluxoids from moving, inhomogeneities known as pinning centers must be
introduced into the lattice. Finally when H > Hc2 , the material returns to the
normal state.
Magnets with superconducting cables must deal with the problem of persistent
currents.[13] As the current in the magnet is ramped up or down, eddy currents11
are induced in the superconductor. These induced shielding currents produce a field
that opposes the change in the field caused by the magnet’s power supply. Because
of the lack of resistance in a superconductor, the decay times for the induced
currents are very long. The persistent currents produce undesired sextupole and
higher multipoles that can be particularly significant at low field values.
11
Eddy currents will be discussed in more detail in Section 10.4.
106 Conductor-dominant transverse fields
in the opposite direction at the symmetrical location on the other side of the pole.
There are several standard end configurations. If the aperture in the end regions
does not need to be open, the simplest configuration is the racetrack coil.[10]
The coil bends around an arc, maintaining the same vertical position as the coils in
the straight part of the magnet. The bedstead end bends the conductors by 90° as
quickly as possible and then crosses over the pole at a fixed z location. This
configuration does keep a clear aperture in the end regions, but the sharp bend in
the conductor may not be acceptable for mechanical reasons. Another type of end
that maintains a clear aperture is the saddle end.[10] In this case, the conductor
turns cross over the pole following an arc that is spread out over z. The magnetic
design of real coil ends is usually done numerically using the Biot-Savart equation.
The end turns introduce additional multipole contributions to the field of
a magnet. However, it is possible to define a new set of multipole coefficients
defined in terms of the field components integrated along the axis of the magnet.
Flux theorems have been developed that can relate these integrated multipoles to
the geometry of the end turns.[14] The z locations where the cross-over for various
conductors begin can be adjusted using spacers to help balance the integrated
multipoles in the magnet.[15]
References
[1] W. Smythe, Static and Dynamic Electricity, 2nd ed., McGraw-Hill, 1950, p. 65.
[2] R. Gupta, Field Calculations and Computations, School at Centre for Advanced
Technology, Indore, India, 1998, p. 32.
[3] N. Schwerg & C. Vollinger, Analytic models for the calculation of the iron yoke
contribution in superconducting accelerator magnets, CERN Accelerator Technology
Department report CERN/AT 2007-33, 2007.
[4] P. Schmuser, Superconducting magnets for particle accelerators, in M. Month &
M. Dienes (eds.), The Physics of Particle Accelerators, AIP Conf. Proc. 249, vol. 2,
1992, p. 1106.
[5] G. Morgan, Shaping of magnetic fields in beam transport magnets, in M. Month &
M. Dienes (eds.), The Physics of Particle Accelerators, AIP Conf. Proc. 249, vol. 2,
1992, p. 1256–1259.
[6] M. Wilson, Superconducting Magnets, Oxford University Press, 1983, p. 32.
[7] J. Coupland, Dipole, quadrupole and higher order fields from simple coils, Nuc. Instr.
Meth. 78:181, 1970.
[8] A. Jain, Basic theory of magnets, in S. Turner (ed.), CERN Accelerator School on
Measurement and Alignment of Accelerator and Detector Magnets, CERN 98–05,
1998, p. 25.
[9] G. Parzen, Random errors in the magnetic field of superconducting dipoles and
quadrupoles, Part. Acc. 6:237, 1975.
[10] M. Wilson, op. cit., p. 27–28.
[11] Ibid., p. 279.
[12] https://nationalmaglab.org/magnet-development/applied-superconductivity-center.
References 107
In this chapter, we continue the discussion of transverse fields that are determined
by the location of the conductors. In the “central” region, far from the magnet ends,
the powerful methods of complex analysis1 can be applied to the calculation of
potentials, magnetic fields, multipoles and forces. Many of the topics in this chapter
are based on a series of important papers by Richard Beth and by Klaus Halbach.
Beginning with the field from a line current, we consider methods for calculating
the fields from current sheets. Then we use the complex form of Green’s theorem to
express the fields of block conductors in terms of contour integrals.
1
A brief summary of some important results from the theory of complex variables is given in Appendix E.
108
5.1 Complex representation of potentials and fields 109
These equations can be put into the form of the Cauchy-Riemann equations by
associating
u ¼ Az
v ¼ μ 0 Vm :
Thus in two dimensions, the vector and scalar potentials are related to each other as
the real and imaginary parts of the complex potential function
WðzÞ ¼ Az þ iμ0 Vm : (5.3)
Now consider the derivative of the complex potential. A complex function with
a continuous derivative is known as an analytic function. A complex derivative
must give the same result independent of the manner that Δz approaches 0. In the
case when Δz ¼ Δx, we have
dW ∂W ∂Az ∂Vm
¼ ¼ þ iμ0 ¼ By iBx :
dz ∂x ∂x ∂x
If we had chosen Δz ¼ iΔy instead, we would obtain the same expression for B.
So in either case we find
dW
i ¼ Bx iBy :
dz
Defining the complex magnetic field as2
we find the relation between the magnetic field and the potential is
dW
B ðzÞ ¼ i ; (5.5)
dz
2
Unfortunately, Beth and Halbach use different definitions for the complex magnetic field H and use different
systems of units, so some care must be exercised in comparing their results. We follow Halbach’s conventions
here in defining the components of H the same way as normal complex variables and using the SI system of
units.
110 Complex analysis of transverse fields
From analytic geometry, this is the condition that indicates that two lines are
perpendicular. Thus the real and imaginary parts of an analytic function describe
orthogonal curves. This indicates in particular that the equipotential lines for Az and
Vm cross at right angles.
We can use Equation 5.8 to write the magnetostatic Maxwell equations in complex
coordinates. The divergence equation
∂x Hx þ ∂y Hy ¼ 0
becomes [1]
∂H ∂H
þ ¼ 0; (5.10)
∂z ∂z
where H ¼ Hx þ iHy is the complex magnetic field intensity. The curl equation
∂x Hy ∂y Hx ¼ Jz ≡ σ
∂2
r2 ¼ 4 : (5.14)
∂z∂z
5.3 Field from a line current 113
Using Equations 5.5 and 5.15 with the derivative operating on the coordinates of
the field point zo , we get the magnetic field
μ0 I 1
Bx i By ¼ i
2π zo z
μI
¼ i 0 eiα (5.18)
2πR
μI
¼ i 0 ðcos α i sin αÞ;
2πR
which agrees with the result in Equation 5.16.
From Equation 5.15, the complex logarithm is
lnðzo zÞ ¼ lnR þ iα:
Using Equations 5.3 and 5.17 we can confirm that the vector potential for a line
current is
μ0 I
Az ¼ lnR (5.19)
2π
3
See Appendix E.
5.4 Field from a current sheet 115
þ
1
dzo ¼ 2πi:
zo z
zðsÞ ¼ z1 þ seiθ :
It follows that dz ¼ eiθ ds and
z2 ¼ z1 þ beiθ : (5.25)
We can write
Substituting for eiθ from Equation 5.25 and evaluating the integral4 gives
μ0 I
Wðzo Þ ¼ ½u2 lnu2 u1 lnu1 u2 þ u1 :
2πðz2 z1 Þ
μ0 I
Wðzo Þ ¼ ½u2 ln u2 u1 ln u1 : (5.26)
2πðz2 z1 Þ
Summing up the contributions to the magnetic field from the field of individual
line currents given in Equation 5.21, we find the field of a current sheet is given by
4
CRC 377.
5.4 Field from a current sheet 117
ð
i KðsÞ
H ðzo Þ ¼ ds; (5.27)
2π zðsÞ zo
where s is the distance along the sheet. If the position along the sheet is specified by
the polar angle ϕ, this can be written as
ð
i dI=dϕ
H ðzo Þ ¼ dϕ: (5.28)
2π zðϕÞ zo
It is possible to determine a unique current distribution dI=dϕ for circular or elliptic
current sheets that can produce any desired two-dimensional field compatible with
Maxwell’s equations in the magnet aperture.[6]
Example 5.2: field due to a circular arc sheet with constant current density
Let us consider a current sheet in the form of a circular arc, as shown in Figure 5.3.
The magnetic field from the sheet is given by Equation 5.28.
ð ϕ2
i dI dϕ
H ðzo Þ ¼ :
2π dϕ ϕ1 zo zðϕÞ
5
GR 2.172.
118 Complex analysis of transverse fields
From Equation 5.28, the field of the arc conductor at the origin is
ð ϕ2
i dI dϕ
H ð0Þ ¼
2π dϕ ϕ1 a eiϕ
1 dI iϕ2
¼ ½e eiϕ1 :
2πa dϕ
If the angular arc completes a full circle, we have a current shell and the integral in
Equation 5.29 becomes
þ
dz
I¼ :
zðz zo Þ
A point where the denominator of the integrand becomes zero is called a pole. If zo is
inside the circle, then the contour integral has simple poles at z = 0 and z ¼ zo .
The residue7 for the pole at z ¼ 0 is
1 1
lim z ¼
z→0 zðz zo Þ zo
and the residue for the pole at z ¼ zo is
1 1
lim ðz zo Þ ¼ :
z→zo zðz zo Þ zo
Therefore, by the residue theorem, the value of the integral is zero and the field inside
the shell vanishes. When zo is outside the shell, the integral only has the pole at z = 0
and the residue theorem gives
6 7
GR 1.622.3. See Appendix E.
5.4 Field from a current sheet 119
1 2πi
I ¼ 2πi ¼ :
zo zo
Therefore, the field outside the shell is
i dI
H ðzo Þ ¼ :
zo dϕ
Since the total current is
dI
I ¼ 2π ;
dϕ
the field can be written as
I
H ðzo Þ ¼ i :
2πzo
This is the same as Equation 5.21 for the field of a line current located at the center of
the shell.
Let us apply the Ampère law, Equation 5.22, for an infinitesimal rectangular
contour across a current sheet, as shown in Figure 5.4. Then we have,
where dI is the current enclosed in the contour. In the limit where the distance
perpendicular to the sheet approaches 0, the path of the observation
points approach the path along the sheet and this results in the “current sheet
theorem.”[7]
dI
H1 ðzÞ H2 ðzÞ ¼ : (5.31)
dz
In addition to determining fields, this result has been used for calculations of
magnetic stored energy and Lorentz forces.[8]
Example 5.3: field from cos ϕ current distribution using contour integration
Consider a closed circular sheet of radius a. We have
dI
¼ I0 cos ϕ
dϕ
z ¼ aeiϕ
eiϕ þ eiϕ z þ z
cos ϕ ¼ ¼ :
2 2a
Substituting into Equation 5.28, we get
þ
i I0 z þ z dz
H ðzo Þ ¼
4πa z zo i z
þ
I0 1 z
¼ þ dz (5.32)
4πa z zo zðz zo Þ
I0
¼ ½I1 þ I2 :
4πa
Using the method of partial fractions,[9] the denominator of the second integral can
be written
1 A B
¼ þ
zðz zo Þ z z zo
1 ¼ ðz zo ÞA þ zB:
Equating powers of z, we find that
1
A¼
zo
1
B¼ :
zo
Then we can write
I2 ¼ I3 þ I4 ;
5.5 cos ϕ current sheets 121
where
þ
1 z
I3 ¼ dz
zo z
ð
ia 2π iϕ
¼ e dϕ ¼ 0
zo 0
and
þ
1 z
I4 ¼ dz: (5.33)
zo z zo
Using z ¼ a2 =z ; we can write this as
þ
a2 dz
I4 ¼
zo zðz zo Þ
þ þ
a2 1 dz 1 dz
¼ þ :
zo zo z zo z zo
For zo inside the contour, the factor in square brackets vanishes because of the
residue theorem and I4 ¼ 0: Then from Equation 5.32,
H ðzo Þ ¼ Hx iHy
I0 I0
¼ ½2πi þ 0 ¼ i :
4πa 2a
From this, we see that the field inside the current sheet is
Hx ¼ 0
I0 (5.34)
Hy ¼ :
2a
The field is only in the vertical direction and has constant magnitude everywhere
inside the sheet in agreement with Equation 4.36.
For zo outside the contour, we have8
þ
z
I4 ¼ dz
zo ðz zo Þ
ð
ia2 2π 1
¼ dϕ
zo 0 a e zo
iϕ
a2
¼ 2πi:
z2o
8
GR 2.313.1.
122 Complex analysis of transverse fields
Substituting these results back into Equation 5.32, the field outside the contour is
I0 a2 I0 a
H ðzo Þ ¼ 0 2 2πi ¼ i 2 :
4πa zo 2zo
If we write zo ¼ x þ iy and multiply the numerator and denominator by ðzo Þ2, we find
that
I0 a I0 a
H ðzo Þ ¼ 4
xy i 4 ðx2 y2 Þ:
r 2r
Equating real and imaginary parts, we find the Cartesian field components outside the
sheet are
I0 a I0 a
Hx ¼ 4
xy ¼ 2 sin 2ϕ
r 2r (5.35)
I0 a 2 I0 a
Hy ¼ 4 ðx y Þ ¼ 2 cos 2ϕ:
2
2r 2r
On the midplane (y = 0), H is positive, along the y direction, and falls off with distance
like 1=x2 .
Example 5.4: field from cos ϕ current distribution using the sheet theorem
Assume again that we have a circular sheet with radius a. The current elements are
located at
z ¼ aeiϕ
dz ¼ izdϕ;
so we have
dI dI dϕ I0
¼ ¼ i cos ϕ
dz dϕ dz z
iϕ
I0 e þ eiϕ
¼ i
z 2
I0 z a
¼ i þ :
z a z
Using the current sheet theorem, Equation 5.31,
I0 I0 a
H1 ðzÞ H2 ðzÞ ¼ i þ :
2a 2z2
The field inside the sheet Hin must be finite at z ¼ 0 and for current in the positive
z direction in the first quadrant of the circle, the field must go in the negative
y direction. Therefore we identify H2 with Hin and get
5.6 Green’s theorems in the complex plane 123
I0
Hin ðzÞ ¼ Hin;x þ iHin;y ¼ i :
2a
Equating real and imaginary parts, we find that
Hin;x ¼ 0
I0 (5.36)
Hin;y ¼
2a
in agreement with Equation 5.34. We identify the field exterior to the current sheet
Hext with H1 in the sheet theorem.
I0 a
Hext;x i Hext;y ¼ i
2z2
and
þ ðð
∂P ∂Q ∂P ∂Q
F dz ¼ i þ þi þ dx dy:
∂x ∂y ∂y ∂x
After substitution, we obtain the second complex Green’s theorem.[1]
ð þ
∂F 1
dS ¼ F dz : (5.39)
∂z 2i
a singularity for the case when zo is inside C, as shown in Figure 5.5. We can isolate
the singularity by constructing a small circular contour C1 around it. Then we can
use Green’s theorem in the region between the two contours to transform the
surface integral into a contour integration. However, this requires evaluation on
both contours C and C1 . Halbach proposed adding a constant term to the function
F in Green’s theorem that is (a) analytic in R and (b) makes the contour integration
around C1 vanish.[1] He assumed that F could be written as the product of two
functions F1 and F2 with the property
∂F ∂F2 ðz Þ
¼ F 1 ðzÞ ; (5.42)
∂z ∂z
where F1 contains the singularity. Then F must have the form
Fðz; z Þ ¼ F1 ðzÞ½F2 ðz Þ F2 ðz o Þ: (5.43)
and
∂F iσ 1
¼ ;
∂z 2π z zo
126 Complex analysis of transverse fields
which is the integrand from Equation 5.41. Applying Green’s theorem, we find that
þ
1 iσ 1
H ¼ ðz z o Þdz;
2i 2π z zo
which simplifies to [1, 10]
þ
σ z z o
H ¼ dz: (5.44)
4π z zo
Let us confirm that Equation 5.44 does indeed vanish for the circular contour C1.
Let
z zo ¼ reiθ :
Then for the contour C1,
þ ð 2π
z z o
dz ¼ ir eiθ dθ ¼ 0:
z zo 0
Thus we can ignore the contours around isolated singularities inside the conductor
region and only evaluate Equation 5.44 on the outer boundary of the conductor.
Other quantities of interest can also be conveniently expressed in terms of
contour integrals. For example, the area A of a current block is given by [10, 11]
þ
1
A¼ z dz: (5.45)
2i
Expressions have also been derived for the stored energy.[1, 12]
z ¼ aeiϕ
zo ¼ reiθ :
Case 1: zo inside the conductor 127
where
2ðπ
dϕ
I1 ¼
aeiϕ reiθ
0
(5.47)
ð
2π
eiϕ
I2 ¼ dϕ:
aeiϕ reiθ
0
The logarithm term cancels because it has the same value at 0 and 2π. Thus we
find that
2π
I1 ¼ :
reiθ
9
GR 2.313.1.
Case 2: zo outside the conductor 129
σa2
H ¼ i
2reiθ
σa2
¼ i ½cos θ isin θ:
2r
Equating real and imaginary parts, we find the field outside the conductor is
σa2
Hx ¼ sin θ
2r
(5.50)
σa2
Hy ¼ cos θ:
2r
The field of an elliptical block conductor has also been found using similar
methods.[12, 13]
where
Δz ¼ znþ1 zn : (5.52)
Since the rectangle has four sides and has to close, we identify z5 ¼ z1. For the side
beginning with vertex 1, we define
Δz aeiθ
β1 ¼ ¼ iθ ¼ e2iθ ;
Δz ae
130 Complex analysis of transverse fields
where a is the length of side 1. The second equation comes from considering z1 as the
origin of a line to z2 in polar coordinates. Note that β is a constant because it is defined
in terms of fixed z locations. In a rectangle, the angles of the other sides with respect
to the x axis increase by 90° at each of the vertices. Thus for the side from vertex 2 to
vertex 3,
Similarly we find,
β3 ¼ β1
β4 ¼ β1 :
Then Equation 5.44 gives
ð z2
σ ½z 1 þ β1 ðz z1 Þ z o
H ¼ dz þ
4π z1 z zo
with similar expressions for the remaining three sides. Defining the constant
αn ¼ z n βn zn z o ;
we can write
σ
H ¼ ½α1 I1 þ β1 I2 þ : (5.53)
4π
For points zo outside the contour,
ð z2
dz
I1 ¼
z1 z zo
¼ ½lnðz zo Þzz21
Case 2: zo outside the conductor 131
and10
ð z2
z
I2 ¼ dz
z1 z zo
¼ ½z þ zo lnðz zo Þzz21 :
Writing out the third term for all four sides gives
z4 z3 ¼ ðz2 z1 Þ;
so this term cancels. Thus we find the field at zo due to the rectangular conductor
block is [5, 10]
σ X4
H ¼ hn
4π n¼1
(5.54)
znþ1 zo
hn ¼ ½ðzn zo Þ βn ðzn zo Þ ln :
zn zo
z ¼ reiϕ
zo ¼ 0
and the contour in Equation 5.44 can be broken into the four parts.
10
GR 2.112.1.
132 Complex analysis of transverse fields
(ð ð ϕ2
σ r2
reiϕ1 iϕ1 r2 eiϕ
H ð0Þ ¼ e dr þ r2 i eiϕ dϕ
4π r1 reiϕ1 ϕ1 r2 eiϕ
ð r1 ð ϕ1 )
reiϕ2 iϕ2 r1 eiϕ
þ e dr þ r1 i eiϕ dϕ :
r2 reiϕ2 ϕ2 r1 eiϕ
σ
¼ ðr2 r1 Þðeiϕ2 eiϕ1 Þ:
2π
Expanding the exponentials, we find that the field of the annular sector conductor is
σ
H ð0Þ ¼ ðr2 r1 Þ½ðcos ϕ2 cos ϕ1 Þ iðsin ϕ2 sin ϕ1 Þ; (5.55)
2π
which agrees with Equation 4.50. Note that the field strength is proportional to the
radial thickness.
I 1
HI ðzo Þ ¼ i : (5.56)
2π zI zo
The location of the image filament is
zI ¼ d þ ðd xÞ þ i y
¼ 2d z :
For a uniform distribution of current, we have
ð
σ dS
HI ðzo Þ ¼ i :
2π 2d z zo
To convert this surface integral to a contour integral, we use the complex Green’s
theorem, Equation 5.39. Choosing
iσ z
F¼ ;
2π 2d z zo
we find that the contribution to the field from the image current in a planar iron
surface is
þ
σ z
HI ðzo Þ ¼ dz : (5.57)
4π 2d z zo
We are also interested in the image currents near a circular iron surface at radius
R, as shown in Figure 5.9. From Chapter 2, we know the image current is in the
same direction as the conductor current. If ρ is the distance of the conductor
filament from the center of the circle, then the image filament is a distance R2 =ρ
from the center. Thus we have
z ¼ ρeiϕ
R2 iϕ R2
zI ¼ e ¼ :
ρ z
For a sheet conductor, the contribution of the images in a circular iron cavity to the
field is a sum over the corresponding image currents. Using Equation 5.27, we
obtain
ð
i KðzÞ
HI ðzo Þ ¼ dz
2π zI zo
ð (5.58)
i Kz
¼ dz;
2π R2 zo z
where KðzÞ ¼ dI=dz.
Example 5.8: image field for cos ϕ current sheet in an iron cavity
Let us examine the image field at the origin for a closed circular sheet with radius
a and a cos ϕ angular current distribution. When zo = 0, the integrand does not have
a singularity. Applying Equation 5.58,
ð 2π
i I0 cos ϕ aeiϕ
HI ð0Þ ¼ izdϕ
2π 0 iz R2
ð
iI0 a 2π cos ϕ
¼ dϕ
2πR2 0 eiϕ
ð
iI0 a 2π eiϕ þ eiϕ
¼ dϕ:
4πR2 0 eiϕ
After the integration, we find the contribution of the image field at the origin is
I0 a
HI ð0Þ ¼ i ¼ HIx iHIy ;
2R2
which gives the field components
HIx ¼ 0
I0 a (5.59)
HIy ¼ :
2R2
The contribution of the image field is in the same direction that we saw in
Equation 5.34 for the field of the conductor. The enhancement of the field due to
the presence of the iron is
H ð0Þ þ HI ð0Þ
Eð0Þ ¼
H ð0Þ
(5.60)
a2
¼ 1þ 2:
R
This shows that the iron cavity can contribute up to a factor of 2 to the field at the
origin.
5.10 Multipole expansion 135
For a block conductor with constant current density, the image current in circular
iron is
ð
iσ z
HI ðzo Þ ¼ dS:
2π R2 zo z
Using Equation 5.39 for Green’s theorem and choosing
iσ zz
F¼ ;
2π R2 zo z
we find that the field due to the image current in the circular iron is [1]
þ
σ zz
HI ðzo Þ ¼ dz : (5.61)
4π R2 zo z
Example 5.9: on-axis image field in circular iron for annular sector conductor
Consider an annular sector conductor extending from radius r1 to r2 inside an iron
cavity of radius R. The image field at the origin is given by Equation 5.61 with zo ¼ 0.
þ
σ
HI ð0Þ ¼ zz dz
4πR2
( ð r2 ð ϕ2 ð r1 ð ϕ1 )
σ
¼ e iϕ1 iϕ
r dr ir2 e dϕ þ e
2 3 iϕ2
r dr ir1 e dϕ :
2 3 iϕ
4πR2 r1 ϕ1 r2 ϕ2
Evaluating the integrals and simplifying gives the field contribution due to the iron.
σ
HI ð0Þ ¼ ðr3 r31 Þ½ðcos ϕ2 cos ϕ1 Þ iðsin ϕ2 sin ϕ1 Þ: (5.62)
6πR2 2
Note that this expression has the same sign and angular dependence as the field from
the conductor given in Equation 5.55. The presence of the iron gives the enhancement
factor at the origin [14]
r22 þ r1 r2 þ r21
Eð0Þ ¼ 1 þ : (5.63)
3R2
The field can also be expressed in terms of the integral in Equation 5.41.
ð
i σ
H ðzo Þ ¼ dS
2π z zo
ð
i σ
¼ dS:
2π z 1 zo
z
Expand the factor in the denominator in a geometric series.
ð
i σ zo zo 2
H ðzo Þ ¼ 1þ þ þ dS
2π z z z
This series converges for observation points inside the magnet aperture up to the
closest conductor. Equating this expression with Equation 5.64 gives
X ð ∞ n1
∞
i σX zo
cn zn1
o ¼ dS
n¼1
2π z n¼1
z
∞ ð n1
i X σ zo
¼ dS:
2π n¼1 z z
The zo factor cancels from both sides of the equation. Then matching term by term,
we find
ð
iσ n
cn ¼ z dS: (5.65)
2π
We can convert this surface integral into a contour integral by using the Green’s
theorem, Equation 5.39, with
iσ z1n
F¼ :
2π 1 n
5.10 Multipole expansion 137
ð r1 ð ϕ1 )
þ ðreiϕ2 Þ1n eiϕ2 dr i ðr1 eiϕ Þ1n r1 eiϕ dϕ :
r2 ϕ2
After performing the integrations and simplifying the algebraic results, we find
that
σ
cn ¼ ðr2n r12n Þðeinϕ2 einϕ1 Þ: (5.68)
2πnð2 nÞ 2
Because of the factor in the denominator, this relation cannot be used when n = 2. For
that case, we return to Equation 5.66 and find
þ
σ
c2 ¼ z1 dz :
4π
138 Complex analysis of transverse fields
For the annular sector, performing the integrals and summing terms, we find the
quadrupole multipole is
σ r2
c2 ¼ ln ðe2iϕ2 e2iϕ1 Þ: (5.69)
4π r1
We can get the n = 1 term from Equation 5.67. The dipole multipole is
σ
c1 ¼ ðr2 r1 Þðeiϕ2 eiϕ1 Þ: (5.70)
2π
Let
d ¼ z2 z1 ¼ jdjeiα
zd ¼ ½ðz1 þ z2 Þ;
where α is the angle between d and the x axis. Substituting, we find
2 3
iI 6
6 d 7
7:
H ðzo Þ ¼
2π 4 25
d
ðzd zo Þ2
4
Recall that the magnetic dipole moment is
m ¼ IA ¼ I l d;
5.11 Field due to a magnetized body 139
where l is a unit distance along the z direction. Define m0 as the magnetic moment
per unit length. Then in the limit as d → 0,
I d→m0
d2
→0
4
and the field at zo due to the doublet at z is
i m0 eiα
H ðzo Þ ¼ :
2π ðz zo Þ2
We now want to express the field in terms of the magnetization M. The direction
of m0 is rotated by π/2 with respect to the direction of d. Let us define β to be the
direction of M with respect to the x axis.
π
β¼α
2
eiα ¼ eiβ eiπ=2 ¼ ieiβ :
Then summing up all the magnetic moments in the magnetized body, we have [18]
ð
1 M
H ðzo Þ ¼ dS: (5.71)
2π ðz zo Þ2
Δzn ¼ znþ1 zn :
Since the triangle is closed, z4 ¼ z1 . The slopes of the sides are
Δz n
βn ¼ ;
Δzn
so we can change the integration variable in Equation 5.72 from z* to z for side
n through the relation
dz ¼ βn dz:
M
H ðzo Þ ¼ fðβ3 β1 Þlnðz1 zo Þ þ ðβ1 β2 Þlnðz2 zo Þ
4πi
þðβ2 β3 Þlnðz3 zo Þg: (5.73)
5.12 Force
The vector force dF on a current filament is
! ! !
dF ¼ I dl B :
If the filament is directed along the z direction, B is in the x-y plane, and so is F.
The force can then be written as the complex variable F ¼ Fx þ iFy , where
dFx ¼ μ0 IHy dz
dFy ¼ μ0 IHx dz:
5.13 Conformal mapping 141
F ¼ μ0 H 2 ;
which gives [1]
þ
iμ0
f ¼ H 2 dz : (5.75)
2
This shows that the transverse force per unit length is proportional to the square of
the magnetic field intensity. Examples of complex force calculations can be found
in references.[1, 19]
half-plane or the interior of the unit circle in the w plane. Once the solution is found
for the problem in the w plane, an inverse mapping z ¼ gðwÞ can be used to obtain
the solution to the original problem. The theory of conformal mapping is a major
subject in its own right. We only have space here to briefly introduce the subject
and present a few examples. Fortunately, approximately half the book by Binns and
Lawrenson is devoted to using conformal mapping in the solution of electric and
magnetic field problems.[20] The interested reader can find many useful examples
there.
The bilinear transformation combines the operations of translation, rotation,
stretching, and inversion.[21]
αz þ β
w¼ ;
γz þ δ
ðw w1 Þðw2 w3 Þ ðz z1 Þðz2 z3 Þ
¼ (5.76)
ðw w3 Þðw2 w1 Þ ðz z3 Þðz2 z1 Þ
This expression can be used to create a transformation that maps three given points
in the z plane to three corresponding points in the w plane. An important bilinear
transformation that maps any point zo in the upper half of the z plane into the
interior of the unit circle in the w plane is given by [22]
iθ0 z zo
w¼e : (5.77)
z z o
The points on the x axis are mapped to the boundary of the circle.
z ¼ i↔w ¼ 0
z ¼ ∞↔w ¼ 1;
v1 þ iðu1 þ 1Þ
w2 ¼ u2 þ iv2 ¼
v1 þ iðu21 þ v21 þ u1 Þ
u1 þ iv1
¼ :
u21 þ v21
Representing w1 and w2 in polar coordinates, we find that
1
r2 ¼
r1
θ2 ¼ θ1 ;
which agrees with the result from the method of images.
Suppose that the boundary of some region in the z plane is made up of a series of
straight line segments, as shown in Figure 5.12. The line segments meet at the
vertices z1 ; z2 , . . . It is possible to map this boundary to the real axis in the w plane
by using the Schwarz-Christoffel transformation,[23] which takes the form of the
differential equation
dz
¼ Gðw u1 Þα1 =π1 ðw u2 Þα2 =π1 ðw un Þαn =π1 ; (5.78)
dw
where G is a complex constant and the αi are the interior angles. The points
u1 ; u2 , . . . on the real axis of the w plane correspond to the vertices in the
z plane. The interior of the figure in the z plane maps to the upper half of the
w plane.
Example 5.13: potential of a line current near the corner of two perpendicular
planes
Consider a line current near the perpendicular intersection of two infinitely permeable
plane surfaces, as shown in Figure 5.13. We solve the problem by using the Schwarz-
Christoffel transformation, which in this case takes the form
5.13 Conformal mapping 145
Figure 5.13 Line current near a corner. ABC and abc lie on infinitely permeable
boundary surfaces. The line current is at z1 .
dz
¼ Gðw u1 Þ1=2 ;
dw
since the vertex angle α1 ¼ π=2. Integrating this equation, we find
pffiffiffiffiffiffiffiffiffiffiffiffiffiffi
z ¼ 2G w u1 þ H;
z2 2H z þ H 2
w u1 ¼ :
4G2
Breaking this equation into real and imaginary parts, leads to
x2 y2 2H x þ H 2 þ 2iðx y H yÞ
u u1 þ iv ¼ : (5.79)
4G2
We choose three points A, B, C on the boundary in the z plane and demand that they
correspond to three points a, b, c along the real axis in the w plane according to the
following prescription:
A : x ¼ 0; y ¼ 1↔a : u ¼ 1; v ¼ 0
B : x ¼ 0; y ¼ 0↔b : u ¼ 0; v ¼ 0
C : x ¼ 1; y ¼ 0↔c : u ¼ 1; v ¼ 0
Applying these constraints to Equation 5.79, we find that
u1 ¼ 0
H¼0
4G2 ¼ 1
and the resulting transformation equation is
w ¼ z2 :
146 Complex analysis of transverse fields
In the w plane, the potential is due to the line current w1 and the image current due to
the plane boundary of the infinitely permeable material, as shown in Figure 5.14.
The potential is given by
μ0 I
Wðwo Þ ¼ ½lnðwo w1 Þ þ lnðwo w 1 Þ:
2π
Transforming back to the z plane, we have
μ0 I
Wðzo Þ ¼ ½lnðz2o z21 Þ þ lnðz2o z 2
1 Þ
2π
μ I
¼ 0 ln½ðzo z1 Þðzo þ z1 Þ þ ln½ðzo z 1 Þðzo þ z 1 Þ
2π
μ I
¼ 0 lnðzo z1 Þ þ lnðzo þ z1 Þ þ lnðzo z 1 Þ þ lnðzo þ z 1 Þ :
2π
This shows that the potential in the z plane is due to the physical line current at z1
together with three image currents,[24] as shown in Figure 5.15. The four currents lie
on a circle centered at the corner of the iron surfaces.
5.14 Integrated potentials 147
We can likewise define Aðr; θ; sÞ as the vector potential describing the same three-
dimensional magnetic configuration. Since the conductors have a finite extent in
the s direction, As also vanishes as s→ ±∞. Define
ð∞
Fðr; θÞ ¼ As ðr; θ; sÞ ds:
∞
so that
ð∞
∂F
Bθ ds ¼ : (5.83)
∞ ∂r
Equating the expressions for the integrated values of Br and Bθ , we find that
1 ∂F ∂E
¼
r ∂θ ∂r
∂F 1 ∂E
¼ :
∂r r ∂θ
These two equations have the same form as the Cauchy-Riemann equations in polar
coordinates. Thus F and E represent the real and imaginary parts of the analytic
potential function
WðzÞ ¼ Fðr; θÞ þ iEðr; θÞ:
This potential can be used to describe the influence of a magnet end on the field
quality of a long magnet.[25]
References
[1] K. Halbach, Fields and first order perturbation effects in two dimensional conductor
dominated magnets, Nuc. Instr. Meth. 78:185, 1970.
[2] K. Miller, Introduction to Advanced Complex Calculus, Dover, 1970, p. 66.
[3] M. Spiegel, Complex Variables, Schaum’s Outline Series, McGraw-Hill, 1964, p. 7,
69–70, 83.
[4] R. Beth, Complex methods for three-dimensional magnetic fields, IEEE Trans. Nuc.
Sci. 18:901, 1971.
[5] R. Beth, Complex representation and computation of two-dimensional magnetic
fields, J. Appl. Phys. 37:2568, 1966.
[6] R. Beth, Elliptical and circular current sheets to produce a prescribed internal field,
IEEE Trans. Nuc. Sci. 14:386, 1967.
[7] R. Beth, Some extensions of complex methods for two-dimensional fields, Proc. 6th
Int. Conf. on High Energy Accelerators, Cambridge Electron Accelerator Lab,
Cambridge, MA, 1967, p. 387.
[8] F. Toral, et al., Further developments on Beth’s current sheet theorem: computation of
magnetic field, energy and mechanical stresses in the cross section of particle accel-
erator magnets, IEEE Trans. Appl. Superconductivity 14:1886, 2004.
[9] M. Protter & C. Morrey, College Calculus with Analytic Geometry, Addison-Wesley,
1964, p. 465–468.
[10] R. Beth, Evaluation of current produced two-dimensional magnetic fields, J. Appl.
Phys. 40:4782, 1969.
[11] M. Spiegel, op. cit., p. 114.
[12] R. Beth, Analytic design of superconducting multipolar magnets, Proc. 1968 Summer
Study on Superconducting Devices and Accelerators, Brookhaven National
Laboratory, p. 843–859.
References 149
[13] R. Beth, An integral formula for two-dimensional fields, J. Appl. Phys. 38:4689, 1967.
[14] J. Blewett, Iron shielding for air core magnets, Proc. 1968 Summer Study on
Superconducting Devices and Accelerators, Brookhaven National Laboratory,
p. 1042–1051.
[15] A. Jain, Basic theory of magnets, in S. Turner (ed.), CERN Accelerator School on
Measurement and Alignment of Accelerator and Detector Magnets, CERN 98–05,
1998, p. 1.
[16] K. Binns & P. Lawrenson, Analysis and Computation of Electric and Magnetic Field
Problems, 2nd ed., Pergamon Press, 1973, p. 48–49.
[17] M. Green, Modeling the behavior of oriented permanent magnet material using
current doublet theory, IEEE Trans. Mag. 24:1528, 1988.
[18] K. Halbach, Design of permanent multipole magnets with oriented rare earth cobalt
material, Nuc. Instr. Meth. 169:1, 1980.
[19] R. Beth, Currents and coil forces as contour integrals in two dimensional magnetic
fields, J. Appl. Phys. 40:2445, 1969.
[20] K. Binns & P. Lawrenson, op. cit., chapters 6–10 and appendix III.
[21] K. Miller, op. cit., p. 202–204.
[22] M. Spiegel, op. cit., p. 203–204, 216.
[23] J. Dettman, Applied Complex Variables, Dover, 1984, p. 260–265.
[24] K. Binns & P. Lawrenson, op. cit., p. 43.
[25] F. Mills & G. Morgan, A flux theorem for the design of magnet coil ends, Part. Acc.
5:227, 1973.
6
Iron-dominant transverse fields
In the previous two chapters, we have been discussing transverse field magnets
where the field shape is controlled by the distribution of the conductors. When
iron was present, it served mainly to enhance the strength of the field produced
by the conductors. In this chapter, we examine transverse field magnets where
the primary roles of the conductor and the iron are reversed. Here the shape of
the field is determined by the shape of the iron surface and the conductors are
used to excite the field in the iron.[1, 2] In addition, the iron reduces the
reluctance in the magnetic circuit, allowing a larger useful field for a given
number of amp-turns from the conductor. These types of magnets typically
have a maximum field less than 2 T, so that iron saturation effects do not
destroy the field quality. We will mainly be concerned with the calculation of
the magnetic fields and do not consider the many engineering considerations
necessary to actually build magnets of this type.
150
6.1 Ideal multipole magnets 151
WðzÞ ¼ cn zn
¼ cn rn einθ ;
where the constant cn gives the strength of the potential. The magnetic field for the
ideal multipole is given by
dW
B ¼ i
dz
¼ i n cn zn1 :
The simplest example of an ideal multipole is the dipole (n = 1), which has the
complex potential
W ¼ c1 z
¼ c1 ðx þ i yÞ:
Figure 6.1 illustrates the iron surface and the lines of magnetic field for a dipole
magnet. The dipole has two poles with opposite polarity. Taking the real and
imaginary parts of W, the vector and scalar potentials are
A z ¼ c1 x
μ0 Vm ¼ c1 y:
where the constant h identifies a particular surface. The vector potential is given by
the equipotential
c1 x ¼ k;
152 Iron-dominant transverse fields
Figure 6.2 Iron surface of an ideal quadrupole, ignoring the asymptotic tails.
where the constant k identifies a particular equipotential line. The magnetic field for
the dipole is
∂Az
By ¼ ¼ c1 ;
∂x
which is a constant.
The iron surfaces for higher order ideal multipoles can be found in a similar
manner. The ideal quadrupole (n = 2) has
W ¼ c2 z2 ¼ c2 ðx2 y2 þ 2 i x yÞ
Az ¼ c2 ðx2 y2 Þ
μ0 Vm ¼ 2c2 x y:
r2 sin 2θ ¼ a2 ;
where a is the radius to the center of a pole and θ is measured from the positive
x axis. There are four poles around the perimeter of the magnet that alternate in
polarity. The surface hyperbola for a normal quadrupole has asymptotes along the
x and y axes. The field components inside the aperture are
Bx ¼ 2c2 y
By ¼ 2c2 x:
6.2 Approximate multipole configurations 153
The vertical field on the midplane varies linearly across the aperture and is an
example of a gradient field.
The symmetry properties of a magnet are determined by the symmetry of the
poles. In an ideal normal 2N-multipole, the poles are located at the azimuthal
angles
π
ϕk ¼ ð2k 1Þ ; k ¼ 1; 2; . . . ; 2N (6.1)
2N
The polarities of the poles alternate in direction. The spacing between the poles
is π/N.
m ¼ N ð2n þ 1Þ; n ¼ 0; 1; 2; . . .
Halbach has described methods for determining the effects on the multipole
coefficients of iron saturation and perturbations in the fabrication or construction of
iron-dominated magnets.[3, 4] These methods involve determining the effect of the
perturbation on the scalar potential associated with the pole surface. Among the
effects he considers are azimuthal and radial displacements of the poles and
modifications in the shape of the pole surface. For example, the addition of an
iron shim with thickness profile hðϕÞ modifies the unperturbed scalar potential
approximately by
B0 h B0 Liron
¼ þ ; (6.2)
μ0 μ
where NI is the number of amp-turns in the coil, B0 is the field on the midplane at
the center of the aperture, h is the gap between the iron boundaries, and Liron is the
path length in the iron. Since μ μ0 and the typical path length in the iron is at
most a few times greater than the gap, we have
h Liron
:
μ0 μ
Figure 6.4 Cross-section through a dipole coil. The aperture of the magnet is
located at negative x in this figure.
For finite permeability, the iron will undergo saturation at high field strengths,
the second term in Equation 6.2 may no longer be negligible, and the field in the gap
will be smaller than indicated by Equation 6.3. The field is lower at the center of the
magnet compared to the field at the edges. This creates a small positive sextupole
component in the field. Increasing the width of the pole beyond the useful aperture
can improve the field quality.
Figure 6.4 shows a cross-section through one of the coils. The force on the
conductor is
ð
! ! !
F ¼ J B dV:
and using
dV ¼ hL dx;
B20 hL
¼ :
2 μ0
156 Iron-dominant transverse fields
50 50
40 40
30 30
20 20
10 10
0 0
0 10 20 30 40 50 60 70 80 90 100
Figure 6.5 POISSON model of a window frame dipole. The dimensions are in
centimeters.
F B20 h
¼
L 2μ0
B20
P¼ :
2μ0
1
This is a model of the 18D72 bending magnet that was built at Brookhaven National Laboratory.
2
We will discuss the POISSON program in more detail in Chapter 11.
6.3 Dipole configurations 157
2.8
2.7
2.6 B0 = 2.588 T
2.5
BY [ T ]
2.4
2.3
2.2
2.1 B0 = 2.065 T
0 5 10 15 20
x [cm]
Figure 6.6 Magnetic field along the midplane aperture for the window frame
dipole.
at any of the iron grid points. This quantity depends on the field strength in the iron
and becomes smaller as the iron saturates at higher fields.
The third column shows the half-width of the good field region, defined here as
the distance at which the field exceeds B0 by more than 103 T. The last column
shows the fractional contribution of the sextupole compared to the dipole contribu-
tion to the field. The strength of the field across the half-aperture is shown in
Figure 6.6 for the two higher values of B0 . The field is smallest at the center of the
aperture and grows as it approaches the conductor.
Another common dipole configuration is the H-dipole,[2, 6, 7] shown in
Figure 6.7. The coils are recessed and hidden from direct view of the useful part
of the magnet aperture. This makes the field less sensitive to errors in the coil
location. The field is not as uniform as that in the window frame configuration.
The error multipoles in a fixed, useful aperture decrease exponentially with
increasing pole width. The iron near the edge of the pole is the first area that
exhibits saturation. The good field region can also be extended by adding or
subtracting material at the outer edges of the pole, rounding the corners, or by
158 Iron-dominant transverse fields
tapering the side edge of the iron. Leakage flux, which circulates around the
conductors, can also cause saturation in the iron near the coils, so the magnetic
field in the pole piece is maximum at the center of the pole. This creates a negative
sextupole component in the field in the aperture.
The effective width of the field is ≈w þ g. This causes the flux from the midplane
to get squeezed going into the poles,
wBpole ≈ ðw þ gÞB0 ;
where B0 is the field on the midplane at the center of the aperture. The field on the
pole is then
g
Bpole ≈ 1 þ B0 :
w
The H-dipole has good field quality and mechanical stability. Simple racetrack-
shaped coils can be used to excite the field in the iron.
The C-dipole, shown in Figure 6.8, is a configuration that allows good
access to the magnet aperture from the side.[2, 6] The field in the iron can
be excited with simple racetrack coils. However, the requirement for accessi-
bility leads to a number of disadvantages. The necessary volume of the iron
yoke is larger than for an H-dipole. There may be considerable leakage flux
surrounding the conductors. The asymmetry in the yoke makes the mechanical
stability worse. At high field levels, the attractive force between the poles can
be quite large. Shims may be required at the edges of the poles to get
acceptable field quality. There is a nonuniform magnetic field across the
aperture, although this may be a desirable feature for applications that require
a gradient field component. The field is smaller on the outside side of the gap
than on the inside. The lack of left-right symmetry allows even harmonics to
also be present in the field between the pole pieces. The fringe fields between
the pole pieces extend outward by about a gap length on both sides of the pole
pieces.
6.4 Quadrupole configurations 159
The effective length of a dipole, taking into account its end windings, is
ð
1 ∞
Leff ¼ BðzÞ dz
B0 ∞ (6.4)
≃ Liron þ h;
where B0 is the dipole strength in the center of the magnet. The quantity h is
a length proportional to the aperture of the dipole, which takes into account the
fringe field extending beyond the iron.
Conformal mapping techniques have been used to improve the modelling of the
fringe field from dipole magnets.[8] Maps were used to transform the field from
single and double-sided pole pieces with a uniform gap to the upper half of the
complex plane.
r2 cos 2θ ¼ a2 ;
160 Iron-dominant transverse fields
Figure 6.9 Cutoff angle θ1 for the iron surface in a quadrupole magnet.
where a is the radius to the center of the pole. The cutoff angle θ1 can be selected to
be ~27° in order to eliminate the first allowed error multipole B6. In this case, the
radius to the cutoff point is given by
r1
≃ 1:12:
a
The choice of the angle θ1 also determines the amount of space available for the
conductor. Excitation of the iron poles by the conductor can be determined by using
the Ampère law around the path shown in Figure 6.10.
ða
BðrÞ B0 Liron
NI ≃ dr þ :
0 μ0 μ
The contribution from the path in the iron may be neglected since μr 1.
The contribution for the path along the x axis vanishes because the field is
perpendicular to the path. Thus we have,
ða
gr
NI ≃ dr
0 μ0
ga2
¼ ;
2μ0
6.4 Quadrupole configurations 161
where NI is the amp-turns around a pole and g is the quadrupole field gradient. Thus
the gradient is given by
2μ0 NI
g¼ (6.5)
a2
and the pole tip field is
Bpole ¼ ga
2μ NI (6.6)
¼ 0 :
a
Saturation in the iron affects the area near the conductor first.
Figure 6.11 shows a POISSON model3 of a quadrupole with hyperbolic pole
pieces. The figure shows 1/8 of the symmetric cross-section. The 45° boundary
splits one of the four poles in half. The region in the vicinity of the origin is the open
aperture. The rectangular box on the side of the pole for x ~14 to 26 cm is the
conductor, which wraps around the pole and returns on the opposite side of the
symmetric half pole piece. The pole piece is part of the iron yoke that provides
a flux path to the symmetric adjacent pole.
High-field dipoles and quadrupoles require pole piece materials with a large
value for the saturation magnetic flux density. A number of soft magnetic materials
with large Bsat are listed in Table 6.2. Also listed are the initial and peak values for
the permeability and the coercivity. The resistivity of the material is important for
considerations of eddy current losses in time-varying operations.
Quadrupoles have also been constructed by approximating the hyperbolic sur-
face with a circular cylinder.[6] Consider a circle tangent to the hyperbola at the
center of the pole, as shown in Figure 6.12. The circular surface is continued out to
a cut-off angle θ1 with respect to the center of the pole and then extends outward
along a radius. Let a be the shortest distance from the center of the magnet to the
pole and R be the radius of curvature of the circle, which is centered at C. Then
R ¼ R sin θ1 þ a sin θ1
a sin θ1
¼ :
1 sin θ1
3
This model is an example file that is part of the POISSON code distribution.
162 Iron-dominant transverse fields
0 5 10 15 20 25 30 35
30
25
20
15
10
0
0 5 10 15 20 25 30 35
One solution to these equations, which makes the first allowed multipole harmonic
B6 ¼ 0; is
θ1 ¼ 31:5
r1
¼ 1:785
a
R
¼ 1:094:
a
Conformal mapping techniques have been used to simplify the design of
quadrupoles and higher order multipole magnets.[10] The desired higher-order
multipole is mapped to a dipole geometry, where it is easier to understand what
effects proposed modifications make to the field in the useful aperture. It is also
possible to numerically determine the field quality in the higher multipole
aperture more accurately by computing the multipole coefficients in the trans-
formed geometry.
References
[1] J. Tanabe, Iron-Dominated Electromagnets, World Scientific, 2005.
[2] T. Zickler, Basic design and engineering of normal-conducting, iron-dominated
electromagnets, in Proc. CERN Accelerator School on Magnets, Bruges, Belgium,
CERN-2010-004, June 2009, p. 65.
[3] K. Halbach, First order perturbation effects in iron dominated two dimensional
symmetrical magnets, Nuc. Instr. Meth. 74:147, 1969.
[4] K. Halbach, Fields and first order perturbation effects in two dimensional conductor
dominated magnets, Nuc. Instr. Meth. 78:185, 1970.
[5] A. Banford, The Transport of Charged Particle Beams, Spon Ltd., 1966, p. 81–84,
100–103.
[6] G. Parzen, Magnetic Fields for Transporting Charged Beams, Brookhaven National
Laboratory Report 50536, 1976.
164 Iron-dominant transverse fields
In this chapter, we consider field configurations that have an axial field component.
In straight channels, these fields are azimuthally symmetric around the system axis
and only have axial and radial components. The basic example of this type of
configuration is the closed circular current loop. Combinations of current loops can
be used to produce desired axial field profiles. The current loop can also be
extended axially to generate an ideal sheet solenoid. We conclude the chapter
with a discussion of bent solenoids. When the bent channel forms a closed ring,
we obtain the toroid configuration.
μ0 Ia2
Bz ¼ ; (7.1)
2fa2 þ z2o g3=2
where I is the current in the loop and zo is the distance of the observation point along
the z axis from the plane of the loop. We now consider the determination of the
vector potential in the case when the observation point P is not restricted to lie
along the z axis, as shown in Figure 7.1. We define a coordinate system where the
x axis lies directly below the observation point P. By symmetry, the vector potential
only has a ϕ component and cannot depend on the azimuthal coordinate ϕ.
The vector potential is given by
þ
μ0 I ds
Aϕ ðρ; zÞ ¼ :
4π R
165
166 Axial field configurations
An arbitrary element of current has the Cartesian coordinates (a cos ϕ; a sin ϕ; 0), so
n o1=2
R ¼ ðρ a cos ϕÞ2 þ ða sin ϕÞ2 þ z2
1=2
¼ ρ2 þ a2 2aρ cos ϕ þ z2 :
The contribution from a current element at ϕ makes the same contribution to the
vector potential as the element at –ϕ. In addition, the contribution of each of these
elements to the vector potential at P is proportional to cos ϕ. Therefore we can write
ð
μ0 I π cos ϕ
Aϕ ðρ; zÞ ¼ a dϕ:
2π 0 fρ þ a 2aρ cos ϕ þ z2 g1=2
2 2
Define
4aρ
k2 ¼ : (7.2)
ða þ ρÞ2 þ z2
Then we have
ð π=2
μ Ia k 2sin 2 θ 1
Aϕ ðρ; zÞ ¼ 0 pffiffiffiffiffiffiffiffi dθ
π 4aρ 0 f1 k2 sin 2 θg
1=2
(7.3)
μ Ia k
¼ 0 pffiffiffiffiffiffiffiffi ½2 I1 I2 ;
π 4aρ
where1
ð π=2
sin 2 θ
I1 ¼ 1=2
dθ
0 f1 k2 sin 2 θg
KðkÞ EðkÞ
¼
k2
and2
ð π=2
1
I2 ¼ 1=2
dθ
0 f1 k2 sin 2 θg
¼ KðkÞ:
The function KðkÞ is the complete elliptic integral of the first kind and EðkÞ is the
complete elliptic integral of the second kind.3 Substituting back into Equation 7.3,
we find
μ0 Ia k 2 2
Aϕ ðρ; zÞ ¼ pffiffiffiffiffiffiffiffi 1 KðkÞ 2 EðkÞ ;
π 4aρ k2 k
1 2 3
GR 8.112.5. GR 8.112.1. Properties of complete elliptic integrals are discussed in Appendix F.
168 Axial field configurations
In order to evaluate these field components, we need the derivatives of the para-
meter k defined in Equation 7.2. We have [1]
∂k zk3
¼
∂z 4aρ
(7.5)
∂k k k3 k3
¼ :
∂ρ 2ρ 4ρ 4a
We also need the derivatives4 of the complete elliptic integrals KðkÞ and EðkÞ with
respect to their parameter k.
∂K E K
¼
∂k k ð1 k2 Þ k
(7.6)
∂E E K
¼ :
∂k k k
Evaluating the derivatives together with a lot of algebra,5 we find that [1, 2]
" #
μ0 I z a2 þ ρ2 þ z2
Bρ ¼ KðkÞ þ EðkÞ (7.7)
2π ρfða þ ρÞ2 þ z2 g1=2 ða ρÞ2 þ z2
and
" #
μ0 I 1 a2 ρ2 z2
Bz ¼ KðkÞ þ EðkÞ : (7.8)
2π fða þ ρÞ2 þ z2 g1=2 ða ρÞ2 þ z2
In the limit ρ → 0, k2 ¼ 0 and the elliptic integrals in Equation 7.8 equal π/2. Then
it is straightforward to show that Bz approaches Equation 7.1 for the axial field on
the axis. Using l’Hopital’s rule and the series expansions
4 5
GR 8.123.2,4. A computer algebra program is really useful here!
7.1 Circular current loop 169
π π 2
EðkÞ ≃ k þ
2 8
π π
KðkÞ ≃ þ k2 þ ;
2 8
it is also possible to show that Equation 7.7 for Bρ approaches 0 on the axis, as it
should.
In the preceding derivation, we have gone through a standard approach of
calculating the vector potential and taking its derivatives to find the field compo-
nents. We have done this to illustrate several useful mathematical properties
involving the use of elliptic integrals. We should note, however, that it is possible
in this case to solve the Biot-Savart equation for the fields directly since the
required integrals are known.[3]
Besides the solution given here in terms of KðkÞ and EðkÞ and cylindrical
coordinates, the problem of the circular current loop has been solved using
a number of alternative methods. The vector potential for the circular loop can be
written in terms of Bessel functions as [4]
ð∞
μ Ia
Aϕ ðρ; zÞ ¼ 0 J1 ðkaÞ J1 ðkρÞ ekj z j dk: (7.9)
2 0
For some applications, it is more convenient to solve for the vector potential of the
current loop in spherical coordinates. Spherical solutions for the vector potential
and field have also been given in terms of elliptic integrals.[5, 6] However in
spherical coordinates, it is sometimes more natural to expand the solutions in
Legendre functions. The vector potential for the current loop for r < a is given in
this case as [1]
μ0 I X
∞
sin α r n 1
Aϕ ðr; θÞ ¼ Pn ðcos αÞ P1n ðcos θÞ; (7.10)
2 n¼1 nðn þ 1Þ a
where P1n is an associated Legendre function and α is the polar angle of the loop. For
r > a, the radial factor in this equation must be replaced with
anþ1
:
r
This type of expansion makes it easier to show that the field of the current
loop approaches that for a magnetic dipole in the limit when r >> a. It is
also possible to solve for the field components directly from Maxwell’s
equations.[7]
170 Axial field configurations
4ab
k2 ¼ :
ða þ bÞ2 þ d2
Since the value of Aϕ is constant for all the points on loop 2, the mutual inductance
is [8]
pffiffiffiffiffi 2 2
M ¼ μ0 ab k KðkÞ EðkÞ : (7.11)
k k
Ω ¼ μ 0 Vm
that does not depend on the coordinate ϕ and can be written as the power series
X
∞
Ωðρ; zÞ ¼ cn ðzÞ ρn :
n¼0
We demand that Ω satisfy Laplace’s equation in the region from the axis up to the
location of the closest coil
r2 Ω ¼ 0;
Substituting the series for Ω into Laplace’s equation and bringing the operator
inside the summation sign, we get
X
n ∂ cn
2
c nρ
n n
2 n2
þρ ¼ 0:
∂z2
In order to satisfy this relation, we need to get cancellations between the cn terms of
order n and second derivative terms two orders higher than n. Therefore let us
demand that
∂ 2 cn
cnþ2 ðn þ 2Þ2 ρn ¼ ρn :
∂z2
Making the substitution n→n 2, we can write the coefficient in terms of the
recursion relation
1 ∂2 cn2
cn ¼ ; n ≥ 2: (7.12)
n2 ∂z2
We know that the radial component of the magnetic field has to vanish on the axis
of the system. Since
∂Ω
Bρ ¼ ;
∂ρ
172 Axial field configurations
we find that
Therefore, we must have c1 ¼ 0, and since Equation 7.12 relates c1 to all the
higher odd terms, the series expansion for Ω can only contain even terms. Thus Ω
has the form
X
∞
Ωðρ; zÞ ¼ c2n ðzÞ ρ2n ;
n¼0
where
∂2n c0
ð1Þn
c2n ðzÞ ¼ ∂z2n : (7.13)
2 ðn!Þ2
2n
The numerical factors in the coefficient can be checked by comparing the values
from Equation 7.13 with the values from recursively using Equation 7.12. Define
the magnetic field on the system axis as
∂Ω ∂c0
B0 ðzÞ ¼ Bz ð0; zÞ ¼ ¼ :
∂z ρ¼0 ∂z
X∞
ð1Þn ∂2n B0 2n
Bz ðρ; zÞ ¼ 2 ∂z2n
ρ (7.14)
n¼0 2 ðn!Þ
2n
In cases involving loops and solenoids, where the on-axis fields are known
analytically, it is possible using this method to achieve high accuracy in computing
the field out to radial distances ~70% of the coil radius.[9]
X
∞
V¼ cn rn Pn ðuÞ;
n¼0
c0 ¼ Vðz0 Þ
1
cn ¼ Bc ; n > 0:
μ0 n rcn1 n1
6
Some important properties of Legendre functions are discussed in Appendix D.
174 Axial field configurations
The unknown quantities are now contained in the coefficients Bcn , which are called the
source terms for the central region. These quantities will be related later to the fields
produced by various current elements, such as circular loops. An important feature of
the zonal harmonic method is that the source terms for a given zo only depend on the
coordinates of the current source. Thus once the source terms have been calculated,
they can be used repeatedly in the series expansions for different field points.
With these definitions, V in the central region is written as [10]
rc X∞
Bcn1 r n
V ¼ Vðzo Þ Pn ðuÞ:
μ0 n¼1 n rc
In order to compute the cylindrical components of the magnetic field, we first make
use of the following relations from Figure 7.3.
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
r ¼ ρ2 þ ðz zo Þ2
∂z r ¼ u
∂ρ r ¼ sin θ
and
z zo
u¼
r
1 u2
∂z u ¼
r
u
∂ρ u ¼ sin θ:
r
7
GR 8.914.2.
7.3 Zonal harmonic expansions 175
to write
Changing the index to m ¼ n 1, the axial field in the central region can be written as
X∞ m
r
Bz ¼ Bmc
Pm ðuÞ: (7.17)
m¼0
rc
Following a similar procedure, the transverse field component in the central region is
X
∞
Bc
Bρ ¼ rc n1
½rn1 sin θ ðuP0n þ nPn Þ:
n¼1
n rnc
and shifting the series index again, we find the transverse field in the central region is
X ∞ m
Bcm r
Bρ ¼ sin θ P0m ðuÞ: (7.18)
m¼0
m þ 1 rc
rr X∞
Brnþ1 rr nþ1
V ¼ V0 ð∞Þ þ Pn ðuÞ:
μ0 n¼1 n þ 1 r
Using the recursion relation Equation 7.16 and shifting the series index, we find the
axial field component in the remote region is
X
∞ r mþ1
r
Bz ¼ Brm Pm ðuÞ: (7.19)
m¼2
r
176 Axial field configurations
and shifting the series index, we find that the transverse field component in the
remote region is
X∞
Brm rr mþ1 0
Bρ ¼ sin θ Pm ðuÞ: (7.20)
m¼2
m r
Now that we have determined the series representations of the field components
in the central and remote regions, we turn to the calculation of the source terms.
Consider a field point on the z axis in the central region with z > zo . In this case,
θ = 0, ρ = 0, u = 1, and r ¼ z zo . From Equation 7.17, we have
X
∞
Bc
B0 ðzÞ ¼ Bz ð0; zÞ ¼ n
ðz zo Þn :
n¼0
rnc
Bc0 ¼ Bz ð0; zo Þ
is the axial magnetic field at the source point. The Taylor series around the point
z0 is
X∞
1 dn B0 n
B0 ðzÞ ¼ n z
ðz zo Þ : (7.21)
n¼0
n! dz o
rnc dn B0
Bcn ¼ ðzo Þ: (7.22)
n! dzn
Suppose that the current source is a circular loop, as shown in Figure 7.4. For
a current loop, we have rc ¼ rr ¼ rL . From Equation 7.1, the axial field at the field
point F is
μ0 Iρ2L
B0 ðzÞ ¼ : (7.23)
2d3
7.3 Zonal harmonic expansions 177
1 X
∞
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ hn Pn ðuÞ;
1 þ h 2hu n¼0
2
where h < 1. Differentiating both sides of this equation with respect to u, we find
1 X
∞
3=2
¼ hn P0nþ1 ðuÞ: (7.24)
f1 þ h2 2hug n¼0
so that
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
d ¼ rL 1 þ h2 2 huL :
1 1X ∞
¼ hn P0nþ1 ðuL Þ:
d3 r3L n¼0
8
GR 8.921.
178 Axial field configurations
Substituting this into Equation 7.23 and comparing with the general series expan-
sion Equation 7.17 for the case of a field point on the axis, we can conclude that the
source term for the circular loop in the central region is [11]
μ0 Iρ2L 0
Bcn ¼ P ðuL Þ: (7.25)
2r3L nþ1
Note that this expression does not depend on any parameters of a field point.
We follow an analogous procedure in the remote region to find that
μ0 Iρ2L 0
Brn ¼ P ðuL Þ: (7.26)
2r3L n1
μ0 I a 2
Bz ð0; zÞ ¼ FðzÞ; (7.27)
2
where
1 1
FðzÞ ¼ : (7.28)
2 3=2 3=2
fa2 þ ðb zÞ g fa2 þ ðb þ zÞ2 g
In the Helmholtz configuration, the spacing 2b between the coils equals the radius
a of the coils. The field of a Helmholtz pair at the origin is
μ0 I 8
Bz ð0; 0Þ ¼ pffiffiffi
a 5 5
(7.29)
μI
≃ 0:7155 0 :
a
7.4 Multiple coil configurations 179
The axial field between the coils falls off slowly with z. At the center of the current
loops the field is
a μ I 1 1
Bz 0; ¼ 0
1 þ pffiffiffi
2 a 2 8
μ0 I
≃ 0:6768 :
a
In the vicinity of the origin, the axial field can be expanded in the Taylor series,
Equation 7.21. The field uniformity is determined by the leading-order terms in this
expansion. Because of the symmetry of the coil arrangement, all the odd power
terms in the series have to vanish. The virtue of the Helmholtz configuration is that
the second derivative term in the expansion also vanishes. Thus the leading
correction in the Taylor series is the fourth order term, which is proportional to
∂ 4 Bz μI
≃ 19:8 05 :
∂z 4 a
Thus in the vicinity of the origin, the axial field is [14]
z 4
μ0 I
Bð0; zÞ ≃ 0:7155 0:825 þ : (7.30)
a a
The field at any point off the axis can be found by using the elliptic integral
solutions for the current loop given in Equations 7.7 and 7.8.[15] In the plane
midway between the coils, the field only has an axial component because of
symmetry. Defining the scaled radius u = ρ/a, the elliptic integral parameter is
16u
k2 ¼
4u2 þ 8u þ 5
180 Axial field configurations
By numerically evaluating this expression, the field is found to be quite uniform in the
vicinity of the axis.[16] It falls off to 99.93% of the on-axis value at a radius of 0.2a.
The Helmholtz pair arrangement has the geometric property that the two coils lie
on the surface of a sphere, as illustrated in Figure 7.6. For this configuration, we have
tan θ ¼ 2
a
sin θ ¼ ;
r
so the radius of the sphere is
a
r¼ ≃ 1:118 a:
sinðtan1 2Þ
The Helmholtz pair also has interesting asymptotic behavior.[17] Expanding the
on-axis field in powers of a/z, the field at large distance is given by
μ0 I a2 3μ0 I a2
Bz ð0; zÞ ≃ þ ð4b2 a2 Þ þ
z3 2z5
The leading term is the magnetic dipole term. However, the next term in the series
vanishes under the Helmholtz condition a = 2b.
An inverse Helmholtz pair has the currents in the two coils flowing in opposite
directions. In this case, the field at the origin vanishes, and the leading multipole
term is the field gradient. The optimum gradient for fixed radius a is [18]
7.4 Multiple coil configurations 181
dBz 48 μ I
¼ pffiffiffi 02
dz 25 5 a
μI
≃ 0:8587 02 :
a
If practical constraints demand it, other gradient optimizations are possible for
fixed spacing b or for fixed radius r2 ¼ a2 þ b2 .[19]
The classic design using three coils is the Maxwell tricoil, shown in Figure 7.7.
The tricoil has a pair of identical coils and a third coil with larger radius in the
symmetry plane between the coil pair. The three coils all lie on the surface of
a sphere. This design improves on the field quality from the Helmholtz pair by also
making the fourth-order term in the Taylor series vanish. Thus the first correction
term is sixth-order. Maxwell’s solution is
rffiffiffi
4
a¼ R
7
rffiffiffi
3
b¼ R
7
49
I¼ I0 :
64
The magnetic field at the origin is
μ0 I
Bz ð0; 0Þ ¼ 60 :
R
An improved three coil design with three circular coils of the same radius has
a larger uniform field region than Maxwell’s design.[15] Another sixth-order
182 Axial field configurations
design used three square coils.[20] An exhaustive study of multi-coil systems using
the methods of zonal harmonics has identified many uniform and gradient field
configurations with higher-order error corrections.[21]
where n is the number of turns per unit length. Performing the integration9 gives
" #
μ0 nI z þ L=2 z L=2
Bz ð0; zÞ ¼ : (7.31)
2 fa2 þ ðz þ L=2Þ2 g1=2 fa2 þ ðz L=2Þ2 g1=2
9
GR 2.264.5.
7.5 Sheet model for the solenoid 183
L
Bz ð0; 0Þ ¼ μ0 n I pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi : (7.33)
L2 þ 4a2
In the limit of an infinitely long current sheet, this expression reduces to the field of
an ideal solenoid, Equation 1.28.
Bz ¼ μ0 n I:
There is a close connection between the derivatives of the on-axis solenoid field
and the on-axis fields of the current loops at the ends of the solenoid.[9, 10] In the
coordinate system in Figure 7.8,
dBSolenoid ð0; zÞ L L
z
¼ n Bz
Loop
0; z þ Bz
Loop
0; z (7.34)
dz 2 2
The off-axis expansion method discussed in Section 7.2 can be used in conjunction
with Equation 7.34 to find the field of a sheet solenoid.[9]
We turn next to calculating the field of a solenoid at any point, including points
off the symmetry axis. We will perform a direct calculation of the field using the
Biot-Savart equation
! !
! μ0 I dl R
dB ¼ :
4π R3
Consider the solenoid geometry shown in Figure 7.9. Because of the azimuthal
symmetry of the current, the field is also azimuthally symmetric. Thus for
mathematical simplicity, we are free to choose the field point F to lie directly above
the x axis. The distance from the source current element to the field point is
!
R ¼ ðρ a cos ϕ0 Þ ^x a sin ϕ0 ^y þ ðz z0 Þ ^z ; (7.35)
and we define
Define the distances from the observation point to the two ends of the solenoids in
terms of the new variables
L
z1 ¼ z
2
L
z2 ¼ z:
2
After doing the integration, we get10
z2
I1 ¼ Ωðz1 Þ;
ðe z2 Þfa2 þ ρ2 2aρ cos ϕ0 þ z22 g
1=2
10
GR 2.264.5.
7.5 Sheet model for the solenoid 185
where we use the symbol Ω here as a shorthand notation that means a second term
similar to the first, but with L replaced by –L, i.e., the other end of the solenoid.
Then we have
ðπ
μ I 0a a ρ cos ϕ0 z2 0
Bz ¼ 0 0 1=2 dϕ Ωðz1 Þ:
2π 0 ða2 þ ρ 2 2aρ cos ϕ Þ a þ ρ þ z2 2aρ cos ϕ
2 2 2 0
cos ϕ0 ¼ 1 þ 2 x2
2 (7.38)
dϕ0 ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffi dx:
1 x2
This gives
ð1
μ0 I 0 a a þ ρ 2ρx2 z2
Bz ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffi dx Ωðz1 Þ:
π 0
2
½ða þ ρÞ 4aρx2 fða þ ρÞ2 þ z2 4aρx2 g1=2 1 x2
2
and
4 aρ
n¼ : (7.40)
ða þ ρÞ2
We find that
μ0 I 0 a z2
Bz ¼ 2
½ða þ ρÞ I2 2ρ I3 Ωðz1 Þ;
πða þ ρÞ fða þ ρÞ2 þ z2 g1=2
2
where11
ð1
dx
I2 ¼
0 ð1 n x2 Þfð1 k2 x2 Þ ð1 x2 Þg1=2
¼ Πðk; nÞ:
The function Πðk; nÞ is the complete elliptic integral of the third kind. The other
integral is
11
GR 8.111.4.
186 Axial field configurations
ð1
x2
I3 ¼ dx:
0 ð1 n x2 Þfð1 k2 x2 Þ ð1 x2 Þg1=2
μ0 I 0 a z2 h aρ i
Bz ¼ KðkÞ þ ðΠðk; nÞ KðkÞÞ Ωðz1 Þ:
π ða þ ρÞfða þ ρÞ2 þ z2 g1=2 2a
2
(7.42)
zL
e
¼ 2
( 2 )1=2 ΩðLÞ;
L
ðe z2 Þ a2 þ ρ2 2 aρcos ϕ0 þ z
2
12 13
GR 8.111.2 and 8.112.1. GR 2.264.6.
7.5 Sheet model for the solenoid 187
where e was defined in Equation 7.37. Substituting these back into the equation for
Bρ and simplifying, we get
ðπ
μ I 0a cos ϕ0
Bρ ¼ 0 dϕ0 Ωðz1 Þ:
fa2 þ ρ2 þ z22 2aρ cos ϕ0 g
2π 1=2
0
KðkÞ EðkÞ
¼ :
k2
Substituting and simplifying, we find the transverse component of the solenoid
field is [22]
1=2
μ I 0 fða þ ρÞ2 þ z22 g
Bρ ¼ 0 ½2ðKðkÞ EðkÞÞ k2 KðkÞ Ωðz1 Þ: (7.43)
4π ρ
Equations 7.42 and 7.43 are exact solutions for the sheet solenoid that are valid
for all points in space, except for observation points at the same radius as the
current sheet. For an arbitrary field point with the cylindrical coordinates (ρ, ϕ, z),
we can write the Cartesian values of the transverse field as
x
Bx ¼ Bρ cos ϕ ¼ Bρ
ρ
(7.44)
y
By ¼ Bρ sin ϕ ¼ Bρ :
ρ
An alternate solution for the field components has also been given in terms of
related functions known as generalized complete elliptic integrals.[23]
The vector potential for the sheet solenoid may also be expressed exactly in
terms of elliptic integrals [22, 24] as
14
GR 3.153.5.
188 Axial field configurations
"
μ I 0 z2 1=2
Aϕ ¼ 0 fða þ ρÞ2 þ z22 g ðKðkÞ EðkÞÞ
4π ρ
# (7.45)
ða ρÞ2
1=2
ðΠðk; nÞ KðkÞÞ Ωðz1 Þ;
fða þ ρÞ2 þ z22 g
where k and n are given by Equations 7.39 and 7.40, respectively. The vector
potential and field for the sheet solenoid may also be expressed as sums of Bessel-
Laplace integrals [25, 26] or in terms of modified Bessel functions.[27]
Example 7.3: mutual inductance and axial force between a solenoid and a loop
Once the vector potential and the field components for the circular current loop and
the sheet solenoid are known, the mutual inductance between a solenoid and a current
loop can be computed as
ð
ΦL 1
MðS; LÞ ¼ ¼ Aϕ ðSÞ dl
IS IS
2πρL
¼ Aϕ ðSÞ;
IS
where the symbols S and L refer to the solenoid and the loop.
7.6 Block model for the solenoid 189
Assume the current is flowing in the same direction in the loop and in the solenoid.
The axial force acting on a current loop due to the magnetic field from the solenoid is
ð ð 2π
! !
Fz ðS; LÞ ¼ IL dlL B ðSÞ ¼ IL ρL Bρ ðSÞ dϕ
0
¼ 2πIL ρL Bρ ðSÞ;
where the minus sign indicates that the force tries to pull the loop and the solenoid
together.
where again Ω is used as a shorthand for the expression in the first term with
L replaced by –L. Performing the integral,15 we find for points along the symmetry
axis
( " 1=2 #)
μ0 J b þ fb2 þ ðz þ L=2Þ2 g
Bz ð0; zÞ ¼ ðz þ L=2Þ ln 1=2
ΩðLÞ: (7.46)
2 a þ fa2 þ ðz þ L=2Þ2 g
The off-axis field from the block solenoid is usually treated by summing over the
fields from a set of current sheets, using one of the methods we have previously
discussed. For example, the block conductor may be simulated using a radial
distribution of current sheets expressed in elliptic integrals.[22] It is also possible
to express the thick solenoid field in terms of a radial expansion of the on-axis field,
[9] as a series of zonal harmonics,[10, 11, 21] or in terms of Bessel-Laplace integrals.
15
GR 2.271.5.
190 Axial field configurations
[25, 26, 28] The good field region calculated for a solenoid from a properly designed
block conductor is frequently larger than that for a sheet solenoid.[21]
The flux leaving a solenoid travels outside and returns through the opposite end.
As a result, the fringe field on the outside of the solenoid can be quite significant.
If the fringe field is unacceptable, it can be reduced by adding supplemental bucking
coils or by using iron shielding. Figure 7.12 shows a POISSON model for the
magnetic field in a typical solenoid. The figure shows 1/4 of a plane projection
through the solenoid. The vertical axis is the centerline of the solenoid. The field is
symmetric on both sides of the vertical axis and both sides of the horizontal axis.
The program used Dirichlet boundary conditions on the left, right, and top borders,
and Neumann boundary conditions on the bottom border. The figure on the left
shows the field from just the coil, while the figure on the right illustrates the reduction
in the exterior field from adding a cylindrical iron return yoke. The current in the coil
was the same for both figures.
(a) (b)
180 180
160 160
140 140
120 120
100 100
80 80
60 60
40 40
20 20
0 0
0 20 40 60 80 100 0 20 40 60 80 100
Figure 7.12 Magnetic field of a solenoid coil (left); field for a coil surrounded by
an iron return path (right).
where (h1 ; h2 ; h3 ) are a set of scale factors. In the Frenet-Serret coordinate system
considered here,[30, 31] the reference curve is a circle and the origin of the unit
vectors moves along the circle. The unit vector s is in the bending plane and tangent
to the circle. The unit vector r is in the bending plane and perpendicular to s.
The unit vector v is perpendicular to the bending plane. The curvilinear scale
factors are
192 Axial field configurations
hr ¼ hv ¼ 1
r
hs ¼ 1 þ ;
ρ
We also include terms in the potential that allow for the possibility of superimposed
transverse fields. The gradient of V is defined as
1
rV ¼ ∂r V ^r þ ∂s V ^s þ ∂v V ^v :
hs
Thus the magnetic field components are
V01 ¼ B1
V10 ¼ A1
V11 ¼ B2
2V20 ¼ A2 :
Br ≃ A1 A2 r þ B 2 v
Bv ≃ B1 þ B2 r 2 V02 v
(7.48)
1
Bs ≃ ðbs A0 1 r þ B0 1 vÞ;
hs
16
In the case where a charged particle has to follow the reference path in the horizontal plane, we must have the
horizontal dipole field A1(s) = 0.
7.8 Toroid 193
1 1
∂r ðhs Br Þ þ ∂s Bs þ ∂v Bv ¼ 0:
hs hs
If no superimposed transverse fields are present, the first-order axial field in the
bent channel is
r
Bs ≃ bs bs : (7.50)
ρ
It is also possible to define an expansion for the field that makes use of “curved
multipoles” that directly correspond to the solution of Laplace’s equation in the
curved coordinate system.[32]
7.8 Toroid
When the bent solenoid channel is extended to form a closed circular ring, we have
a toroid, as shown in Figure 7.14. The current loops from the solenoid are centered
on the circular system axis and the plane of the loops lie in the ρ-z plane.
^ depend on the azimuthal location around
The direction of the unit vectors ^ρ and ϕ
the toroid. Since the coils are closer together on the side nearer to the center of
194 Axial field configurations
curvature, we expect that the field inside the bent solenoid will have a gradient with
respect to the system axis, in agreement with Equation 7.50.
Because of the symmetry of the configuration, all of the field components must
be independent of the azimuthal angle ϕ. Let the mean radius of the toroid equal
b and the radius of the current loops equal a. Then a simple application of
the Ampère law on the midplane (z = 0) shows that Bϕ ¼ 0 for ρ < b a and for
ρ > b þ a since no net current is enclosed in a circular path in those regions.
However, applying the Ampère law on a circular path on the midplane, we find
the field inside the toroid is
μ0 N I
Bϕ ¼ ; (7.51)
2πρ
where N is the number of conductor turns around the circumference and ρ is the
radius of the path. This shows that the field varies like 1=ρ inside the toroid.
A cross-section of the toroid at some azimuthal angle ϕ is shown in Figure 7.15a.
The angle α gives the location of an element of the current loop. The current
element has the Cartesian coordinates
xl ¼ ðb þ a cos αÞ cos ϕ
yl ¼ ðb þ a cos αÞ sin ϕ
zl ¼ a sin α
(a) (b)
Figure 7.15 (a) Cross-section of the toroid at the azimuthal location ϕ; (b) location
of the field observation point F.
The location of the field observation point (ρ, 0, z) is shown in Figure 7.15b. The
resulting distance vector is
!
R ¼ ½ρ ðb þ a cos αÞcos ϕ ^x ðb þ a cos αÞsin ϕ ^y þ ðz a sin αÞ ^z :
Applying the Biot-Savart law to any point inside the toroid shows that Bz and Bρ
vanish. It follows that the field inside the toroid has to have the form
B ¼ Bϕ ðρ; zÞ:
The analytic results for Bϕ are complicated [33] expressions defined in terms of
integrals of elliptic integrals. Alternatively, one could examine the field inside the
toroid by evaluating one of the integrals in terms of complete elliptic integrals and
performing the other integral numerically.
References
[1] W. Smythe, Static and Dynamic Electricity, 2nd ed., McGraw-Hill, 1950, p. 270–275.
[2] J. Stratton, Electromagnetic Theory, McGraw-Hill, 1941, p. 263.
[3] R. Good, Elliptic integrals, the forgotten functions, Eur. J. Phys. 22:119, 2001.
[4] W. Panofsky & M. Phillips, Classical Electricity and Magnetism, 2nd ed., Addison-
Wesley, 1962, p. 156.
[5] R. Schill, General relation for the vector magnetic field of a circular current loop:
a closer look, IEEE Trans. Magnetics 39:961, 2003.
[6] J. D. Jackson, Classical Electrodynamics, Wiley, 1962, p. 142.
[7] D. Redzic, The magnetic field of a static current loop: a new derivation, Eur. J. Phys.
27:N9, 2006.
[8] G. Harnwell, Principles of Electricity and Magnetism, 2nd ed., McGraw-Hill, 1949,
p. 329.
196 Axial field configurations
In the previous chapters, we considered magnets that for the most part had
continuous transverse or longitudinal fields along some system axis. In this chapter,
we look instead at magnetic channels where the on-axis field is periodic. Periodic
field configurations are used for focusing charged particle beams and for produc-
tion of radiation at light sources. We begin by considering the field produced by
helical conductor windings. Then we examine several examples where we demand
some desired field configuration along the axis and then find off-axis field compo-
nents that satisfy Maxwell’s equations.
since z progresses by a distance λ as the azimuthal angle goes once around the
circumference. With this definition, α = 0 corresponds to the limiting case when the
helix reduces to a circular loop. It is convenient to write λ as a function of a.
197
198 Periodic magnetic channels
λ ¼ 2πa tan α
The azimuthal angle and axial distance are connected through the helix constraint
z ¼ a ϕ tan α: (8.1)
Consider the cross-section through the helix shown in Figure 8.2. A single
conductor follows a helical path in a current sheet of radius a. Let the observation
point F have cylindrical coordinates (ρ, ϕ, z) and the current element at the
location S on the conductor have coordinates (a; ϕ0 ; z0 ). The current element is
given by
!
dl ¼ a sin ϕ0 dϕ0 ^x þ a cos ϕ0 dϕ0 ^y þ a tan ϕ0 dϕ0 ^z
Taking into account the constraint between z0 and ϕ0 , the direct evaluation of
B using the Biot-Savart equation only requires an integration over ϕ0 : The
8.1 Field from a helical conductor 199
integration limits for a winding of finite length can be found using Equation 8.1.
Although the resulting integral for the general case is complicated, the solution for
observation points along the axis of the helix is fairly straightforward.[1]
The vector potential for an infinitely long helix is given by
ð !
! μ K ða; ϕ0 ; z0 Þ 0
A ðρ; ϕ; zÞ ¼ 0 dS : (8.2)
4π R
The sheet current density only has components in the ϕ0 and z0 directions. The pitch
angle α for the helical winding can be written as
Kz 0
tan α ¼ ;
K ϕ0
Kϕ0 ¼ ka Kz0 :
! I
K ða; ϕ0 ; z0 Þ ¼ ð^z 0 þ kaϕ^0 Þ δðϕ0 kz0 εÞ; (8.3)
a
where ε is the azimuthal angle of the winding at z0 ¼ 0. The Dirac delta function
enforces the constraint between changes in z0 and changes in ϕ0 .
We can write the periodic delta function in Equation 8.3 as a Fourier series. Let
τ ¼ ε þ kz0 . Then
X
∞
f ðϕ0 Þ ¼ δðϕ0 τÞ ¼ a0 þ ½an cos nϕ0 þ bn sin nϕ0 :
n¼1
1 1X ∞
δðϕ0 τÞ ¼ þ ½cos nτ cos nϕ0 þ sin nτ sin nϕ0
2π π n¼1
" # (8.4)
1 X∞
0
¼ 1 þ 2 cos nðϕ τÞ :
2π n¼1
The distance from the current element to the observation point can be written as
1=2
R ¼ fa2 þ ρ2 2 aρ cosðϕ ϕ0 Þ þ ðz z0 Þ2 g : (8.5)
We need to express the unit vectors in Equation 8.3 in terms of unit vectors in the
coordinate system for the observation point. The axial unit vectors are identical,
z^0 ¼ ^z . The azimuthal unit vector is given by
Thus there are in general nonvanishing components of the vector potential in all
three dimensions.
Substituting Equations 8.3–8.6 into Equation 8.2, the axial component of the
vector potential is given by
" #
X
ðð 1þ2
∞
cos nðϕ0 kz0 εÞ
μ0 I
a dϕ0 dz0 :
n¼1
Az ðρ; ϕ; zÞ ¼ 1=2
4π 2πa 0 2
fa2 þ ρ2 2 aρ cos ðϕ ϕ Þ þ ðz z0 Þ g
(8.7)
μ0 I μ IX∞
Az ðρ; ϕ; zÞ ¼ ln a þ 0 Kn ðnkaÞ In ðnk ρÞ cos nðϕ kz εÞ (8.8)
2π π n¼1
μ0 I μ IX∞
Az ðρ; ϕ; zÞ ¼ ln ρ þ 0 Kn ðnk ρÞ In ðnkaÞ cos nðϕ kz εÞ : (8.9)
2π π n¼1
μ0 Ika X
∞
Aρ ðρ; ϕ; zÞ ¼ ½Knþ1 ðnkaÞ Inþ1 ðnkρÞ Kn1 ðnkaÞ In1 ðnkρÞ
2π n¼1 (8.11)
sin nðϕ kz εÞ
μ0 Ika X
∞
Aρ ðρ; ϕ; zÞ ¼ ½Knþ1 ðnkρÞ Inþ1 ðnkaÞ Kn1 ðnkρÞ In1 ðnkaÞ
2π n¼1
(8.12)
sin nðϕ kz εÞ :
(8.13)
μ0 I kρ
Aϕ ðρ; ϕ; zÞ ¼
4π
μ0 I ka X
∞
þ ½Knþ1 ðnkaÞ Inþ1 ðnkρÞ þ Kn1 ðnkaÞ In1 ðnkρÞcos nðϕ kz εÞ
2π n¼1
(8.14)
202 Periodic magnetic channels
μ0 I ka2
Aϕ ðρ; ϕ; zÞ ¼
4πρ
μ0 I ka X
∞
þ ½Knþ1 ðnkρÞ Inþ1 ðnkaÞ þ Kn1 ðnkρÞ In1 ðnkaÞcos nðϕ kz εÞ :
2π n¼1
(8.15)
The solution for the magnetic field components can be found by taking the curl
of A. For the case ρ < a the field components are [2, 3, 4]
μ0 Ik2 a X
∞
Bρ ðρ; ϕ; zÞ ¼ nK0 ðnkaÞ In0 ðnkρÞ sinðnðϕ kz εÞÞ
π n¼1 n
μ0 Ika X
∞
In ðnkρÞ
Bϕ ðρ; ϕ; zÞ ¼ nKn0 ðnkaÞ cosðnðϕ kz εÞÞ (8.16)
π n¼1 ρ
μ0 Ik μ0 Ik2 a X
∞
Bz ðρ; ϕ; zÞ ¼ nK0 ðnkaÞ In ðnkρÞ cosðnðϕ kz εÞÞ;
2π π n¼1 n
μ0 Ik2 a X
∞
Bρ ðρ; ϕ; zÞ ¼ nK0 ðnkρÞ In0 ðnkaÞ sinðnðϕ kz εÞÞ
π n¼1 n
μ0 I μ0 Ika X
∞
Kn ðnkρÞ
Bϕ ðρ; ϕ; zÞ ¼ þ nI 0 ðnkaÞ cosðnðϕ kz εÞÞ (8.17)
2πρ π n¼1 n ρ
μ0 Ik2 a X
∞
Bz ðρ; ϕ; zÞ ¼ nI 0 ðnkaÞ Kn ðnkρÞ cosðnðϕ kz εÞÞ:
π n¼1 n
(a)
(b)
Figure 8.3 (a) The dependence of Bz along one period of the helix; (b) the depen-
dence of Bρ (solid) and Bϕ (dashed) versus z. The parameters used here were
λ = 10 cm, a = 10 cm, ρ = 5 cm, ϕ = 0, ε = 0, I ¼ 105 A, and N = 40 terms in the
sums.
μ0 Ik2 a 0
Bρ ð0; 0; zÞ ¼ K1 ðkaÞ
2π
Bϕ ð0; 0; zÞ ¼ 0
μIk
Bz ð0; 0; zÞ ¼ 0 :
2π
204 Periodic magnetic channels
There is a close relationship between these results and our previous results for
the field of a solenoid. In practice, the conductor in a solenoid is wound in many
helical layers over a cylindrical form. The helical pitch length λ for a solenoid is
very small. In the previous chapter, the field for a solenoid was derived by assuming
that the field came from a longitudinal distribution of parallel infinitesimal current
loops. In the limit that k → ∞, it can be shown by taking the asymptotic limits for
the Bessel functions that Equation 8.16 approaches the on-axis field of an infinitely
long solenoid
Bz ð0; 0; zÞ ¼ μ0 nI;
where n = 1/λ is the number of turns per unit length.[6] In the opposite limit where
k → 0, the helical fields approach that of a straight conductor.
In order to satisfy this equation for all x, y, and z, Bx and Bz must have the form
where f and g are unknown factors. Substituting these field expressions back into
the divergence equation gives the constraint
f α þ βB0 þ gγ ¼ 0: (8.19)
γB0
g¼ ;
β
αB0
f ¼ :
β
Substituting these values for f and g into Equation 8.19, we find the wave number
constraint
γ2 ¼ α2 þ β2 : (8.20)
4π2
λ2 ¼ :
β2 α2
For a periodic solution, we have the additional constraint that β > α. The solution
for the planar transverse field is then [7]
αB0
Bx ¼ sinðαxÞ sinhðβyÞ cosðγz φÞ
β
By ¼ B0 cosðαxÞ coshðβyÞ cosðγz φÞ (8.21)
γB0
Bz ¼ cosðαxÞ sinhðβyÞ sinðγz φÞ:
β
The other two solutions consistent with Equation 8.18 have the transverse
dependences
coshðαxÞ cosðβyÞ
206 Periodic magnetic channels
and
coshðαxÞ coshðβyÞ:
These solutions can be derived following the same procedure used above. Each
solution has a unique relation among the wavenumbers.[7]
It’s important to keep in mind that this type of derivation only represents part of
the problem. What we have shown is that our desired on-axis field profile is a valid
solution of Maxwell’s equations. However, what we have not considered is
a configuration of conductors that actually produces that desired field.
The obvious trial solution here would be a series of transverse permanent magnet
or electrically excited magnetic poles that alternate in direction along the system
axis. Oftentimes the field from a realistic coil distribution can only approximate the
desired field. We define the problem of finding a current distribution that produces
a specified magnetic field as an inverse problem to distinguish it from the situation
encountered using the Biot-Savart formula, where we find the magnetic field
produced by a given current distribution. Finding a suitable current distribution
frequently involves using numerical optimization methods.
Bρ0 ¼ B0 cosðγz ϕ εÞ
(8.22)
Bϕ0 ¼ Bz0 ¼ 0;
where ε is an initial phase shift. For points off the axis, we look for a solution for Bρ
with the form
Using Equations 8.23–8.26, we can write the div B = 0 equation directly in terms
of G.
1 1
½ρ ∂2ρ ðρGÞ þ ∂ρ ðρGÞ G γ2 ρ G ¼ 0: (8.27)
ρ ρ
∂2 ðγ ρ GÞ ∂ðγ ρ GÞ
γ2 ρ2 þγρ ½1 þ ðγ ρÞ2 ðγ ρ GÞ ¼ 0: (8.28)
∂ðγ ρÞ 2 ∂ðγ ρÞ
This is the differential equation for the modified Bessel1 function I1 , where the
unknown variable is γ ρG and the argument of the Bessel function is γ ρ. Thus we
have
γ ρ GðρÞ ¼ I1 ðγ ρÞ
I1 ðγ ρÞ
GðρÞ ¼ : (8.29)
γρ
1
Some properties of the modified Bessel functions are described in Appendix C.
208 Periodic magnetic channels
where I0 is the modified Bessel function of order 0. Substituting this back into
Equation 8.23, we find
1
Bρ ¼ C I0 ðγ ρÞ I1 ðγ ρÞ cosðγz ϕ εÞ:
γρ
1 2
I0 ðuÞ ≃ 1 þ u þ (8.31)
4
and4
1 1 3
I1 ðuÞ ≃ uþ u þ (8.32)
2 16
Thus near the axis, we take the leading terms for I0 and I1 and find that
C
Bρ ≃ cosðγz ϕ εÞ:
2
Comparing this with Equation 8.23 gives C ¼ 2B0 . The solution for the helical
transverse field is then [8]
1
Bρ ¼ 2B0 I0 ðγ ρÞ I1 ðγ ρÞ cosðγz ϕ εÞ
γρ
I1 ðγ ρÞ (8.33)
Bϕ ¼ 2B0 sinðγz ϕ εÞ
γρ
Bz ¼ 2B0 I1 ð γ ρÞ sinðγz ϕ εÞ:
Note that the argument for the Bessel functions involves the longitudinal wave-
number γ.
2 3 4
AS 9.6.26. AS 9.6.12. AS 9.6.10.
8.4 Axial fields 209
Thus
γ C
Bρ ¼ B0 I1 ðαρÞ sinðγz εÞ þ :
α ρ
5
AS 11.3.25.
210 Periodic magnetic channels
γ2 ∂I0 ðαρÞ
I1 ðαρÞ ¼ 0:
α ∂ρ
Using6
∂I0 ðαρÞ
¼ α I1 ðαρÞ;
∂ρ
we find that α = γ. Thus the solution for the periodic axial field is
Bρ ¼ B0 I1 ðγρÞ sinðγz εÞ
Bϕ ¼ 0 (8.34)
Bz ¼ B0 I0 ðγρÞ cosðγz εÞ:
References
[1] W. Smythe, Static and Dynamic Electricity, 2nd ed., McGraw-Hill, 1950, p. 276–278.
[2] T. Tominaka, Vector potential for a single helical current conductor, Nuc. Instr. Meth. A.
523:1, 2004.
[3] T. Tominaka, M. Okamura & T. Katayama, Analytic field calculation of helical coils,
Nuc. Instr. Meth. A. 459:398, 2001.
[4] R. Hagel, L. Gong & R. Unbehauen, On the magnetic field of an infinitely long helical
line current, IEEE Trans. Mag. 30:80, 1994.
[5] S. Park, J. Baird, R. Smith & J. Hirshfield, Exact magnetic field of a helical wiggler,
J. Appl. Phys. 53:1320, 1982.
[6] T. Tominaka, Magnetic field calculation of an infinitely long solenoid, Eur. J. Phys.
27:1399, 2006.
[7] D. Sagan, J. Crittenden, D. Rubin & E. Forest, A magnetic field model for wigglers and
undulators, Proc. 2003 Part. Accel. Conf., Portland, OR, p. 1023.
[8] J. Blewett & R. Chasman, Orbits and fields in the helical wiggler, J. Appl. Phys.
48:2692, 1977.
[9] R. Fernow & R. Palmer, Solenoidal ionization cooling lattices, Phys. Rev. Special
Topics – Accel. and Beams 10:064001, 2007.
6
AS 9.6.27.
9
Permanent magnets
Since M is uniform here, the second term vanishes. In addition, the first term
vanishes on the flat end faces. Thus A only depends on the surface contributions
around the sides of the cylinder. We assume the sources for the potential are
azimuthally directed Amperian currents with current density
! ¼! M ^n
Ka
^
¼ M ϕ: (9.1)
Thus the field in a bar magnet is analogous to the field in a solenoid, where the
Amperian currents here take the place of the conduction currents in a solenoid.
The magnetic flux density in the axial direction is then given by Equation 7.32
211
212 Permanent magnets
μ0 nI
Bz ð0; zÞ ¼ ðcos β2 cos β1 Þ;
2
where βi are the angles from the observation point along z to the two outer edges at
the ends of the cylinder. Making the substitution Ka ¼ n I; we get
μ0 K a
B¼ ðcos β2 cos β1 Þ
2
μM
¼ 0 ðcos β2 cos β1 Þ (9.2)
2
BR
¼ ðcos β2 cos β1 Þ;
2
where BR is the remanent field for the permanent magnet. For a point P inside the
magnet, cos β1 < 0 and cos β2 > 0. Thus B points along the positive z direction.
From Gauss’s law for a pillbox on one of the end faces, B must also be directed
along +z outside the magnet. We can find the magnetic intensity H from
μ0 H ¼ B μ0 M:
Thus [2]
Ka
HP ¼ ðcos β2 cos β1 Þ Ka
2
(9.3)
Ka
¼ ðcos β2 cos β1 2Þ:
2
Since the two cosine terms are both smaller than one, this expression is negative.
Thus H inside the magnet points in the opposite direction from M and B.
Now let us consider the situation at the point Q on the end face of the magnet.
In this case cos β2 ¼ 0 and
9.1 Bar magnets 213
μ0 K a
BQ ¼ cos β1 : (9.4)
2
On the inside surface of the end face, we find from Equation 9.3 that
Ka
HQin ¼ ðcos β1 þ 2Þ;
2
which points along −z. On the outer surface of the end face, M = 0 and we find from
Equation 9.4 that
Ka
HQout ¼ cos β1 ;
2
which points along +z. The behavior for H is similar to the case of the electric field
from a surface distribution of charge. It is sometimes convenient when modeling
bar magnets to assume that the end faces contain a distribution of fictitious
magnetic charges or “poles.” In terms of magnetic charges, we can express the
surface and volume charge densities as [3]
!
σm ¼ M ·^n
!
ρm ¼ r·M (9.5)
and the scalar potential as
ð ð
1 σm 1 ρm
Vm ¼ dS þ dV: (9.6)
4π R 4π R
Thus in a bar magnet, we can assume that B comes from the Amperian currents
flowing azimuthally around the cylinder and that H comes from magnetic charges
on the flat end faces. In a real bar magnet, M is typically weaker near the end faces
than it is in the central region. In this case, there will also be contributions to the
field from the volume integrals above.
Consider a bar magnet with radius a that is much smaller than its length L. This is
sometimes referred to as a “magnetic needle.”[4] We look at the field on the axis
outside the magnet at location z. We have
( 2 )1=2
a
cos β2 ¼ 1 þ :
z L=2
For a L,
2
1 a
cos β2 ≃ 1 þ
2 z L=2
214 Permanent magnets
and we have a similar expression for cos β1 with L → –L. Then from Equation 9.2,
μ0 M a 2 μ0 M a2
BðzÞ ¼ 2
:
4ðz L=2Þ 4ðz þ L=2Þ2
We interpret this as the field due to a positive magnetic charge at face 2 and a
negative charge at face 1. The field strength falls off like the inverse square distance
from the charge.
around the circuit since there are no conduction currents. Assuming the gap is small
and that the leakage flux is negligible, this gives
Hg Lg Hm Lm þ Hi Li ¼ 0;
where the subscripts g, m, and i refer to gap, magnet, and iron. The quantity Hm Lm
takes the place of NI in magnetic circuits powered by conductors. Neglecting
leakage, the flux is constant around the loop, so we have
ΦB ¼ Bi Ai ¼ Bg Ag ¼ Bm Am ;
Wg ¼ ½ μ0 B2g Ag Lg
¼ ½ H g L g Bg Ag :
Hm Lm ≃ Hg Lg
Thus the stored energy in the gap is proportional to the BH product and the volume
of the permanent magnet.
Since B and H point in opposite directions, permanent magnets operate in
the second quadrant of the hysteresis curve. To maximize the energy stored in the
gap, it is desirable to operate at a point where the BH product is maximum.
Rewriting Equation 9.8 as
216 Permanent magnets
3.5
2.5
1.5
0.5
0
0 1 2 3 4 5
Bm Lm
¼ ;
H m Am R
the quantities on the right-hand side can be adjusted to get the B-H “load line” for
the magnetic circuit to pass through the point where the BH product is maximum.
Figure 9.3 shows a simple POISSON model1 of a magnetic circuit energized by
a pair of permanent magnets. The figure shows one quarter of a symmetric circuit.
The x and y axes are symmetry planes in this figure. The use of a 2D program such
as POISSON assumes that the configuration is uniform over a large distance in the
third dimension (into the figure).
The permanent magnet material is located in the rectangle near the origin.
The magnetization was oriented in the vertical direction. Considerable field fring-
ing is evident in the air gap for this simple geometry. The amount of fringing
depends strongly on the type of permanent magnet material that is used.[1]
1
The calculation was actually made with the PANDIRA program in the POISSON distribution.
9.4 Model for rare earth materials 217
a high remanent field but is easily demagnetized and has a large leakage flux.[1]
Ferrite contains iron oxides (Fe2 O3 ). It has a low BR, but it is a cheap material.
Ceramic ferrites are compounds of barium or strontium ferrite. Samarium and
neodymium are rare earth elements.[8] They have large values of ðBHÞmax , can
produce a large field from a compact magnet, are very resistant to demagnification,
and have small leakage flux.
BR ≃ μ 0 H C :
218 Permanent magnets
In this case, the fields from assemblies of blocks can be linearly superimposed.
If the magnetization is anisotropic, the relation between B and H can be written as
Bk ¼ μk Hk þ BR
(9.10)
B? ¼ μ ? H ? ;
where k refers to the easy direction in the material and ? refers to the direction
perpendicular to it. For convenience, we define the reluctivity γ = 1/μ. Then from
Equation 9.10, we have
BR
Hk ¼ γk Bk
μk (9.11)
≃ γk Bk HC :
and
! ! ! !
H ¼ γk B k þ γ? B ? H C : (9.13)
The vectors BR and HC are directed in opposite directions along the easy axis.
Taking the divergence of Equation 9.12, we find that
! ! !
r·ðμk H k þ μ? H ? Þ ¼ r· B R ≡ ρm ; (9.14)
which relates H to the density of magnetic charges. Taking the curl of Equation 9.13
in a region with no conduction currents, we get
! ! ! !
r ðγk B k þ γ? B ? Þ ¼ r H C ≡ Jm ; (9.15)
which relates B to the Amperian currents. Thus the material can be treated
!
magnetically as vacuum together with either a charge density r · B R or
!
a current density r H C : For homogeneous materials, these charges or currents
vanish everywhere except on the surface.
The scalar potential for the permanent magnet material can be written as
ð
1 ρm
μ 0 Vm ¼ dV
4π R
ð !
1 ðr· B R Þ
¼ dV; (9.16)
4π R
9.5 Rare earth model in two dimensions 219
where R is the distance between the source point and the field point. If the material
is inhomogeneous, we can define G = 1/R and make use of Equation B.3
! ! !
r·ðG B R Þ ¼ G r· B R þ B R ·rG
which vanishes because we can choose the surface to lie outside the material where
BR ¼ 0. We also have
!
R
rG ¼ 3 :
R
Thus the scalar potential for inhomogeneous material is
ð! !
1 B R· R
μ 0 Vm ¼ dV: (9.17)
4π R3
For homogeneous material, we can use the divergence theorem in Equation 9.16
and find that
ð!
1 B R ·^
n
μ 0 Vm ¼ dS
4π R
! ð
B R n^
¼ · dS; (9.18)
4π R
BR ¼ BRx þ i BRy :
∂x BRy ∂y BRx ¼ μ0 σ:
1
B ðzo Þ ¼ ½I1 I2 ; (9.20)
2πi
where
ðð
∂x BRy
I1 ¼ dx dy
zo x i y
ðð
∂y BRx
I2 ¼ dx dy:
zo x i y
which gives
ð" ð #
BR y BR y
I1 ¼ dx dy
zo x i y ðzo x i yÞ2
þ ðð
BR y BR y
¼ dy dx dy:
zo x i y ðzo x i yÞ2
The first term vanishes for a line integral evaluated outside the material where
BR ¼ 0. Thus we have
ðð
BR y
I1 ¼ dx dy:
ðzo x i yÞ2
9.6 Multipole expansion for continuously distributed material 221
BR 1
F¼ :
2π zo z
The field for a homogeneous block can then be written in terms of the contour
integral
þ
BR dz
B ðzo Þ ¼ : (9.22)
4π i zo z
It is also possible by direct integration to write this in the alternate forms [9]
þ
BR dx
B ðzo Þ ¼
2π i zo z
þ
BR dy
¼ : (9.23)
2π zo z
1 1 X∞
zm
¼
¼ o
:
zo z z 1 zzo m¼0
zmþ1
If we let n = m + 1, then
1 X∞ n1
zo
¼ : (9.24)
zo z n¼1
zn
1 X
∞
n zn1
¼ o
: (9.25)
ðzo zÞ2 n¼1
znþ1
Now let us consider the design of an ideal multipole magnet. We use a polar
coordinate system with z ¼ r eiϕ . We assume the material is located in an
annular region between the radii r1 and r2 . Assume that we have the freedom
to specify the direction of the easy axis everywhere in the material. Let the easy
axis in the material located at angle ϕ be rotated through an angle βðϕÞ with
9.6 Multipole expansion for continuously distributed material 223
respect to its direction at ϕ = 0. The contribution to the field from the material at the
angle ϕ is
The first term in the square brackets comes from the easy axis distribution
and the second term comes from the azimuthal dependence in zn+1 in
Equation 9.25. We can make the multipole as large as possible by making
the quantity in square brackets equal to 0. This determines the constraint on
the easy axis angles
βðϕÞ ¼ ðn þ 1Þ ϕ: (9.30)
where
ð 2π
I¼ ei ðNnÞϕ dϕ
0
2π for n ¼ N
¼
0 for n ≠ N:
B ðzo Þ ¼ BN zN1
o :
separated geometrically by the angle π/4 and the easy axis rotates by π/2 from block
to block.
Let Cn be the contribution to the multipole Bn for some reference block. Then the
contribution of block m to Bn is
2π 2π
Cn exp i m ðN þ 1Þ exp i m ðn þ 1Þ ;
M M
where the first exponential comes from the rotation of the easy axis and the second
exponential comes from the definition of Bn . The sum of all M blocks gives
X
M1
Nn
Bn ¼ C n exp i 2πm ;
m¼0
M
where ν is an index to the set of allowed multipoles. The multipole order that
corresponds to a given ν is given by [9]
n ¼ N þ νM: (9.34)
Figure 9.6 shows the geometry of a trapezoidal block. Only the two segments
marked 1 and 2 contribute to the contour integral. On path 1 we have
y ¼ x tan α;
y ¼ x tan α:
226 Permanent magnets
We can express
e inα
ð1 i tan αÞn ¼ :
cosn α
Writing the exponential term as a sine function, we find that Cn for the trapezoidal
block is
" n1 #
BR r1
Cn ¼ 1 cosn α sin nα: (9.35)
π ðn 1Þr1n1 r2
Substituting this back into Equation 9.33, we find that the field inside the aperture is
given by
" n1 #
X ∞ n1
zo n r1
B ðzo Þ ¼ BR 1 Kn ; (9.36)
ν¼0
r 1 n 1 r2
References 227
References
[1] J. Bahrdt, Permanent Magnets Including Undulators and Wigglers, Proc. CERN
Accelerator School on Magnets, Bruges, Belgium, June 2009, CERN-2010-004,
p. 185.
[2] D. Tomboulian, Electric and Magnetic Fields, Harcourt, Brace & World, 1965, p. 240.
[3] L. Eyges, The Classical Electromagnetic Field, Dover, 1980, p. 139.
228 Permanent magnets
[4] W. Smythe, Static and Dynamic Electricity, 2nd ed., McGraw-Hill, 1950, p. 434.
[5] P. Lorrain & D. Corson, Electromagnetic Fields and Waves, 2nd ed., Freeman, 1970,
p. 409.
[6] High Performance Permanent Magnets, Magnet Sales and Manufacturing, Inc.,
Culver City, CA, 1993.
[7] A. Chao & M. Tigner, Handbook of Accelerator Physics and Engineering, World
Scientific, 1999, p. 366.
[8] J. Becker, Permanent magnets, Sci. Am. 223:92, 1970.
[9] K. Halbach, Design of permanent multipole magnets with oriented rare earth cobalt
material, Nuc. Instr. Meth. 169:1, 1980.
[10] A. Chao & M. Tigner, op. cit., p. 468.
[11] Q. Peng, S. McMurry & J. Coey, Cylindrical permanent magnet structures using
images in an iron shield, IEEE Trans. Mag. 39:1983, 2003.
10
Time-varying fields
Until this point, we have only examined magnetic effects due to steady currents or
magnetic materials in stationary configurations. In this chapter, we will partially
relax this constraint by considering phenomena where there are slow variations in
current or magnetic flux. By slow, we mean slow enough that we can ignore all
effects of electromagnetic radiation. We begin with a discussion of Faraday’s law,
which presents another connection between electric and magnetic phenomena.
This is followed by a more detailed discussion of the energy associated with
a magnetic field, including the energy loss from the hysteresis cycle in ferromag-
netic materials. We find that Faraday’s law leads to the production of eddy currents
in some materials, while the skin effect can restrict currents to a layer near the
surface. We introduce the displacement current, which finally allows us to give
a complete set of Maxwell’s equations for stationary media. We conclude the
chapter with a brief discussion of magnetic measurements.
229
230 Time-varying fields
þ ð
!! ∂ !
ℰ ¼ E ·dl ¼ B ·^n dS; (10.2)
∂t
where E is the electric field intensity and the surface S is bounded by the closed
circuit. The field E acts on a distance element dl in its rest frame. Because of the
tangential boundary condition on E, it follows that C can refer to any closed loop in
space, not just a physical circuit.[1] The line integral1 on the left side of
Equation 10.2 is called the electromotance ε. If the flux links a coil with N turns,
the electromotance must be multiplied by N. Since the contour integral is non-zero,
the induced electric field in this case is nonconservative, i.e., work is done on
a charge going around the contour.
! μ ^
B ðtÞ ¼ 0 IðtÞ ϕ:
2πρ
1
Historically, ε has also been referred to as an emf or electromotive force.
10.1 Faraday’s law 231
Then since the surface S is arbitrary, we find the differential form of Faraday’s law is
!
! ∂B
r E ¼ : (10.3)
∂t
This equation is valid at any point in space. We can relate this to the vector potential by
! ∂ !
r E ¼ ðr A Þ
∂t
!
∂A
¼ r ;
∂t
so that
!!
! ∂A
r E þ ¼ 0:
∂t
Since its curl vanishes, the quantity in parentheses must be the gradient of a scalar
function, which we denote Ve .
!
! ∂A
rVe ¼ E þ :
∂t
L I ¼ N ΦB :
Taking the time derivative of both sides, we find that an alternative definition
of L is
ε
L¼ : (10.5)
dI=dt
across the current element and the source must supply the power
dP ¼ I dVe
!
¼ J dA ðrVe ·dl Þ;
where A is the cross-sectional area of the conductor. Since J and dl are parallel, we
can use Equation 10.4 to write
!!
! ∂A !
dP ¼ E þ · J dτ;
∂t
The first term on the right side is the power used to compensate for energy losses
from heating in the conductor. The second term is the power used to set up the
magnetic field associated with the increasing current, which is the subject of
interest here. If we let W represent the energy stored in the magnetic field, then
10.2 Energy in the magnetic field 233
ð !
dW ∂A !
¼ · J dτ: (10.6)
dt ∂t
into Equation 10.6, we find the energy stored in the magnetic field is
ð
!!
W ¼ ½ J · A dτ: (10.7)
! !
If the current distribution is a current loop, we let J dτ→I dl and Equation 10.7
becomes
ð
!!
W ¼ ½ I A ·dl :
W ¼ ½ I ΦB : (10.8)
Returning again to Equation 10.7, we can use the curl H = J equation to write the
energy as
ð
! !
W ¼ ½ ðr H Þ· A dτ:
Rewriting the first term in terms of B and using the divergence theorem in
the second, we get
2
J.D. Jackson, Classical Electrodynamics, Wiley, 1962, p. 176.
234 Time-varying fields
ð ð
!! ! !
W ¼ ½ H · B dτ þ ½ H A ·^n dS:
Looking at the surface integral, we know that the field from a conductor element falls
off like 1=R2 and the vector potential falls off like 1=R, while the surface area only
grows like R2 . By evaluating at a sufficiently large distance, the second integral
vanishes. Thus the energy stored in the magnetic field from conduction currents is
ð
!!
W ¼ ½ B ·H dτ (10.9)
dW dΦB
¼ NI
dt dt
dB
¼ NIA
dt
NI dB
¼ Al ;
l dt
where N is the number of conductor turns, A is the cross-sectional area of the
sample, l is the mean circumference of the ring, and B is the average flux density
inside the sample. Using the Ampère law, this can be written
dW dB
¼ HV ;
dt dt
where V is the volume of the sample.
Now consider the hysteresis loop shown in Figure 10.2. The energy supplied by
the source in moving from point a to point b along the loop is
10.4 Eddy currents 235
ðb
Wab ¼ V H dB:
a
Since dB is the independent variable, the value of this integral is the area projected
on the B (vertical) axis in the figure. Going from point b to point c along the loop,
I is in the same direction, but is decreasing. Thus the electromotance changes sign
and some energy is returned to the source.
ðc
Wbc ¼ V H dB:
b
The sum of these two integrals is the area inside the hysteresis loop in the first
quadrant. If we continue this analysis for a complete cycle, we find that the net
energy lost in the ferromagnetic material per cycle is [4]
þ
W ¼ V H dB: (10.12)
The range of current densities over which this linear relation holds depends on the
material. Thus we have
! ∂ !
r ðr J e Þ ¼ σ μ ðr H Þ;
∂t
where Je is the eddy current density. We can use the vector identity B.7 on the left
side of this equation and the curl H = J equation on the right side to get
!
! ! ∂J e
rðr· J e Þ r2 J e ¼ σ μ :
∂t
Since the divergence term on the left side vanishes, we find that [6]
!
2! ∂J e
r J e ¼ σμ : (10.14)
∂t
This is a form of the diffusion equation. The rate of build-up of the eddy currents is
controlled by the factor σμ.
If instead, we begin by taking the curl of the Ampère law, we find
! !
r ðr H Þ ¼ r J e
!
¼σr E:
Applying Equation B.7 on the left-hand side and Faraday’s law on the right, we
obtain the equation
!
2! ∂H
r H ¼ σμ : (10.15)
∂t
10.4 Eddy currents 237
Thus the magnetic field associated with the eddy currents also satisfies a diffusion
equation with the same characteristic constant. If one specifies the time depen-
dence for H and the geometry of the configuration, the diffusion equation can be
solved for the spatial and time dependence of the magnetic field due to eddy
currents.[5] This may lead to a series of terms, each with its own characteristic
time dependences.
Example 10.2: time constant for eddy currents in a solid iron core
Consider a long H-dipole with a solid iron yoke, as shown in Figure 10.3. For slow
time changes, eddy currents can flow throughout the volume of the iron yoke
surrounding the coil.[7] The magnetic flux from the eddy currents is not symmetric
with the flux from the coils, which causes the iron saturation to vary with transverse
position.
Consider a path through the iron yoke at the midplane in the region 0 ≤ x ≤ d.
Assuming there is no leakage flux, all of the return flux from the conductor has to
pass across this path. Assume the current in the conductor is changing with time.
Then the magnetic field in the vicinity of the path is in the y direction, the induced
electric field is in the z direction, and on the midplane both fields are only functions
of x and t.
! ¼ B ðx; tÞ ^y
B y
!
E ¼ Ez ðx; tÞ ^z :
The eddy currents flow parallel to the z axis until they reach the magnet end,
where they reverse direction and flow back at the symmetric (x, y) position on
the other side of the magnet. From the Ampère and Ohm’s laws, we have
∂By ðx; tÞ
¼ σ μ Ez ðx; tÞ;
∂x
238 Time-varying fields
and
where p is the variable conjugate to t in the Laplace transform. Taking the derivative
of Equation 10.16 with respect to x and substituting Equation 10.17, we get
Defining k 2 ¼ σμp, the solution for the magnetic field consistent with the boundary
conditions is [7]
The field across the return yoke is asymmetric and is larger on the side nearer the coil.
At the edge of the path closest to the conductor, we have
pffiffiffiffiffiffiffiffiffiffi
kd ¼ d σμp
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
¼ σ μ p d2
pffiffiffiffiffiffi
¼ p τ;
where [7]
τ ¼ σμd2 : (10.18)
The variable τ has the dimensions of time. It gives a characteristic time for eddy
current effects in this configuration. Note that it depends quadratically on the width
d of the return yoke.
Eddy currents can be suppressed by restricting the rate of change of the desired
field or by constructing the magnet in such a way that potential eddy current loops
are minimized. Magnet yokes are frequently constructed by assembling thin iron
laminations for this reason.
10.5 Skin effect 239
Now assume that the current is flowing along a conducting slab that occupies the
space y ≤ 0. Then the component of J flowing in the z direction, for example, is
Jz ¼ Jz0 eiωt :
Applying Equation 10.19 to this, we find
d2 Jz0
¼ ζ 2 Jz0 :
dy2
Jz0 ¼ JS eζ j y j ;
where JS is the spatial dependence of the current density on the surface of the slab.
Using the relation
i ¼ ½ ð1 þ 2 i 1Þ ¼ ½ ð1 þ iÞ2 ;
we can write
1 pffiffiffiffiffiffiffiffiffi
ζ ¼ pffiffiffi ð1 þ iÞ ωσμ:
2
The boundary condition for large | y | eliminates the negative solution for ζ. Thus
rffiffiffiffiffiffiffiffiffi
ωσμ
ζ ¼ ð1 þ iÞ :
2
240 Time-varying fields
We see that the current density decreases exponentially with distance into the
surface. In addition, there is a phase shift of the current flowing inside the material
with respect to the current flowing on the surface. These effects scale with the skin
depth parameter δ. For a copper conductor with current varying at 1 kHz, the skin
depth is ~2.1 mm.
The vector D is called the electric flux density3 and is related to the electric field
intensity by
! !
D ¼ε E (10.24)
for linear materials, where ε is the permittivity. Taking the time derivative of
Equation 10.23, we get
!
∂ρ ∂D
¼ r· : (10.25)
∂t ∂t
3
Historically, the vector D is also known as the electric displacement.
10.7 Rotating coil measurements 241
!
Comparing Equations 10.22 and 10.25, we see that the quantity ∂ D =∂t acts
like an additional kind of current. Thus we define the displacement current
density as
!
! ∂D
Jd¼ : (10.26)
∂t
Taking this into account, the Ampère law must then be modified as [11]
!
! ! ∂D
r H ¼ J þ : (10.27)
∂t
This shows that a magnetic field can also be produced by a time-varying electric
field.
At this point, we can summarize the complete set of Maxwell’s equations for
stationary media in Table 10.1. It is important to keep in mind that writing the
equations in this form assumes the validity of the constitutive relations
! ¼μ!H
B
! !
J ¼σ E
! !
D ¼ε E:
One technique, which relates directly with our previous discussions of the multi-
pole content of fields, involves measurements in long magnets with large aperture
using a rotating coil.[14] The azimuthal component of the field can be measured
using a radial coil, the principle of which is shown in Figure 10.4. The flux through
the wire loop with N turns is
ð r2
ΦB ðθÞ ¼ N L Bθ ðr; θÞ dr
r1
X
∞ ð r2
¼ N L ðAn sin nθ þ Bn cos nθÞ rn1 dr
n¼1 r1
X∞ n
r2 rn1
¼ N L ðAn sin nθ þ Bn cos nθÞ ;
n¼1
n
where we have used Equation 4.8 to express the azimuthal field in terms
of multipole field components. If the coil rotates at a constant rate, we have θ ¼
ω t and
Performing a Fourier analysis on the voltage signal allows the multipole coeffi-
cients to be determined.
10.7 Rotating coil measurements 243
ð 2π X
∞ ð 2π
VðθÞ sin mθ dθ ¼ ωNL ðr2 r1 Þ An cos nθ sin mθ dθ
n n
0 n¼1 0
ð 2π
Bn sin nθ sin mθ dθ ¼ ωNLðrn2 rn1 Þ Bn π δmn :
0
Thus
ð 2π
1
Bm ¼ VðθÞ sin mθ dθ: (10.29)
π ωN L ðrm
2 r1 Þ
m
0
and
ð 2π
1
Am ¼ VðθÞ sin mθ dθ: (10.33)
2πωNLRm sin mδ 0
References
[1] W. Panofsky & M. Phillips, Classical Electricity and Magnetism, 2nd ed., Addison-
Wesley, 1962, p. 159.
[2] P. Lorrain & D. Corson, Electromagnetic Fields and Waves, 2nd ed., Freeman, 1970,
p. 351–356.
[3] J. Stratton, Electromagnetic Theory, McGraw-Hill, 1941, p. 126–128.
[4] P. Lorrain & D. Corson, op. cit., p. 398–400.
[5] G. Moritz, Eddy currents in accelerator magnets, in D. Brandt (ed.), Proc. CERN
Accelerator School on Magnets, Bruges, Belgium, June 2009, CERN-2010-004,
p. 103.
[6] G. Harnwell, Principles of Electricity and Magnetism, 2nd ed., McGraw-Hill, 1949,
p. 340–341.
[7] K. Halbach, Some eddy current effects in solid core magnets, Nuc Instr. Meth.
107:529, 1973.
[8] E. Kreyszig, Advanced Engineering Mathematics, Wiley, 1962, chapter 4.
[9] J. Marion, Classical Electromagnetic Radiation, Academic Press, 1965, p. 149.
[10] P. Lorrain & D. Corson, op. cit., p. 105–106.
[11] J. D. Jackson, Classical Electrodynamics, Wiley, 1962, p. 177–178.
[12] C. Reymond, Magnetic resonance techniques, in S. Turner (ed.), CERN Accelerator
School on Measurement and Alignment of Accelerator and Detector Magnets, CERN
98-05, 1998, p. 219.
[13] J. Kvitkovic, Hall generators, in S. Turner (ed.), ibid., p. 233.
[14] A. Jain, Harmonic coils, in S. Turner (ed.), ibid., p. 175.
11
Numerical methods
Until now, we have mostly considered magnetostatic problems that have analytic
solutions. In practice, this has usually restricted our choice of problems to situa-
tions where one of the problem boundaries coincides with a coordinate system axis
and to solutions that can be written as products of functions of a single coordinate.
Frequently, more complicated problems can only be solved using numerical meth-
ods.[1] A number of commercial and freeware programs are available for solving
magnetic problems. A lot of effort has been devoted to making many of these
programs accurate, efficient, and user-friendly. If such a program is available and
can address the problem under consideration, it is often the best choice to use it.
However, there are occasions when new code must be written to solve a problem.
It is also important to have some basic understanding about the methods involved in
obtaining these solutions. In this chapter, we will examine three numerical methods
that have been used for solving boundary value problems involving the Poisson
equation: finite differences, finite elements, and integral equations. In each case,
the analytical equation or its solution is approximated in some way that leads to
a matrix equation for the unknown potential or field. We conclude the chapter with
a discussion of inverse problems and optimization techniques.
245
246 Numerical methods
in the greatly simplified situation shown in Figure 11.1. The line is discretized into 7
nodes. Let u refer to the unknown quantity, which we assume satisfies the Dirichlet
boundary conditions u0 ¼ 0 and u6 ¼ 0. Assume that the source function f has the
value f3 ¼ v at the center of the line and is 0 otherwise. The values of u at the five
interior nodes are the unknown quantities. Using Equation 11.5, each interior node
satisfies the equation
We can write the equations for the five unknowns in the form of a matrix equation
Cu ¼ g: (11.6)
For a square grid in two dimensions, let us designate the node under consideration
as node 0 and its four nearest neighbors as nodes 1–4, as shown in Figure 11.2.
We can write the Laplacian operator in terms of the values at the five nodes as [2]
1
r2 u ≃ ½uðx h; yÞ þ uðx þ h; yÞ þ uðx; y hÞ þ uðx; y þ hÞ 4uðx; yÞ:
h2
(11.7)
2 u1 u2 u3 u4 1 1
r2 u ≃ þ þ þ þ u0 ; (11.8)
h2 pð p þ rÞ qðq þ sÞ rð p þ rÞ sðq þ sÞ pr qs
where p, q, r, s are dimensionless scaling factors for the spacings from node 0 to its
nearest neighbors. Using this expression for problems where the physical boundaries
of conductors and iron are parallel to the x and y axes, it is possible to set up the
equations for the interior nodes together with nodes coinciding with the boundaries.
Higher-order difference equations for the Laplacian are also possible.[4]
A complication arises in setting up the node equations for nodes adjacent to
boundaries that do not align exactly with the grid spacing, for example nodes next
to circular boundaries in a rectangular grid. For Dirichlet boundary conditions, we
can make use of the fact that we know the value of uðx; yÞ on the boundary.
Consider a node 0 adjacent to the boundary shown in Figure 11.3. The two-
dimensional Laplacian operator acting at u0 can be approximated as [5]
2 uA uB u3 u4 1 1
r u≃ 2
2
þ þ þ þ u0 ; (11.9)
h sð1 þ sÞ tð1 þ tÞ 1 þ s 1 þ t s t
two regions (a) and (b) with different permeabilities. Let us consider an arbitrary
point 0 along the boundary. Since point 0 is part of region a, the Laplace equation is
where we use A for the unknown function here. For region a, node 1 is fictitious and
must be eliminated from the final difference equation. Node 0 is also a part of
region b, so we have
Substituting for Aa1 from Equation 11.10 and Ab3 from Equation 11.11, we find the
difference equation at boundary point 0 is
2μa 2μb
4A0 Ab1 A2 Aa3 A4 ¼ 0: (11.12)
μa þ μb μa þ μb
x x1
s¼
x2 x1
y y1
t¼ ;
y2 y1
This expression varies continuously in x and y and reduces correctly to the node
values at the corners of the rectangle.
The discrepancy between the result from using the difference approximation and
the exact result from solving the differential equation is known as the truncation
error. The error can be estimated by examining the first term in the Taylor series that
was neglected in deriving the difference formula under consideration. For a square
mesh, the error on the second derivative goes like
2h2 ∂4 u
∼ :
4! ∂x4 0
The error is proportional to h2, so one method of improving the accuracy in a finite
difference calculation is to reduce the mesh spacing. We can monitor the improve-
ment in accuracy by finding the maximum absolute value for the difference
where the superscript refers to the mesh spacing used for the solution and ði; jÞ
refers to nodes common to both mesh spacings. This approach is ultimately limited
by the growth in the size of the coefficient matrix and by rounding errors in the
numerical calculations. An alternative approach for increasing the accuracy of the
calculation is to use higher-order difference equations.
The quality of a solution can be monitored by calculating the residual for each of
the interior nodes. For a general node for the Poisson equation, the residual is
defined as
11.2 Example solution using finite differences 251
Ri; j ¼ 4 ui; j ui; jþ1 ui; j1 uiþ1; j ui1; j h2 μi; j Ji; j : (11.15)
If the difference equation exactly satisfies Poisson’s equation, the residual should
be 0.
For problems using iteration algorithms, we can compute
for the unknown function at the node ði; jÞ, where the superscript refers to the
iteration number. For these methods, one can estimate the quality of the solution by
calculating the difference at all the nodes. Let M refer to the absolute value of the
largest difference in the mesh.
M ¼ max jeij j
For the 5-point Laplacian operator in Equation 11.7, the error ε between the exact
solution of the difference equation and the approximate solution after n iterations is
bounded by [8]
Mρ2
ε ≤ ; (11.16)
4h2
where ρ is the radius of the smallest circle that encompasses the entire field region.
In problems where iron saturation is a consideration, the permeability of the iron
can be made a variable at each of the nodes in the iron regions.[9] The permeabil-
ities are stored on a separate mesh. After each iteration of the potential, the field in
the iron region is updated. A table of B-μ values can be used to relate the perme-
ability to the field at the node. The mesh of permeability values is then updated
using, for example, an under-relaxation algorithm.
In applying the Poisson equation here, four types of node patterns are required,
as shown in Figure 11.7. In each case, node 0 refers to the node we are currently
evaluating. For a general interior node where all the neighbor nodes are in the same
region, we apply the pattern (a), which results in the relation
4A0 ¼ A1 þ A2 þ A3 þ A4 þ f ; (11.17)
where f ¼ h2 μ J for nodes inside the conductor and 0 otherwise. For nodes on the
symmetry plane, pattern (b) gives
4A0 ¼ A1 þ 2A2 þ A3 þ f :
For the left side of the iron sheet, we can use Equation 11.12 for pattern (c) with
μa ¼ 1 and μb ¼ 100. For the right side of the sheet, we use Equation 11.12 for
pattern (d) with μa ¼ 100 and μb ¼ 1.
11.2 Example solution using finite differences 253
For problems with a very large number of unknown nodes, it is not practical to
solve the matrix equation using direct methods. Instead iterative methods must be
used. A common method is to use the Successive Overrelaxation (SOR) algorithm.
[10, 11] Let us define Anj;k to be the value of the potential at the interior node at
location ðj; kÞ after n iterations. On the next iteration, we update the value of the
potential according to the prescription
j;k ¼ ð1 αÞ Aj;k þ α Aj;k ;
Anþ1 n
(11.18)
Thus the updated value of the potential has two contributions. The first term is an
adjustable fraction of its value on the previous iteration. The second term is
a fraction of the Poisson equation residual at the node, calculated from the values
of the potential at the neighbor nodes. Note that the calculation of the residual uses
values for two nodes that have already been updated for a given iteration and values
for two nodes from the previous iteration. The iterations continue until
Anþ1 An
j;k j;k
max ≤τ
Anj;k
over all the interior nodes.1 The tolerance τ ¼ 105 was used in this example. This
criterion was satisfied after 7,511 iterations.
The magnetic field was calculated at the center of every square formed by four
neighbor nodes, as shown in Figure 11.8.
1
Bx ¼ ∂ y A z ¼ ðA1 þ A2 þ A3 A4 Þ
2h
1
By ¼ ∂x Az ¼ ðA1 þ A2 A3 A4 Þ:
2h
The results of the calculations for the magnetic field are shown in Figure 11.9.
1
Or just the difference in values if the potential is 0.
254 Numerical methods
Figure 11.9 Magnetic field pattern for the example finite difference problem.
The line integral in Equation 11.19 vanishes since we require that the potential
satisfy either Dirichlet or Neumann boundary conditions everywhere on the
boundary. If the reluctivity is constant over an element, we can perform the
integration over b to get the simplified energy functional
ðð ( " 2 2 # )
γ ∂A ∂A
F¼ þ J A dx dy: (11.21)
2 ∂x ∂y
If we write this expression for each of the three nodes, we have three equations that
can be solved for the three unknown coefficients ci in terms of the potentials and
coordinates at the nodes. Substituting the result back into Equation 11.22, we find
ðx2 y3 x3 y2 Þ þ ðy2 y3 Þx þ ðx3 x2 Þy
A ¼ A1
2S
ðx3 y1 x1 y3 Þ þ ðy3 y1 Þx þ ðx1 x3 Þy
þA2 (11.23)
2S
ðx1 y2 x2 y1 Þ þ ðy1 y2 Þx þ ðx2 x1 Þy
þ A3 ;
2S
where S is the area of the triangle.
S ¼ ½½ðx2 y3 x3 y2 Þ þ ðy2 y3 Þx1 þ ðx3 x2 Þy1 (11.24)
The coefficients of the node potentials in this equation are known as shape functions,
ζ.[14] Thus we can also write the interpolation function for the potential as
A ¼ ζ 1 A1 þ ζ 2 A2 þ ζ 3 A3 : (11.25)
The shape function ζ 1 has the properties that
ζ 1 ðx1 ; y1 Þ ¼ 1
ζ 1 ðx2 ; y2 Þ ¼ 0
ζ 1 ðx3 ; y3 Þ ¼ 0
∂F γ h i
¼ ðy2 y3 Þ2 A1 þ ðy3 y1 Þðy2 y3 ÞA2 þ ðy1 y2 Þðy2 y3 ÞA3
∂A1 4S
γ h 2
i
þ ðx3 x2 Þ A1 þ ðx1 x3 Þðx3 x2 ÞA2 þ ðx2 x1 Þðx3 x2 ÞA3
4S
ðð
½ðx2 y3 x3 y2 Þ þ ðy2 y3 Þx þ ðx3 x2 Þy
J dx dy ¼ 0
2S
with analogous expressions for the derivatives with respect to A2 and A3. For
elements containing current, the second integral can be evaluated by assuming
that J is constant and that x and y are evaluated at the centroid of the triangle.
x1 þ x2 þ x3
xc ¼
3
y1 þ y2 þ y3
yc ¼ :
3
In this case, the numerator in the last term is 2S/3, so the integral has the value JS/3.
Thus minimization of the functional over the triangular element leads to the matrix
equation
2 32 3 2 3
C C12 C13 A1 1
γ 4 11 J S
C21 C22 C23 54 A2 5 ¼ 4 1 5: (11.27)
4S 3
C31 C32 C33 A3 1
equation. A major advantage of the integral equation method over methods based
on the solution of a differential equation is that the mesh only needs to encompass
the iron region (and possibly the conductor region if the current density is not
uniform).[16, 17] The boundary condition at infinity follows naturally and does not
have to be imposed at the edge of a mesh. An important disadvantage is that the
resulting matrix equation is dense, so the solution time grows rapidly as the number
of elements is increased. Also the flux density computed near the iron elements can
be strongly affected by the discretization. There are many ways to formulate
a solution of Poisson’s equation using integral equations.[17, 18, 19] In addition,
it is also possible to formulate procedures which combine differential and integral
equation techniques.[20]
We describe here an integral equation approach that uses the magnetization of
iron elements as the unknown function.2 Recall that the magnetization is related to
the magnetic field intensity by
!
! B ! !
M ¼ H ¼ χðHÞ H ; (11.29)
μ0
where χ is the susceptibility. The field intensity has contributions from both
conductor currents and from the magnetization in the iron. The total field inten-
sity is
! ! !
H ¼ H c þ H m;
2
This procedure was adopted by a group at the Rutherford High Energy Laboratory in the development of the
GFUN program.
260 Numerical methods
Break the iron region into N elements and assume the magnetization is constant
over the area of each element. The magnetization at some element i will depend on
the field at i due to the conductors and on the field at i due to the magnetization from
all the elements. Thus we have
" ðX ! #
N
! ! 1 ! R ij
M i ¼ χi H ci ri M j · 3 dVj : (11.32)
4π j¼1
Rij
It is important to note that the components of G depends on the directions of the unit
vectors used to define M, but do not depend on the magnitude of M. We can then
rewrite Equation 11.32 as
" #
! ! X
N
!
M i ¼ χi H ci þ Gij M j :
j¼1
Two of the terms vanish because the derivative in field coordinates acts on the unit
vector along the magnetization in source coordinates. Another term vanishes because
of the curl operator acting on the linear vector R. Thus only the third term on the
right-hand side remains. Writing this out in terms of components, we have
! !
ðM^ j ·ri Þ R i j
R3i j
2 3
ðx x Þ^
x þ ðy y Þ^
y þ ðz z Þ^
z
¼ M ^ jx ∂ix þ M ^ jz ∂iz 4
^ jy ∂iy þ M i j i j i j 5
2 2 2 3=2
fðxi xj Þ þ ðyi yj Þ þ ðzi zj Þ g
^ jx h
M !i M^ jy h !i M ^ jz h !i
¼ R 2
i j ^
x 3ðx i x j Þ R þ R 2
i j ^
y 3ðy i y j Þ R þ R 2
i j^z 3ðzi z j Þ R :
R5i j R5i j R5i j
After inserting this expression into Equation 11.33, we can identify the three-
dimensional coupling constants
ð 2 2
1 Rij 3ðxi xj Þ
Gix; jx ¼ dVj
4π R5ij
ð
1 3ðxi xj Þðyi yj Þ
Gix; jy ¼ dVj
4π R5ij
and similarly for the other components.[21, 22] The G coupling matrix is sym-
metric. There are constraints on the sum of the diagonal elements.[21]
0 if i ≠ j
Gix; jx þ Giy; jy þ Giz; jz ¼
1 if i ¼ j
where dSj is the two-dimensional area element. This can be written in the form
ð
1 ^ ^ jy I2 dSj ;
Gij ¼ ½M jx I1 þ M (11.35)
4π
where
ð∞
R2i j ^x 3ðxi xj Þ½ðxi xj Þ ^x þ ðyi yj Þ ^y zj ^z
I1 ¼ dzj
∞ R5i j
and
ð∞
R2i j ^y 3ðyi yj Þ½ðxi xj Þ^x þ ðyi yj Þ ^y zj ^z
I2 ¼ dzj :
∞ R5i j
Let rij be the distance between the observation point and the centroid of the iron
element in the x-y plane. Then the integral I1 can be broken into the three simpler
integrals3
ð∞
dzj 2
3=2
¼ 2
∞ frij þ zj g
2 2 rij
ð∞
dzj 2
5=2
¼ 4
∞ fr2
ij þ zj g
2 3rij
ð∞
zj
dzj ¼ 0
∞ fr2 þ z 2 g5=2
ij j
2 4ðyi yj Þ !
I2 ¼ 2
^y r ij :
rij r4ij
Inserting these results into Equation 11.35, we find the two-dimensional coupling
constants are the dimensionless, geometric factors
3
GR 2.271.5, 2.263.3, 2.271.7.
11.4 Integral equation method 263
ð
1 ðxi xj Þ2 ðyi yj Þ2
Gix; jx ¼ dSj
2π r4ij
ð
1 2ðxi xj Þðyi yj Þ
Gix; jy ¼ Giy; jx ¼ dSj (11.36)
2π r4ij
ð
1 ðyi yj Þ2 ðxi xj Þ2
Giy; jy ¼ dSj :
2π r4ij
Once we know the coupling constants G, we can solve Equation 11.34 to find the
magnetization in each of the iron elements. Then the field at any position can be
found from the sum of the fields due to all the current elements together with the
sum of all the fields due to the iron magnetizations. In applications where saturation
in the iron is important, the permeability of all the iron elements must be recom-
puted using the magnetizations and a μ–H table for the iron material. The process is
then iterated until the maximum change in permeability in any element is less than
some tolerance value.
Assume, for example, that there are two iron elements in the first quadrant. We can
find the right-hand side of Equation 11.34 using the techniques in Chapters 4 and 5.
We evaluate the field at the centroids of each of the iron elements. The magnetization
has a constant magnitude and direction in each element. The coupling constants for
the field due to the magnetizations can be found from Equation 5.72.
þ
M dz
H ðzo Þ ¼ :
4π i z zo
We can determine the coupling constants Gxx and Gyx by evaluating H* with M = 1
and the coupling constants Gxy and Gyy by evaluating H* with M ¼ i.
The field components due to the magnetization in each iron element can be found
from
Hmx ¼ Gxx Mx þ Gxy My
(11.37)
Hmy ¼ Gyx Mx þ Gyy My :
This gives the field anywhere outside the iron block. However, when the observation
point is inside the block, the numerical procedure must ensure that H and M point in
opposite directions, as they must inside a magnetic material. For the dipole config-
uration, we define the matrix coefficients C as sums over the coupling constants in the
four quadrants. For example,
where the minus signs take into account the reversal in the sign of Mx in the second
and fourth quadrants. The other three coefficients are similarly defined. The matrix
Equation 11.34 can be written for the case of two iron elements as
11.5 The POISSON code 265
ðs1 C1x1x ÞM1x C1x1y M1y C1x2x M2x C1x2y M2y ¼ Hc1x
C1y1x M1x þ ðs1 C1y1y ÞM1y C1y2x M2x C1y2y M2y ¼ Hc1y
C2x1x M1x C2x1y M1y þ ðs2 C2x2x ÞM2x C2x2y M2y ¼ Hc2x
C2y1x M1x C2y1y M1y C2y2x M2x þ ðs2 C2y2y ÞM2y ¼ Hc2y ;
where si ¼ 1=ðμri 1Þ and the numeral subscripts refer to the two iron elements in
the first quadrant.
After solving the matrix equation, M is known for all the iron blocks. The
contribution of the iron to the field at any location can be found using
Equation 11.37, where G is now evaluated for the desired field point.
4
http://laacg.lanl.gov/laacg/services
266 Numerical methods
þ ð
∂A ∂A ! !
γðBÞ ^x ^y ·dl ¼ μ0 J ·^n dS:
∂y ∂x
The contour around each mesh point follows a twelve-sided path through the
interior of the six surrounding triangles. After a lengthy calculation,[24]
a difference equation for the potential at node 0 can be derived in terms of the
potentials at the six neighbor nodes as
X6
μ X6
Ai wi þ 0 Ji Si
i¼1
3 i¼1
A0 ¼ :
X 6
wi
i¼1
In this equation, Si is the area of the triangle. The wi are coupling coefficients that
involve the parameters of the triangles on adjacent sides of the line connecting A0 to
neighboring node i. Looking at the diagram in Figure 11.13,
w1 ¼ ½ðγ1 cot θ1 þ γ2 cot θ4 Þ
with similar expressions for the other five couplings.
The vector potential varies linearly inside any triangle. As a result, the magnetic
field is constant over the area of the triangle. Values for the potential are updated
using the successive over-relaxation algorithm. The new values of the field are then
used to estimate new values for γ and for the couplings w. An under-relaxation
algorithm is used to update the final values of the couplings for each iteration
winþ1 ¼ ð1 αÞwin þ αwinew ;
Table 11.1. The first REG command defines the problem domain and specifies the
mesh size and the boundary conditions. The PO commands define points around
the boundary of regions. The second region defines the conductor and specifies the
current. The third region defines the iron sheet. POISSON sets up a triangular mesh
using this information and then solves the Poisson equation using the SOR algo-
rithm. For this example, the program used 112,896 mesh points, converged in 1,160
iterations, and had an average residual of 5 107 . The resulting field distribution,
shown in Figure 11.14, agrees qualitatively well with the finite difference result in
Figure 11.9. The shielding effects of the iron sheet are clearly apparent in the figure.
0 20 40 60 80 100
100
90
80
70
60
50
40
30
20
10
0
0 20 40 60 80 100
by constraining the coil geometry. For example, by specifying that the unknown
currents lie on a cylindrical current sheet, it is possible to use Fourier-Bessel
transforms to find the azimuthal and longitudinal current components that produce
a specified target field inside a magnet aperture.[25] Another interesting approach
used a numerical variational process to modify the contours of a uniform current
density block conductor.[26] The target field was expressed in terms of a multipole
expansion of the transverse field in the aperture. Higher-order multipoles were
minimized by varying the geometry of the outer boundary of initial circular or
elliptical current blocks.
A powerful technique for solving inverse problems is to make use of numerical
optimization methods. Let us consider in more detail the numerical solution for two
interesting inverse problems. As the first example, assume we have a solenoid
channel with a constant axial field B1 and that we need to design an interface region
to a second solenoid channel with constant axial field B2 . Assume the interface has
length L measured from the center of the last magnet in the first channel to the
center of the first magnet of the second channel. Assume in addition that the
11.6 Inverse problems and optimization 269
transition region and second channel have to accept the full magnetic flux present in
the first channel. Then the desired field profile in the interface region must satisfy
the four constraints
Bz ð0Þ ¼ B1
Bz ðLÞ ¼ B2
dBz
ð0Þ ¼ 0
dz
dBz
ðLÞ ¼ 0:
dz
A model field profile that satisfies these constraints is
B1
Bz ðzÞ ¼ ; (11.38)
1 þ c z2 þ d z3
where
3 ðB1 B2 Þ
c¼
B2 L2
2 ðB1 B2 Þ
d¼ :
B2 L3
If r1 is the inner radius of the coils in the first channel, then the requirement for
constant flux puts an additional constraint on the allowed inner radius of the
downstream coils.
sffiffiffiffiffiffiffiffiffiffiffi
B1
rðzÞ ≥ r1 :
Bz ðzÞ
For example, let the coil C1 be the last solenoid in a 10 T channel with a fixed
inner radius of 10 cm and coil C14 be the first solenoid in a 2 T channel. Assume the
transition region is 7 m long and contains 12 solenoids that are 45 cm long,
separated by 5 cm, and have adjustable inner radius, radial thickness, and current
density. The axial field for each solenoid uses Equation 7.46. The merit function
f for the minimizer compares the desired value of the field at N locations zi from
Equation 11.38 with the calculated sum of the fields from all the coils, each with
a set of parameters aj .
" ! #2
X
N X
14
f ¼ Bz ðzi ; aj Þ Bz ðzi Þ :
i¼1 j¼1
270 Numerical methods
Figure 11.15 The optimized axial field (dots) and the desired field profile (line) in
the transition region between two solenoid channels.
The calculation shown here used N = 36. Minimization of this function used
methods that do not require calculation of the derivatives.[27] The initial mini-
mization was done using a simplex algorithm. The most useful parameters to adjust
were the current densities and the inner radii of the coils. The axial parameters are
severely constrained here by the chosen geometry for the transition region. After
a preliminary solution had been found, the Powell direction-set method was used
for the final minimization. The optimized axial field is compared with the desired
field profile in Figure 11.15.
As a second example of optimization, let us consider the design of the central
section of a long dipole magnet with a circular cross-section. Assume that field
quality in the dipole aperture is the matter of concern and that we want to minimize
the strength of the first four allowed harmonics of the dipole field. We saw in
Chapters 4 and 5 that the multipole coefficients depend on the limiting angles of
annular conductor blocks. In order to eliminate four multipoles we will need to use
at least three blocks. We choose here a conductor design with two radial layers,
each of which has two conductor blocks, as shown in Figure 11.16. The contribu-
tion to the multipoles from an annular conductor block with constant current
density was given in Equations 5.68–5.70. We again use a minimization algorithm,
where the merit function is now given by
" ! #2
X
4 X
16
f ¼ wi bn ðni ; aj Þ e
b n ðni Þ :
i¼1 j¼1
11.6 Inverse problems and optimization 271
The index i sums over desired multipole orders, the index j sums over coils, and the
normalized multipole ratio is defined as
Bn ½T=mn1 rn1
ref ½m
n1
bn ¼ :
B1 ½T
The eb n factors are the desired values of the multipole ratios, which we take here as
0. The aj are the set of parameters that describe conductor block j. The parameters
of the coils in quadrants 2–4 are related to the parameters of the coils in quadrant 1
by the dipole symmetry. For this calculation the adjustable parameters are the end
angles of the blocks nearest the midplane and the start and end angles of the second
block in each layer. The start angle of the two blocks nearest the midplane are made
as close to 0 as possible to maximize the dipole field. The wi are weights that
determine the importance of satisfying the constraint on multipole ni . The reference
radius used for the multipole calculations was 2/3 of the magnet aperture. After
minimization, the allowed multipole ratios b3 ; b5 ; b7 , and b9 have strengths ~104 .
In the design of actual magnets,[28] a need for high precision field quality
may require that allowed multipoles higher than b9 are also minimized.
In addition, the conductor may have to be described in terms of individual
turns of the cable separated by the appropriate insulation thickness, instead of
the continuous conductor blocks used here. This introduces the additional
constraint that there must be an integral number of turns in a conductor block.
In addition, if the coils are surrounded by an iron shell, saturation effects, which
cause the multipole strength to vary with the excitation current, may have to be
taken into account.
272 Numerical methods
References
[1] A. Frisiani, G. Molinari & A. Viviani, Introduction, in M. Chari & P. Silvester
(eds.), Finite Elements in Electrical and Magnetic Field Problems, Wiley, 1980,
p. 1–10.
[2] N. Gershenfeld, The Nature of Mathematical Modeling, Cambridge University Press,
1999, p. 86–91.
[3] K. Binns & P. Lawrenson, Analysis and Computation of Electric and Magnetic Field
Problems, Pergamon Press, 2nd ed., 1973, p. 246.
[4] D. Jones, Methods in Electromagnetic Wave Propagation, vol. 1, Oxford University
Press, 1987, p. 113.
[5] G. Smith, Numerical Solution of Partial Differential Equations: Finite Difference
Methods, 2nd ed., Oxford University Press, 1978, p. 216.
[6] K. Binns & P. Lawrenson, op. cit., p. 268.
[7] W. Press, S. Teukolsky, W. Vetterling & B. Flannery, Numerical Recipes in Fortran,
2nd ed., Cambridge University Press, 1992, p. 117.
[8] W. Milne, Numerical Solution of Differential Equations, Dover, 1970, p. 217.
[9] G. Parzen & K. Jellett, Computation of high field magnets, Part. Acc. 2:169, 1971.
[10] W. Press, et al., op. cit., p. 857–860.
[11] K. Binns & P. Lawrenson, op. cit., p. 260–265.
[12] N. Gershenfeld, op. cit., p. 93–101.
[13] M. Chari, Finite element solution of magnetic and electric field problems in electrical
machines and devices, in M. Chari & P. Silvester (eds.), Finite Elements in Electrical
and Magnetic Field Problems, Wiley, 1980, p. 87–107.
[14] R. Gallagher, Shape functions, in M. Chari & P. Silvester (eds.), ibid., p. 49–67.
[15] P. Silvester & M. Chari, Finite element solution of saturable magnetic field problems,
IEEE Trans. on Power Apparatus and Systems 89:1642, 1970.
[16] J. Simkin & C. Trowbridge, Three dimensional nonlinear electromagnetic field
computations using scalar potentials, IEE Proc. 127:368, 1980.
[17] J. Simkin & C. Trowbridge, Magnetostatic fields computed using an integral equation
derived from Green’s theorems, Proc. Compumag, Rutherford Appleton Laboratory,
Oxford, 1976, p. 5.
[18] C. Trowbridge, Applications of integral equation methods to the numerical solution of
magnetostatic and eddy current problems, in M. Chari & P. Silvester (eds.), op. cit.,
p. 191–213.
[19] C. Trowbridge, Progress in magnet design by computer, Proc. 4th Int. Conf. Magnet
Technology, Brookhaven National Laboratory, 1972, p. 555.
[20] B. McDonald & A. Wexler, Mutually constrained partial differential and integral
equation field formulations, in M. Chari & P. Silvester (eds.), op. cit., p. 161–190.
[21] M. Newman, C. Trowbridge & L. Turner, GFUN: an interactive program as an aid to
magnet design, Proc. 4th Int. Conf. Mag. Tech., Brookhaven National Laboratory,
1972, p. 617.
[22] F. Mingwu, S. Hanguang & W. Jingguo, Some experiences of using integral equation
method to calculate magnetostatic fields, IEEE Trans. Mag. 21:2185, 1985.
[23] A. Winslow, Numerical solution of the quasilinear Poisson equation in a nonuniform
triangle mesh, J.Comp. Phys. 1:149, 1966.
[24] J. Billen & L. Young, Poisson-Superfish, Los Alamos National Laboratory report LA-
UR-96–1834, 2006, p. 577–596.
[25] R. Turner, A target field approach to optimal coil design, J. Phys. D: Appl. Phys. 19:
L147, 1986.
References 273
[26] G. Morgan, Two dimensional, uniform current density, air core coil configurations for
the production of specified magnetic fields, Proc. 1969 Part. Accel. Conf., Washington
DC, p. 768.
[27] W. Press, et al., op. cit., p. 402–413.
[28] E. Bleser, et al., Superconducting magnets for the CBA Project, Nuc. Instr. Meth.
Phys. Res. A 235:435, 1985.
Appendices
A
Symbols and SI units [1, 2]
References
[1] E. R. Cohen (ed.), The Physics Quick Reference Guide, American Institute of Physics,
1996, p. 37–47.
[2] D. Halliday & R. Resnick, Physics for Students of Science and Engineering, Wiley,
1962, appendix G.
274
B
Vector analysis
Vector analysis plays an essential role in describing the theory of magnetic phe-
nomena.[1, 2] A vector V is a quantity that has both a magnitude and a direction.
A scalar S is a quantity that only has an associated magnitude. Vector fields are
functions that describe a physical quantity at every point in space.
The vector differential operator (del) is
∂ ∂ ∂
r¼ ^x þ ^y þ ^z :
∂x ∂y ∂z
When r is applied to a scalar function, it results in a vector known as the gradient.
∂S ∂S ∂S
rS ¼ ^x þ ^y þ ^z :
∂x ∂y ∂z
The gradient gives a measure of the rate of change of a vector. The dot product of r
with a vector forms a scalar known as the divergence.
∂2 S ∂2 S ∂2 S
r2 S ¼ r · rS ¼ þ þ :
∂x2 ∂y2 ∂z2
275
276 Vector analysis
The cross product of r with a vector forms another vector known as the curl.
^x ^y ^z
!
r V ¼ ∂x ∂y ∂z :
Vx Vy Vz
The curl gives a measure of the tendency of the vector to circulate around some
source. According to Helmholtz’s theorem,[3] a vector function that is bounded at
infinity can be uniquely defined by specifying its divergence and its curl.
If we consider a volume of space V enclosed by a surface S, then we find that any
changes in a vector W inside the volume must match the flux of W through the
boundary surface. This is the basis for an important result known as Gauss’s diver-
gence theorem.[4]
ð ð
! !
r·W dV ¼ W ·^n dS;
where n is the normal vector to the surface. If, on the other hand, we break up the
surface S into a number of smaller areas and look at the net result of the circulation
in all the subareas, we find that the circulations cancel in the interior of the region
and only give a net result around the perimeter of S. The result is known as Stokes’s
theorem.[5]
ð þ
! !!
ðr W Þ·^n dS ¼ W ·dl
Some other integral relations involving the gradient, divergence, and curl are less
common, but still useful.[6]
ð ð
rP dV ¼ P ^n dS
ð þ
!
^
n rP dS ¼ P dl
ð ð
! !
r W dV ¼ W ^n dS;
where P is a scalar function and S is the surface that bounds the volume V.
The differential vector operators for cylindrical and spherical coordinate systems
are given in Table B1.
Some important vector identities are
! ! ! ! !! ! !!
A ðB C Þ ¼ B ðA ·C Þ C ðA ·B Þ (B.1)
Vector analysis 277
Cylindrical Spherical
1 ^ þ ∂z S ^z 1 1
rS ∂ρ S ^ρ þ ∂ϕ S ϕ ∂r S ^r þ ∂θ S ^θ þ ^
∂ϕ S ϕ
ρ r r sin θ
! 1 1 1
r· V ∂ρ ðρ Vρ Þ þ ∂ϕ Vϕ þ ∂z Vz ∂r ðr2 Vr Þ
ρ ρ r2
1
þ ½∂θ ðVθ sin θÞ þ ∂ϕ Vϕ
rsin θ
!
r V 1 ^ 1
∂ϕ Vz ∂z Vϕ ^ρ þ ð∂z Vρ ∂ρ Vz Þ ϕ ½∂θ ðVϕ sin θÞ ∂ϕ Vθ ^r
ρ rsin θ
1 1
þ ½∂ρ ðρ Vϕ Þ ∂ϕ Vρ ^z þ ½∂ϕ Vr sin θ∂r ðr Vϕ Þ^θ
ρ rsin θ
1 ^
þ ½∂r ðr Vθ Þ ∂θ Vr ϕ
r
1 1 1
r2 S ∂ρ ðρ∂ρ SÞ þ 2 ∂2ϕ S þ ∂2z S ∂r ðr2 ∂r SÞ
ρ ρ r2
1
þ 2 ∂θ ðsin θ ∂θ SÞ
r sin θ
1
þ ∂2ϕ S
r2 sin 2 θ
!! ! ! ! ! ! ! ! !
rð A · B Þ ¼ A ðr B Þ þ B ðr A Þ þ ð A ·rÞ B þ ð B ·rÞ A (B.2)
! ! !
r·ðS V Þ ¼ V ·rS þ S r· V (B.3)
! ! ! ! ! !
r·ð A B Þ ¼ B ·ðr A Þ A ·ðr B Þ (B.4)
!
r·ðr V Þ ¼ 0 (B.5)
! ! !
r ðS V Þ ¼ rS V þ S r V (B.6)
! ! !
r r V ¼ rðr· V Þ r2 V (B.7)
r rS ¼ 0 (B.8)
! ! ! ! ! ! ! ! ! !
r ð A B Þ ¼ A ðr· B Þ B ðr· A Þ þ ð B ·rÞ A ð A ·rÞ B (B.9)
278 Vector analysis
References
[1] G. Harnwell, Principles of Electricity and Magnetism, 2nd ed., McGraw-Hill, 1949,
p. 636–649.
[2] J. Marion, Classical Electromagnetic Radiation, Academic Press, 1965, p. 447–456.
[3] L. Eyges, The Classical Electromagetic Field, Dover, 1980, p. 387.
[4] P. Lorrain & D. Corson, Electromagnetic Fields and Waves, 2nd ed., Freeman, 1970,
p. 13–16.
[5] Ibid., p. 21–22.
[6] E. R. Cohen (ed.), The Physics Quick Reference Guide, American Institute of Physics,
1996, p. 162.
C
Bessel functions
Using the method of separation of variables for the Laplace equation in cylindrical
coordinates gives rise to Bessel’s equation.[1, 2]
d dR
ρ ρ þ ðk2 ρ2 n2 Þ R ¼ 0:
dρ dρ
In this equation, R ¼ RðρÞ and k and n are separation constants. The parameter
n must be an integer to keep the azimuthal dependence of the solution single-
valued, i.e., we must have
Bessel’s equation is a second order differential equation that has two indepen-
dent classes of solution. One class involves Bessel functions of the first kind,[3]
RðρÞ ¼ Jn ðkρÞ: The behavior of the first three Bessel functions Jn are shown as
a function of kρ in Figure C1. All functions of this type are well-behaved at ρ = 0.
They are oscillatory with a decreasing amplitude that approaches zero as kρ → ∞.
The first root of the function J0 ðxÞ occurs at x = 2.405, where x = kρ. The first root of
J1 ðxÞ occurs at x = 3.832. The series expansion is
X
∞
ð1Þk xnþ2k
Jn ðxÞ ¼ :
k ¼ 0
k! ðn þ kÞ! 2
279
280 Bessel functions
dJ0 ðxÞ
¼ J1 ðxÞ:
dx
The other class of solutions to Bessel’s equation are the Bessel functions of
the second kind,[4] RðρÞ ¼ Nn ðkρÞ: The behavior of the first three Bessel functions
Nn are shown as a function of kρ in Figure C2. These solutions are also oscillatory
with decreasing amplitude that approach zero as kρ → ∞. However, they diverge at
ρ = 0, so they cannot be used in magnetostatics for any region that contains the
Bessel functions 281
origin. The Nn ðxÞ functions satisfy the same recurrence relations as Jn ðxÞ. The
derivative of N0 is given by
The solutions of this equation are known as modified Bessel functions. This same
equation can be produced by replacing k with i k in the ordinary Bessel equation.
One class of radial solutions involves the modified Bessel function In ðkρÞ.[5]
The behavior of the first three modified Bessel functions In are shown in Figure C3.
All functions of this type are well behaved at ρ = 0. They are related to the ordinary
Bessel functions by
Iν ðxÞ ¼ iν Jν ði xÞ:
The other class of solutions for the modified Bessel’s equation are the functions
Kn ðkρÞ: The behavior of the first three modified Bessel functions Kn are shown in
Figure C4. These solutions diverge at ρ = 0, so they cannot be used in any region
that contains the origin. The functions Kn satisfy the recursion relations2
2n
Kn ðxÞ ¼ Kn1 ðxÞ Knþ1 ðxÞ
x
2K 0 n ðxÞ ¼ Kn1 ðxÞ þ Knþ1 ðxÞ:
References
[1] F. Bowman, Introduction to Bessel Functions, Dover Publications, 1958.
[2] M. Abramowitz & I. Stegun (eds.), Handbook of Mathematical Functions, Dover
Publications, 1972, chapter 9.
[3] G. Arfken, Mathematical Methods for Physicists, 3rd ed., Academic Press, 1985,
p. 573–584.
[4] Ibid., p. 596–601.
[5] Ibid., p. 610–619.
1 2
GR 8.486.1,8.486.2. GR 8.486.10,8.486.11.
D
Legendre functions
Separation of variables for the Laplace equation in spherical coordinates gives the
partial differential equation
1 1
∂θ ðsin θ ∂θ YÞ þ ∂2ϕ Y þ lðl þ 1Þ Y ¼ 0
sin θ sin θ
2
for the angular dependence. The solution of this equation is given in terms of the
spherical harmonic functions Ylm [1, 2]
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
ð2l þ 1Þ ðl mÞ! m
Yl m ðθ; ϕÞ ¼ P ðcos θÞ eimϕ ;
4π ðl þ mÞ! l
where x = cos θ. The solutions of this equation are called associated Legendre
functions of the first and second kind,
fPm
l ðxÞ; Ql ðxÞg:
m
283
284 Legendre functions
l m Ylm
0 0 1
pffiffiffiffiffi
4π
rffiffiffiffiffi
3
1 0 cos θ
4π
rffiffiffiffiffi
3
1 1 sin θ ei ϕ
8π
rffiffiffiffiffiffiffiffi
5
2 0 ð3 cos2 θ 1Þ
16π
rffiffiffiffiffi
15
2 1 sin θ cos θ ei ϕ
8π
rffiffiffiffiffiffiffiffi
15
2 2 sin2 θ e2i ϕ
32π
Only the functions of the first kind have convergent power series over the complete
range 0 ≤ x ≤ 1, so we choose
ΘðθÞ ¼ Pm
l ðcos θÞ:
ð1Þm dlþm
l ðxÞ
Pm ¼ l ð1 x2 Þm=2 lþm ðx2 1Þl :
2 l! dx
The associated Legendre functions for l ≤ 3 and m > 0 are given in Table D2.
In problems with azimuthal symmetry, we have m = 0. In this case, the associated
Legendre functions reduce to the ordinary Legendre polynomials.
l m Plm
1 1 sin θ
2 1 3 cos θ sin θ
2 2 3 sin2 θ
3
3 1 ð5 cos2 θ 1Þ sin θ
2
3 2 15 cos θ sin2 θ
3 3 15 sin3 θ
The Legendre polynomials form a complete set of orthogonal functions over the
interval 1 ≤ cos θ ≤ 1. The behavior of the Legendre polynomials for l ≤ 4 are
shown in Figure D1 as a function of x. Legendre polynomials satisfy the recurrence
relation
References
[1] G. Arfken, Mathematical Methods for Physicists, 3rd ed., Academic Press, 1985,
chapter 12.
[2] M. Abramowitz & I. Stegun (eds.), Handbook of Mathematical Functions, Dover
Publications, 1972, chapter 8.
E
Complex variable analysis
Complex variables
In a Cartesian coordinate system, the complex variable z is given by
z ¼ x þ i y;
pffiffiffiffiffiffiffi
where x is called the real part of z, y is called the imaginary part of z, and i ¼ 1.
In polar coordinates, z can be written in the form
z ¼ r eiθ
¼ r ðcos θ þ i sin θÞ;
where r is called the modulus of z and θ is the argument of z. The De Moivre
formula is useful for evaluating powers of z.
z ¼ x i y:
287
288 Complex variable analysis
Care is required in working with the complex counterparts of some real func-
tions. An important function in magnetostatics is the complex logarithm function.
This is defined as
w ¼ ln z
¼ lnðr eiθ Þ
¼ ln r þ i ðθ þ 2πnÞ;
where n ¼ 0; 1; 2; …. This function has multiple branches of angular
width 2π, depending on the value of n.[3] We customarily compute this function
using the principal branch where n = 0 and where θ is in the range –π < θ ≤ π. In this
case, the function changes discontinuously when crossing the negative x axis,
which is called a branch cut.
Complex differentiation
The derivative of the complex function F is defined as [4]
lim ðz z0 Þn FðzÞ
z→z0
For a region not including the origin, the Cauchy-Riemann equations can be written
in polar coordinates as
∂u 1 ∂v
¼
∂r r ∂θ
1 ∂u ∂v
¼ :
r ∂θ ∂r
Series
Theorem E.2 (power series) [8] Let f(z) be analytic on a domain G and let zo
be an arbitrary point of G. Let d ¼ dðz0 Þ be the distance between zo and the
boundary of G. Then there exists a power series
X
∞
f ðzÞ ¼ cn ðz zo Þn
n¼0
Theorem E.3 (Taylor series) [9] Let f(z) be analytic and single-valued in an
open region G. Let a be any point in G and let C be a circle with center at a, which
together with its interior lies entirely in G. Then at every point z in C, the series
Theorem E.4 (Laurent series) [10] Let f(z) be analytic for the annular region
G: R1 < jz zo j < R2
and let C be any simple closed contour lying inside G and having zo in its interior.
Then for points z in G, the function f(z) may be expanded in the series
290 Complex variable analysis
X
∞
f ðzÞ ¼ ck ðz zo Þk ;
k¼∞
where
þ
1 f ðzÞ
ck ¼ dz
2πi ðz z0 Þkþ1
The Laurent series is valid in a region surrounding, but not including, a singularity.
Note that this series includes negative values of k. The coefficient c1 has special
significance and is known as the residue.
Complex integration
If instead of f(z), we consider the contour integral of f(z) / (z – zo), then we have the
following theorem.
Theorem E.7 (Cauchy’s Integral Formula) [13] Let f(z) be analytic within and
on a simple closed contour C. Then, if zo is a point inside C,
þ
1 f ðzÞ
f ðzo Þ ¼ dz:
2πi ðz zo Þ
This gives the value of f ðzo Þ at the singularity zo inside a region in terms of the
contour integral around the boundary C.
Complex variable analysis 291
Theorem E.8 [14] Let f(z) be analytic within and on a simple closed contour
C. Then all derivatives of f ðzo Þ exist at a point zo inside C and are given by
þ
ðnÞ n! f ðzÞ
f ðzo Þ ¼ dz:
2πi ðz zo Þnþ1
Theorem E.9 (Residue theorem) [15] Let f(z) be analytic within and on a simple
closed contour C, except for a finite number of isolated singularities inside C. Let σ be
the sum of the residues at the singular points of f(z) that lie inside C. Then
þ
1
f ðzÞ dz ¼ σ:
2πi
In other words, the value of the contour integral is 2πi times the sum of the
residues for the enclosed singularities. For a pole of order n, the residue can be
found as [16]
1 dn1
a1 ¼ lim n1 ½ðz aÞn f ðzÞ:
ðn 1Þ! z→a dz
Conformal mapping
We can define a function F that maps a complex variable z into a variable w in
another two-dimensional space.
w ¼ FðzÞ:
Assume that two curves that cross at a point zo in the z space are separated by an
angle θ. A mapping w ¼ FðzÞ is conformal, or angle preserving, if the mapped
curves in the w space cross at the point wo ¼ Fðzo Þ with the same angle θ.
dw
¼ Aðz x1 Þα1 =π1 ðz x2 Þα2 =π1 ðz xn Þαn =π1 ;
dz
where A is a complex constant, maps the interior of the polygon in the w plane onto
the upper half of the z plane and maps the boundary of the polygon onto the real axis
of the z plane.
References
[1] W. Smythe, Static and Dynamic Electricity, 2nd ed., McGraw-Hill, 1950, p. 72–94.
[2] W. Panofsky & M. Phillips, Classical Electricity and Magnetism, 2nd ed., Addison-
Wesley, 1962, p. 61–72.
[3] J. Dettman, Applied Complex Variables, Dover, 1984, p. 58–60.
[4] K. Miller, Introduction to Advanced Complex Calculus, Dover, 1970, p. 48–52.
[5] M. Spiegel, Complex Variables, Schaum’s Outline Series, McGraw-Hill, 1964, p. 63.
[6] Ibid., p. 67.
[7] W. Kaplan, Advanced Calculus, Addison-Wesley, 1952, p. 510–512.
[8] R. Silverman, Introductory Complex Analysis, Dover, 1972, p. 204–205.
[9] K. Miller, op. cit., p. 114–115.
[10] J. Dettman, op. cit., p. 163–164.
[11] Ibid., p. 81–82.
[12] W. Kaplan, op. cit., p. 516.
[13] J. Dettman, op. cit., p. 113.
[14] Ibid., p. 114.
[15] K. Miller, op. cit., p. 132–133.
[16] M. Spiegel, op. cit., p. 172–173.
[17] E. Kreyszig, Advanced Engineering Mathematics, Wiley, 1962, p. 719.
[18] J. Dettman, op. cit., p. 256–259.
[19] M. Spiegel, op. cit., p. 204, 218–219.
F
Complete elliptic integrals
Each of these integrals depends on a parameter k called the modulus that satisfies
k2 ≤ 1. The third type of integral also depends on a second parameter n called the
characteristic.[2] When the upper limit of integration is
π
θ0 ¼ ;
2
these functions define the complete elliptic integrals of the first, second, and third1
kinds.
1
One should be aware that the complete elliptic integral of the third kind is sometimes defined with a negative
sign before the factor n in the denominator.
293
294 Complete elliptic integrals
ð π=2
dθ
KðkÞ ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
0 1 k2 sin2 θ
ð π=2 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
EðkÞ ¼ 1 k2 sin2 θ dθ
0
ð π=2
dθ
∏ðk; nÞ ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
0 ð1 þ n sin θÞ 1 k2 sin2 θ
2
where k2 < 1.[2] Efficient numerical algorithms have been developed to calculate
the complete elliptic integrals.[3]
The dependences of the complete elliptic integrals of the first and second kinds are
shown as a function of k in Figure F1. Both functions have the value π/2 for k = 0.
The function EðkÞ has the value 1 for k = 1, while KðkÞ approaches ∞ as k → 1.
The behavior of the complete elliptic integral of the third kind for several values of n is
shown as a function of k in Figure F2. The function ∏ðk; nÞ increases as k increases
for all values of n. For a given value of k, the function increases as n becomes more
negative.
If the vector potential is defined in terms of complete elliptic integrals, we
need to take derivatives to find the magnetic field. In this case, we need to know
Complete elliptic integrals 295
5.0 n = –0.7
4.5
4.0
n = –0.3
3.5
Π (k,n)
3.0 n=0
2.5
n = 0.3
2.0
1.5
2
GR 8.123.2, GR 8.123.4.
296 Complete elliptic integrals
For the complete elliptic integral of the third kind, the derivatives are given by
[4, 5]
∂∏ðk; nÞ k
¼ ½EðkÞ ð1 k2 Þ∏ðk; nÞ
∂k ð1 k Þðk2 nÞ
2
∂∏ðk; nÞ 1 k2 n k 2 n2
¼ EðkÞ þ KðkÞ ∏ðk; nÞ :
∂n 2ðn 1Þðk2 nÞ n n
References
[1] CRC Standard Mathematical Tables, 14th ed., Chemical Rubber Co., 1965, p. 282.
[2] M. Abramowitz & I. Stegun (eds.), Handbook of Mathematical Functions, Dover,
1965, p. 590–591.
[3] W. Press, S. Teukolsky, W. Vetterling & B. Flannery, Numerical Recipes in Fortran,
2nd ed., Cambridge University Press, 1992, p. 254–261.
[4] National Institute of Standards and Technology, Digital Library of Mathematical
Functions, http://dlmf.nist.gov, equation 19.4.4.
[5] http://functions.wolfram.com/EllipticIntegrals/EllipticPi/introductions/CompleteElliptic
Integrals/05
Index
297
298 Index
Diffusion equation 236 energy loss 235
Dirac delta function 13, 199
Dipole magnets Images method
C-dipole 158 circular boundary 37, 79, 115, 144
excitation 154 plane boundary 34, 133
H-dipole 157, 237 Inductance
symmetries 90, 263 coaxial solenoids 20
window frame 154 ideal solenoid 19
Dirichlet problem 59, 248 mutual 19, 170, 188
Displacement current 241 self 19, 232
Divergence operator 275 Integral equation method 68, 258
Divergence theorem 276 Integrated potentials 147
Domains 28 Intrinsic coercivity 30, 217
Doublet 138 Inverse Helmholtz coils 180
Inverse problem 206, 267
Easy axis 217 Iron enhancement 36, 79, 86, 134
Eddy currents 235
Effective length 159 Langevin equation 28
Electric field intensity (E) 230 Laplace equation 41, 53
Electric flux density (D) 240 cylindrical coordinates 55
Electric scalar potential 231 polar coordinates 71
Electromotance 230 rectangular coordinates 54
Elliptic integrals 293 spherical coordinates 57
End fields 105 Laplacian operator 247, 275
Energy stored in magnetic field 20, 215, 232 Laurent series 289
Legendre functions 283
Faraday’s law 229 Legendre polynomials 177, 284
Ferromagnetism 28 Lenz’s law 229
Finite difference method 245 Linear materials 27
Finite element method 254 Lines of induction 3, 63, 66
Finite element (triangular) 255
Flux linkage 19 Magnetic charges 53, 213
Force Magnetic circuit 31
between current loops 4 C-shaped electromagnet 32
between two parallel wires 8 with permanent magnets 214
between solenoid and loop 189 Magnetic field
complex 140 block 93, 124
in terms of inductance 20 circulation 13, 26
magnetized body 50 complex 109
on charged particle 5 current sheet 80, 116
on conductor 155 divergence 11
on localized current 49 flux density (B) 3
Frenet-Serret coordinates 191 harmonic expansion 72
intensity (H) 26
Gauss’ law for magnetism 11 radial expansion 170
Gradient field 161 Magnetic field from conductors
Gradient operator 275 annular sector 97, 131
Green’s function 68 circular loop 10, 168
Green’s theorem 67 circular arc sheet 118
complex 124 cos(ϕ) sheet 121
cylindrical block 14, 126
Harmonic functions 41 dipole 33, 154, 158
Helical wiggler fields 208 finite plane sheet 115
Helix 197 helix 202
relation to solenoid 204 Helmholtz coils 178
Helmholtz coils 178 ideal multipole sheet 84
Helmholtz theorem 276 infinite plane sheet 8
Hysteresis 30 line current 7, 114
Index 299
Maxwell tricoil 181 Paramagnetism 28
overlapping circular 95 Periodic axial fields 209
pair of loops 178 Permanent magnet
parallel current sheets 82 assemblies 224
quadrupole 90, 161 circuit excitation 214
rectangular block 97, 129 materials 216
solenoid 15, 182, 189 rare earth model 217
toroid 31, 194 Permeability 4, 27, 29, 31, 162
Magnetic field from image currents Permittivity 240
annular sector in iron cavity 135 Persistent currents 105
block in circular iron 135 Perturbations 153
cos(ϕ) sheet in iron cavity 134 Planar wiggler fields 204
line current 132 Poisson equation 41, 255
line current in plane boundary 35, 133 energy functional 255
sheet in circular boundary 134 POISSON program
Magnetic field from magnetized bodies algorithm 265
bar magnet 212 conductor near iron shield 266
dipole moment 47 permanent magnet 216
doublet 139 quadrupole 161
magnetized block 25, 104, 139, 264 sample input file 267
permeable sphere 66 solenoid 190
rare earth materials 218, 221, 224 window frame dipole 156
triangular block 140 Poles (complex) 288, 291
Magnetic flux 3, 40 Potential
Magnetic measurements 241 complex 108
Magnetic moment 6, 47 scalar 51
Magnetic needle 213 vector 39
Magnetic units 274
Magnetization 23 Quadrupole magnets
Maxwell’s equations 241 circular 161
complex 112 excitation 160
Meissner effect 104 hyperbolic 152
Modified Bessel functions 202, multipole 152
208, 281 permanent magnet assembly 227
Multipole field components
annular sector 99, 137 Racetrack coil 106
approximations 90, 153 Radial axial field expansion 170
complex 137 Refraction at current sheet 18
dependence on current distribution 88 Reluctance 31, 215
dimensionless coefficients 75, 271 Reluctivity 218
errors 91, 93, 106, 138 Remanence 30, 217
inside magnet aperture 74 Residual 250
normal 89 Residue 290
rare earth materials 221 Residue theorem 291
rotating coil 243 Resistivity 30, 161
skew 89 Rotating coil measurements 242
symmetries 89, 93 Rowland ring 32
trapezoidal block 226
Saddle coil 106
Neumann problem 59, 248 Saturation
Non-allowed harmonics 93 cause 29
effects 31, 86, 155, 157, 162, 237
Ohm’s law 236 modelling 251, 258, 263
Optimization Scalar potential 51
matching solenoid channels 268 Scalar potential examples
minimizing higher multipoles 270 cos(mθ) sheet in circular cavity 87
Orthogonal equipotentials 111 dipole 151
300 Index
Scalar potential examples (cont.) Toroid 193
inside magnet aperture 73 Torque 5
line current 114
magnetized body 52, 102, 213 Vector identities 276
permeable sphere 66 Vector potential
quadrupole 152 block conductors 94
rare earth magnets 219 circulation 39
Schwarz-Christoffel transformation current sheet 81
292 definition 39
line current near perpendicular divergence 41
planes 144 Vector potential examples
Segmentation efficiency 227 bifilar conductors 43
Segmented magnet assemblies circular loop 167, 169
224 cos(mθ) sheet in circular cavity 85
Separation of variables 53 dipole 151
Shape functions 256 helical conductor 200
Sheet theorem 119 ideal multipole sheet 84
Skin depth 240 ideal multipole block 100
Skin effect 239 inside magnet aperture 72
Soft magnetic alloys 30 line current 42, 77
Solenoid magnet line current in circular cavity 77
block model 189 line current in slot 59
derivative of field 183 localized current 47
sheet model 182 magnetized body 49
Spherical harmonics 283 quadrupole 152
Successive overrelaxation 253 sheet solenoid 187
Successive underrelaxation 266
Superconductors 104 Wiggler 204, 206
Susceptibility 27
Stokes Theorem 276 Zonal harmonic expansion 172