Magnetic forces in and on a magnet

Given the shape of a magnet and its magnetization, point by point, which force does it exert on itself, also point by point? We explain what 'force' means in such a context and how to define it by using the Virtual Power Principle. Mathematically speaking, this force is a vector-valued distribution, with Dirac-like concentrations on surfaces across which the magnetization is discontinuous, i.e., material interfaces. To find these concentrations, we express the force as the divergence of a (symmetric) 2-tensor which generalizes a little the classical Maxwell tensor.


ALAIN BOSSAVIT
with respect to the mechanical configuration (here denoted by u all along). The main issue there is to link the magnetization vector M, the one that appears in the constitutive law B = µ 0 (H + M ), with the material properties, in a way that conforms with the intuitive idea of 'permanent magnet'. We do that (Section 3) by introducing a 'magnetostriction characteristic' m(x, β), supposed to be provided by measurements, where β is a linear map that characterizes the deformation of matter about point x. Then M (x) = m(x, 1), and the partial derivative ∂ β m (evaluated at β = 1) will play a role in the description of magnetostrictive forces. An important property of m is assumed, namely r m(x, 1) = m(x, r) for all rotations r, in compliance with the (strangely named... [12]) "Material Frame Indifference Principle". Thanks essentially to this, we can distinguish four contributions to the total force, namely Laplace force, inhomogeneity force, anisotropy force and magnetostrictive force, two of which are not plain vector fields but vector-valued distributions, in the Schwartzian sense of the word.
The non-vectorial part of the total force distribution exactly corresponds to the already mentioned force concentration at material interfaces. (It's a vector-valued field, but its domain is restricted to such interfaces.) After a reinterpretation of the classical 'Maxwell stress tensor' (Section 4), we use the so-called 'pillbox trick' [1], [7], to generate analytic formulas for these interface forces (Section 5).

On notations.
We work in the familiar environment of classical physics, the real three-dimensional Euclidean space E 3 , equipped with a dot product ' · ' and with the metric this generates (the norm (X · X) 1/2 of a vector X is denoted by |X|). Assume an orthonormal basis for the few definitions that follow.
Given two vector fields X and Y, we define X ·Y as the scalar field x → X(x)·Y (x) and X · ∇Y as the vector field whose components are, using Einstein's convention, X j ∂ j Y i . For ∇Y · X, the components are ∂ i Y j X j . One has as easily checked. A 2-tensor T can be described as a field of linear maps T (x), acting on vectors anchored at x. By X · T and T · Y , we mean the vector fields X i T ij and T ij Y j . By . The (right-handed) divergence div T of a 2-tensor T is the vector field ∂ j T ij . The 2-tensor whose components are δ ij (Kronecker symbol) is called δ, properly defined by X · δ · Y = X · Y for all X, Y. Given a function g, one has div (g δ) = ∇g. Given 2-tensors S and T , their scalar product is defined as the scalar field (S : Integrals such as f · v should be understood as (f (x) · v(x)) dx, with the whole space as integration domain (otherwise, the domain is indicated) and dx as the Euclidean measure of volumes. We indulge in a few notational abuses, most of them innocuous, with one exception: When f is a vector-valued distribution, not merely a vector field, its effect on the test-field v (a smooth, compactly supported vector field) is denoted by f · v, as if f were a vector field. (See (10) for a typical example.) 'Inverted commas', like this, intend to signal the few concepts and expressions familiar in Electromagnetism that may not be more generally known.
2. The concept of magnetic force.
2.1. The basic electromagnetic equations. We shall work on the so-called "Magnetoquasistatics with motion" model, in which displacement currents and electric charges are neglected (implying that Coulomb forces are ignored too), so the Maxwell equations take the form Matter moves at the velocity v, so the particle that was at point x at time 0 reaches at time t the point x + u(t, x), where u(t, x) = t 0 v(s, x)ds, and ∂ t u = v. We call 'configuration parameter', in what follows, the vector field x → u(t, x), denoted by u. (Beware u and v are data in (3)(4), like J s , not unknowns.) It is convenient to assume that v is defined over all space, not only at points occupied by matter, and that this extension is smooth with compact support. This way, the map x → x + u(t, x) is a diffeomorphism for t small enough, which is all right since we are only interested in the force field at time 0 in the 'reference configuration' u = 0. (There is of course an infinity of such possible extensions to the vacuum or to regions occupied by non-conductive dielectrics. But they all lead to the same force [2].) In (3), J s represents a source current, maintained by some external agency, and E + v × B is the electric field in the comoving frame (again, with a slight departure from Relativity [11], as if the Lorentz factor γ was equal to 1), hence the induced current J = η u (E + v × B) by Ohm's law. The conductivity η u changes as matter moves, which the subscript u recalls. (For simplicity, we assume that u = 0 and η = 0 on the support of J s .) We need a constitutive law to link the 'magnetic induction' B to the 'magnetic strength' H in (3)(4). It will be B = µ 0 (H + M (u)), where the 'magnetization' M only depends on the configuration u, not on the local values of B or H. Equivalently (we shall often prefer this variant), H = ν 0 B − M (u), where the 'reluctivity' ν 0 is the inverse of the (vacuum's) 'permeability' µ 0 . A 'magnet' is then, in still vague terms, a connected part of the support of M (u), and its dependence on only u makes it 'permanent', in a sense compatible with the common acception of this word. We assume that no other magnetic materials than such magnets are around. How M does depend on u will be described later: The idea is that "magnetization follows the matter in its motion." We also need to provide initial conditions for (3)(4). To this effect, set where E is the 'electric field' and A 0 a vector field whose curl is the value of B at time zero. Setting B = rot A, then, satisfies (4). Last, replace E by −∂ t A in (3). This results in a parabolic equation in A, with initial condition a familiar model (similar to convection-diffusion of a vector quantity A). The first term on the left is minus the 'induced current' J. There is no uniqueness of A, but B = rot A is uniquely determined, as well as the induced current J = rot H − J s , and hence E and A, at places and times where η u = 0.

2.2.
Magnetic energy. We now introduce 'magnetic energy'. First, suppose v = 0 in (5) for all t. Dot-multiply both sides by ∂ t A ≡ −E, and integrate over all space. This gives, taking B = rot A and d t M (u) = 0 into account, Since J = −η u ∂ t A, the first term on the left is the integral over all space of J 2 /η u , in which one recognizes the totality of Joule losses. On the right-hand side, − ∫ J s ·E represents the power brought to the system by the network that provides the sourcecurrent J s . (Cf. [2] for a study of these points.) The braced term on the left must then be, up to an additive constant, the energy stored in the magnetic field, hence our definition of magnetic energy as Implicitly, this fixes the constant at 0, so Ψ(u, 0) = 0, which is satisfying (null magnetic field B = rot A, zero stored energy), but of course an arbitrary choice. Now, come back to (5) with a nonzero velocity. One sees, by differentiating Ψ(u, A) with respect to A, that (5) can be written as Dot-multiply by ∂ t A and integrate. By the chain rule, the rate of variation in time of the magnetic energy is d The Joule losses being the energetic balance is, finally, 2.3. Definition of magnetic force. Formula (7) is valid whatever the motion 3. Computing the magnetic force.
3.1. The directional derivative trick. In need of the derivative of M (u) with respect to u at u = 0, we look for the directional derivative lim t=0 (M (tv)−M (0))/t, where v is a fixed smooth vector field. This limit is ∂ t (M (tv)) at t = 0. In more physical mood, we consider, separately for each v, a virtual evolution t → tv of the configuration, during which the magnetization field is M (tv). Its value at point x and time t is thus M (tv)(x), which we'll find more convenient to denote by where M is the magnetization field in the reference state (the one that appears in the law B = µ 0 (H + M )). This being settled, what we look for is the time-derivative ∂ t M v at t = 0, and this for all v.

3.2.
Magnetostriction characteristic and 'frame indifference'. Let us introduce a function m(x, β), vector-valued, with β a linear map of dimension 3, positive definite, acting on vectors anchored at point x. The interpretation is as follows. First (with symbol 1 for the identity map), m(x, 1) is the magnetization M (x) in the reference state. Next, think of a small sphere of radius centered at x, "small" meaning that the average magnetization inside it is equal to m(x, 1) up to terms in and higher. The image of the sphere under β is a volume with the same material content, but a different shape, hence a different average magnetization, whose value is what we denote here by m(x, β). (The name distortion for β, and the symbol β itself, come from [6].) The role of m is thus to assign to a point x a magnetization that may depend not only on the nature of the matter occupying this point, but also on the distortion of the matter around it. (The fact that distortions at neighboring points may not be compatible is irrelevant.) Applying polar decomposition to β, one has β = rs, where s is a symmetric 3 × 3 matrix, positive definite (called strain, or deformation) and r a rotation. Consider a chunk of matter around point x, with characteristic m, under strain s, so it bears a magnetization M = m(x, s). Turn it by r, without modifying the strain. The magnetization, which must follow the rotation, becomes rM , so one has-one must have: r m(x, s) = m(x, rs).
This kind of relation is called 'conjugation' in mathematics: The characteristics m for a piece of matter and for the same piece r-rotated are conjugate by the rotation r. Property (8) is much less commonplace than it may seem: Replace the group of rotations by another subgroup, isomorphic to SO (3), of the linear group, and it's no more true that ρ m(x, s) = m(x, ρs) for all ρ of the new group. (The same observation can be made about many constitutive laws. Explaining why this is often called 'material frame indifference', or 'objectivity' ( [8], p. 8, [12], p. 140) would make a long story.) An important consequence of (8) is that the knowledge of m(x, s) for symmetric distortions s (in coordinates, s ji = s ij ) suffices to know m(x, β) for all distortions. We shall take advantage of that in a moment. For clarity, denote bym the restriction of m to the set of such symmetric distortions. We shall have use for the partial derivative ∂ sm (x, s), which is a symmetric vector-valued 2-tensor.
3.3. The derivative ∂ u M (u)| u=0 . Let's return to the virtual motion x → x+tv(x). The particle located at x + tv(x) at time t was at x at time 0, and has undergone the distortion β v (t) = 1 + t∇v(x) in this process, so M v (t, x + tv(x)) = m(x, β v (t)).
By polar decomposition, β v (t) = r v (t) s v (t), and since t is meant to go to zero, we may use the approximations where ω v = (rot v)/2 and ∇ sym v is the symmetrized gradient, (∂ j v i + ∂ i v j )/2 in Euclidean coordinates. (The compound ω v × denotes a linear map operating on 3-dimensional vectors.) Now, for t small enough, and neglecting higher order terms, (This relation is called the 'Lin constraint' in some parts of the literature ( [5], [9]).) Differentiating in chain with respect to t and letting t go to 0, one obtains where " : " is ad-hoc notation here for this process of chain derivation: Being a derivative with respect to s, ∂ sm has a slot that can be filled by an object of the same type as s, such as ∇ sym v. Remark that ∂ sm is evaluated for s = 1, in the reference configuration.
(It may help at this stage to momentarily return, as a kind of reality check, to Euclidean coordinates. Using Einstein convention, v · ∇M is v j ∂ j M i , a vector, ∇ sym v is the 2-tensor (∂ j v i + ∂ i v j )/2, of the same type as s, and ∂ sm is a vectorvalued 2-tensor, symmetric, whose entries we denote bym k i j . Then ∂ sm : According to (9), which becomes, after an integration by parts over all space, To streamline the notation a little, introduce the 'magnetostrictive stress tensor' to be taken, let's stress that again, in the sense of distributions. We shall refer to the last term in (11) as 'magnetostrictive force'. Let us sum up. Assuming the law B = µ 0 (H + M ) all over (with M = 0 at places, possibly), the force exerted by the magnetic compartment of the coupled system on its elastic compartment appears as being the sum of four contributions: • The Laplace force, (rot H) × B, • The inhomogeneity force, −∇M · B, • The anisotropy force, − 1 / 2 rot(H × B), • The magnetostrictive force, div(σ ms (B)).
The terminology thus proposed can be justified as follows. Laplace force we met earlier, cf. (7). 'Inhomogeneity' refers to the fact that this term would vanish if the magnetization M was spatially uniform. This of course cannot happen, but we are warned that discontinuities of M across material interfaces will result in forces, concentrated at these boundaries. Isotropy, requiring parallel B and H, would result from the invariance of magnetic energy density under rotations, which cannot hold, as a rule, at places where M = 0. (See [3] for an example of anisotropy due to the law B = µH when µ is not a scalar but a 3 × 3 matrix. Then H × B = 0. ) Last, the term div(σ ms (B)) will vanish where and when ∂ sm = 0, that is, when the B-H law does not depend on local strain. Such a dependence seems to be the hallmark of 'magnetostriction', a word used with various meanings in the literature, but one than might benefit from a precise definition. This partition of the force distribution into four components is generic: They appear (some of them null, possibly) whatever the B-H law [3]. Moreover, Laplace force and anisotropy force always take the above form, and a stress tensor similar to σ ms will pop up in case of magnetostriction.
Our next task: Find what to do with expressions such as rot(H × B) or ∇M · B, which make no sense at points belonging to material interfaces, because of discontinuities of M , B and H. Enter a new actor, the 'Maxwell tensor'.

The classical Maxwell tensor. This is the 2-tensor (with M standing for Maxwell, not for magnetization!)
. It's one of the workhorses of Electrical Engineering, used mainly to compute the total electromagnetic force exerted on a definite volume of matter.
To describe this procedure, let us assume the simplest case where B = µ 0 H all over (no magnets, in particular). Neglect electrostatic phenomena, so the only force is the Laplace one, whose density is (rot H)×B. We assume that a preliminary computation has resolved the field equations, hence, at some instant, a 'magnetic pair' (B, H) such that B = µ 0 H and div B = 0. We want the total force on the matter contained in some bounded domain D with boundary S, that is to say, the integral F = ∫ D (rot H) × B. Now, , and this last term is null.
So, by Gauss's theorem, where n is the outgoing unit normal on surface S: The Maxwell tensor allows one to replace a volume integral by a surface integral, a priori a less expensive task, and with a smoother integrand, namely M T · n, than (rot H) × B. Indeed, using (12), one has M T · n = HB n − 1 / 2 (H n B n + H τ · B τ )n, where H n and H τ are the normal and tangential parts of H on S, with the same convention for B. Since H = ν 0 B, this simplifies as M T · n = ν 0 B τ B n + 1 / 2 ν 0 (|B n | 2 − |B τ | 2 )n, with a separation between tangential part and normal part of M T · n which may be welcome in some finite-element computations (cf., e.g., [13]). One must be attentive to choose surface S in such a way that, for almost all its points x, the traces B n and H τ are well defined. The standard recommendation, "take for S a surface that lies entirely in the air and encloses the piece you are interested in", guarantees this, but is way too restrictive: S can traverse matter, as we shall see, provided one uses the 'right' Maxwell-like tensor. Thus described, Maxwell's tensor appears as just an auxiliary, a tool, in force computations. The literature, however, suggests otherwise, for many authors seem to consider this object as more fundamental than that: To the point where one could derive the force field from the Maxwell tensor, whose expression would thus take place among the basic postulates of Electromagnetics, instead of being derived from them. We take the opposite stand: First, derive a formula for the force field ((11), for instance), then find a Maxwell-like tensor of which this force, a priori a rather unsmooth distribution, would be the divergence.

4.2.
A Maxwell-like tensor for the case B = µ 0 (H + M ). Consider which reduces to the M T of (12) when M = 0 all over. This time, not yet the force, but close! Indeed, a comparison with (11) shows what is lacking: The magnetostrictive force, which doesn't make a problem since we already know it as the divergence of a 2-tensor (the σ ms (B) of (11)), and the rot(H × B) term, for which we use the easily proven identity rot(H × B) = div(H ⊗ B − B ⊗ H) with the appropriate −1/2 factor. All in all, Proposition 3. In the B = µ 0 (H + M ) case, the total force field (Laplace force plus magnetic force of (11)) is where M + T is the following extension of the classical Maxwell tensor: Notice the symmetry of M + T (in coordinates, M + T ij = M + T ji ). One may feel frustrated by the limitation to B-H laws of the form B = µ 0 (H + M ), but a close examination of the term 1 / 2 ν 0 |B| 2 δ in (15) will suggest a generalization: One has 1 / 2 ν 0 |B| 2 = 1 / 2 µ 0 |H + M | 2 , which is what is called the coenergy density, ϕ(H) in the notation of [3]. (To check this, remark that the derivative in B not used in the present paper, but for which one may refer to [3]. To sum up: Starting from an expression such as the left-hand side of (14) for the force, one may find a 2-tensor, M + T in this case, whose divergence is this expression, hence a way to compute the total force on the volume enclosed by a surface S, provided the trace M + T · n on S is well defined. This forbids to include in it a part of a material interface, but leaves enough leeway to, for instance, "let S traverse iron" [10], thus avoiding the recourse to hazardous techniques such as the creation of 'virtual airgaps' [4].

4.3.
The pillbox trick. Assume we know the force field as a distribution, as in (14). This doesn't tell us yet the density of force on a material interface S, because of discontinuities of the magnetization M -and as a consequence, of the fields B and H-across S. The pillbox trick allows one to compute this density.
Equip S with a field of unit normals n. (Doing this locally, near a point of interest x ∈ S, for instance, will be enough.) The pillbox is a volume Σ containing a part Σ ∩ S of S and all points x + λn, with x ∈ Σ ∩ S and |λ| ≤ l/2. In other words, Σ is an extrusion of its own equatorial section Σ ∩ S (Fig. 1). The thickness l of the pillbox, taken small enough to avoid the crossing of normals, is destined to tend to 0.
We have a Maxwell-like tensor (the M + T of (15)) which, integrated over the surface of the pillbox, will give the total force over it, ∫ ∂Σ M + T · n Σ . If now we let the height l of the box go to zero, the limit of the integral will be what we are looking for, the force exerted on Σ ∩ S. The passage to the limit will rely on the following (fairly obvious) 'pillbox Lemmas':  The objective in what follows is to express the interface force in terms of such well-defined quantities.
5. Forces at material interfaces: Results. It will be convenient to separate the force field (14) into three components, namely (rot H) × B−∇M · B, the anisotropy force − 1 / 2 rot (H × B), and the magnetostrictive force div(σ ms (B)).

5.1.
Interface force due to (rot H) × B−∇M · B. This term's contribution to the force over the pillbox is 5.2. Interface force due to − 1 / 2 rot(H × B). The part affecting the pillbox is which thanks to the second pillbox Lemma tends to 1 5.3. Interface force due to magnetostriction. This force is just the jump of σ ms (B) · n. It would be worth a deeper study, but that would require more information on the characteristic m.
6. Conclusion. We were able to determine the force inside a magnet (alone, or plunged into the field created by a given current, our J s ), in the form of a vectorvalued distribution. The latter was presented as the sum of a standard vector field and of Dirac-like vector fields borne by material interfaces, for which we could give analytical formulas. Let us conclude with a few remarks, which may let one perceive what remains to be done: • It should be clear that knowing the electromagnetic field is not enough to determine the force: The latter also depends on the characteristic m that connects magnetization and matter deformation.
• Anyway, the very notion of magnetic force is not so limpid when one deals with coupled problems. The definition of what is force depends on where one places the boundary, in modelling, between the magnetic compartment and the elastic compartment of the system. We deliberately wrote Eq. 5 in a way that should make this clear: The source current J s and the Amperians rot M are interchangeable to some extent. Exploiting this possibility would change the magnetic force, without of course changing the solution of the coupled problem.
• Our treatment of magnetostriction was minimal, and insufficient in many respects. Instead of relying on the characteristic m, one might describe the dependence of M on strain by the enunciation of explicit laws: For instance, request that the induction flux due to M through a small material surface does not vary with deformation. Alternatively, enforce this rule about line integrals of amperians, rather than induction fluxes. Force formulas based on both these laws were derived