    # The magnetic field from a homogeneously magnetized cylindrical tile

The magnetic field of a homogeneously magnetized cylindrical tile geometry, i.e. an angular section of a finite hollow cylinder, is found. The field is expressed as the product between a tensor field describing the geometrical part of the problem and a column vector holding the magnetization of the tile. Outside the tile, the tensor is identical to the demagnetization tensor. We find that four components of the tensor, N_xy,N_xz,N_yz and N_zy, can be expressed fully analytically, while the five remaining components, N_xx,N_yx,N_yy,N_zx and N_zz, contain integrals that have to be evaluated numerically. When evaluated numerically the tensor is symmetric. A comparison between the found solution, implemented in the open source magnetic framework MagTense, and a finite element calculation of the magnetic flux density of a cylindrical tile shows excellent agreement.

## Authors

##### This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

## 1 Introduction

A magnetic field is used in a very large variety of applications. Calculating the magnetic field generated by a homogeneously magnetized magnet of a specific shape and size is relevant for several applications, including permanent magnet motors and generators, magnetic bearings and nuclear magnetic resonance (NMR) systems. The magnetic field generated by a permanent magnet or a system of these is often calculated using a finite element framework (FEM) approach, which is limited in precision and speed by the fineness of the finite element mesh used in the calculation. However, the magnetic field generated by a few specific magnet geometries can be expressed analytically if the magnet is uniformly magnetized. This is for example the case for an ellipsoid or a cylinder of infinite lengthOsborn_1945 . The magnetic field is often calculated through the use of the demagnetization tensor, which expresses the relation between the magnetic field and the magnetization

 H=Happl−Nd⋅M (1)

where is the magnetic field, is an externally applied field, is the magnetization and is the demagnetization tensor field, which is a rank-2 symmetric tensor Moskowitz_1966 .

The magnetic field anywhere in space can be calculated analytically through the demagnetization tensor field for a few specific geometries with uniform magnetization, e.g. for a rectangular prism Joseph_1965 ; Smith_2010 or a hollow sphere Prat-Camps2016

. If analytical expressions for the tensor are not known, the internal field within the magnet can be estimated using an average demagnetization factor

Aharoni_1998 ; Joseph_1966 , which can also be computed for non-solid magnetic samples, i.e. powder samples Breit_1922 ; Bleaney_1941 ; Bjoerk_2013 ; Arzbacher_2015 .

In this work we consider the magnetic field generated by a solid cylindrical tile. The cylindrical tile is an angular section of a hollow finite cylinder, as shown in Fig. 1. This geometry is important in several types of applications, for example as segments of a permanent magnet magnetic resonance imaging (MRI) scanner, as magnetic pieces of a generator/motor or as segmented pieces of a Halbach cylinder. For the cylindrical tile magnet geometry there has been a previous attempt to calculate the magnetic field Ravaud_2009 , but no closed-form analytical solution could be obtained for any of the components. The derived expressions for the components of the magnetic field all included definite integrals that have to be evaluated numerically. The cylindrical tile geometry can be seen as a continuation of a previous work, which calculated the magnetic field from a 3D permanent magnetic ring Ravaud_2008 . The magnetic field generated by a whole cylinder is also known Kraus_1973 ; Chen_1991 ; Chen_2006 ; Caciagli_2018 . The magnetic field generated by the cylindrical tile has also been considered in the less complex case of a cylinder of infinite length, i.e. a two dimensional system. In two dimensions, the magnetic field produced by a cylindrical tile has been studied extensively, and the field has been calculated numerous times, with focus on different applications such as magnetic gears Lubin_2010 , motors with Halbach arrays Markovic_2009 , segmented two dimensional Halbach cylinder Shi_2012 and brushless permanent magnet motors Zhu_1993 .

In this work we will derive analytical expressions of the components of the magnetic flux density formulated as a tensor field for a cylindrical tile, for as many of the components as possible. The magnetic field can then easily be computed from the magnetic flux density and the magnetization. The computer code of this model (written in Fortran with a Matlab interface) is publicly available as a part of the MagTense code MagTense .

## 2 The magnetic field of a cylindrical tile

We consider a cylindrical tile in 3D composed of six surfaces. The cylindrical tile is completely specified by six parameters, namely the angular span of the tile from to , the height of the cylinder from to , and the internal, , and external radii of the cylinder, , as illustrated in Fig. 1. We assume the cylinder tile to be homogeneously magnetized with the magnetization . Our objective is then to find the magnetic flux density, , at any point , be it inside or outside the cylinder tile. Once the magnetic flux density is known, the magnetic field can easily be calculated from , where the vacuum permeability is denoted .

We are only interested in the magnetostatic limit and so no free volume currents are flowing. Following Griffiths we may write the magnetic vector potential assuming the Coulomb gauge, , at the position as:

 A(\mathbcalr)=μ04π∫V′∇′×M(\mathbcalr′)|\mathbcalr−\mathbcalr′|dV′+μ04π∫SM(\mathbcalr′)×^n(\mathbcalr′)|\mathbcalr−\mathbcalr′|da′. (2)

Note that there are two sets of coordinates. The coordinates marked with a are the coordinates of the structure that creates the magnetic field, whereas the non-marked coordinates are to the point at which the field is evaluated. It is seen that the first term in the equation is the zero-vector since the magnetization is assumed constant. This also means that the marked coordinates thus refer to the surface of the object that creates the magnetic field. The unit vector normal to the surface with surface element is .

Given we now have:

 (3)

Note that the curl-operator is in non-marked coordinates. As the curl-operator and the integration are with respect to different sets of coordinates, the order of these are interchangable.

We now consider cylindrical coordinates, where the unit vectors can be expressed in the Cartesian base as:

 ^n1 = ^r=(cosθ′,sinθ′,0) ^n2 = ^θ=(−sinθ′,cosθ′,0) ^n3 = ^z=(0,0,1), (4)

noting that they are given in marked coordinates. Note however that the base coordinate system is the same regardless of whether marked coordinates are used or not, and therefore the mark on e.g. or on will be omitted.

With these the integral can be expressed as

 B(\mathbcalr)=μ04π∇×∫SD(\mathbcalr−\mathbcalr′)M(\mathbcalr′)×^n(\mathbcalr′)da′. (5)

with

 D = 1∣∣\mathbcalr−\mathbcalr′∣∣=1√(x−r′cosθ′)2+(y−r′sinθ′)2+(z−z′)2 r′ = √x′2+y′2 (6)

The components of the magnetic flux density are computed by integrating along the entire closed surface of the cylindrical tile. Here the six surfaces of the cylindrical tile are most easily integrated in cylindrical coordinates. In this orthonormal basis, the cylindrical tile is fully symmetric with respect to all axes, which reduces the number of integrals to three.

### 2.1 Surface 1 - the arc surface

We first consider the surface for which is constant, i.e. the arc surface. Here attains the values or . We consider the surface for which the normal points outwards, i.e. the outer arc surface for which . In order to find the flux density we need to integrate the un-marked curl of the surface current density:

 ∇×(D(−^xMzsinθ′+^yMzcosθ′+^z(Mxsinθ′−Mycosθ′))) (7) = ^x((Mxsinθ′−Mycosθ′)∂D∂y−Mzcosθ′∂D∂z) − ^y((Mxsinθ′−Mycosθ′)∂D∂x+Mzsinθ′∂D∂z) + ^z(Mzcosθ′∂D∂x+Mzsinθ′∂D∂y),

The flux density for surface 1 then becomes:

 B1x(\mathbcalr) = μ0r24π∫z2z1∫θ2θ1(Mxsinθ′−Mycosθ′)∂D∂y −Mzcosθ′∂D∂zdz′dθ′∣∣∣r′=r2 B1y(\mathbcalr) = −μ0r24π∫z2z1∫θ2θ1(Mxsinθ′−Mycosθ′)∂D∂x +Mzsinθ′∂D∂zdz′dθ′∣∣∣r′=r2 B1z(\mathbcalr) = μ0r24π∫z2z1∫θ2θ1Mz(cosθ′∂D∂x+sinθ′∂D∂y)dz′dθ′∣∣∣r′=r2. (8)

For the inner arc surface, the normal points inwards and a sign change is required as well as changing where is evaluated from to .

### 2.2 Surface 2 - the vertical surface

The second surface is defined for constant and thus lies in the -plane. We consider the surface normal along the positive azimuthal unit-vector, , and obtain:

 ∇×(D(−^xMzcosθ′−^yMzsinθ′ (9) +^z(Mxcosθ′+Mysinθ′))) = ^x((Mxcosθ′+Mysinθ′)∂D∂y+Mzsinθ′∂D∂z) − ^y((Mxcosθ′+Mysinθ′)∂D∂x+Mzcosθ′∂D∂z) + ^z(Mzcosθ′∂D∂y−Mzsinθ′∂D∂x),

which gives the following contribution to the flux density:

 B2x(\mathbcalr) = +Mzsinθ′∂D∂zdr′dz′∣∣∣θ′=θ2 B2y(\mathbcalr) = −μ04π∫r2r1∫z2z1(Mxcosθ′+Mysinθ′)∂D∂x +Mzcosθ′∂D∂zdr′dz′∣∣∣θ′=θ2 B2z(\mathbcalr) = μ04π∫r2r1∫z2z1Mz(cosθ′∂D∂y −sinθ′∂D∂x)dr′dz′∣∣∣θ′=θ2

The surface with the normal anti-parallel to the azimuthal unit-vector is found by a change of sign on the flux density and replacing where is evaluated from to .

### 2.3 Surface 3 - the horizontal surface

Surface 3 is defined for constant and lies in the -plane. Considering the normal parallel to the positive direction, i.e. we get:

 ∇×(D(^xMy−^yMx)) = ^xMx∂D∂z+^yMy∂D∂z−^z(Mx∂D∂x+My∂D∂y),

 B3x(\mathbcalr) = μ04π∫r2r1∫θ2θ1Mx∂D∂zr′dr′dθ′∣∣∣z′=z2 B3y(\mathbcalr) = μ04π∫r2r1∫θ2θ1My∂D∂zr′dr′dθ′∣∣∣z′=z2 B3z(\mathbcalr) = −μ04π∫r2r1∫θ2θ1(Mx∂D∂x+My∂D∂y)r′dr′dθ′∣∣∣z′=z2

The surface with its normal vector anti-parallel to the direction requires a change of sign and replacing where is evaluated from to .

### 2.4 Tensor components

From the above calculations it is seen that the local flux density at the point from a homogeneously magnetized cylindrical tile can be written as the product between a tensor and the magnetization vector:

 B(\mathbcalr)=μ04πN(\mathbcalr−\mathbcalrp)⋅M(\mathbcalrp), (11)

We explicitly remark that the magnetization vector has been assumed to be constant. The center of the tile is denoted .

Outside the tile, the magnetic field is equal to the magnetic flux density given above with a factor of . Therefore the tensor is identical, with a factor of , to the demagnetization tensor defined in Eq. 1 outside the tile. Within the tile, the demagnetization tensor as defined in Eq. 1 and the magnetic field can be obtained from the relation , as the above expression give both within and outside the tile.

The total flux density may be found by summing the contributions from all the surfaces of the cylindrical tile as given above. We introduce the following nomenclature for simplicity:

 A(r′) = ∫z2z1∫θ2θ1sinθ′∂D∂ydz′dθ′ B(θ′) = ∫r2r1∫z2z1∂D∂ydr′dz′ C(z′) = ∫r2r1∫θ2θ1∂D∂zr′dr′dθ′ D(r′) = ∫z2z1∫θ2θ1cosθ′∂D∂ydz′dθ′ E(r′) = ∫z2z1∫θ2θ1cosθ′∂D∂zdz′dθ′ F(θ′) = ∫r2r1∫z2z1∂D∂zdr′dz′ G(r′) = ∫z2z1∫θ2θ1sinθ′∂D∂xdz′dθ′ H(θ′) = ∫r2r1∫z2z1∂D∂xdr′dz′ I(r′) = ∫z2z1∫θ2θ1cosθ′∂D∂xdz′dθ′ J(r′) = ∫z2z1∫θ2θ1sinθ′∂D∂zdz′dθ′ K(z′) = ∫r2r1∫θ2θ1∂D∂xr′dr′dθ′ L(z′) = ∫r2r1∫θ2θ1∂D∂yr′dr′dθ′

The integrals above are all functions of the unmarked coordinates as well as a constant marked coordinate. In order to shorten the notation in the following, the integrals are given as a function of the constant marked coordinate only.

Using the integral expressions, the components of the tensor, , are given by

 Nxx = r2A(r2)−r1A(r1)+cosθ2B(θ2)−cosθ1B(θ1) +C(z2)−C(z1) Nxy = r1D(r1)−r2D(r2)+sinθ2B(θ2)−sinθ1B(θ1) Nxz = r1E(r1)−r2E(r2)+sinθ2F(θ2)−sinθ1F(θ1) Nyx = r1G(r1)−r2G(r2)+cosθ1H(θ1)−cosθ2H(θ2) Nyy = r2I(r2)−r1I(r1)+sinθ1H(θ1)−sinθ2H(θ2) +C(z2)−C(z1) Nyz = r1J(r1)−r2J(r2)+cosθ1F(θ1)−cosθ2F(θ2) Nzx = K(z1)−K(z2) Nzy = L(z1)−L(z2) Nzz = r2A(r2)−r1A(r1)+r2I(r2)−r1I(r1)+cosθ2B(θ2) (12) −cosθ1B(θ1)+sinθ1H(θ1)−sinθ2H(θ2)

### 2.5 Evaluating the integrals

The integrals given above can more easily be evaluated by introducing a rotation and translation trick following Varga1998 . The goal is to express the integrals as a function of and the primed coordinates only, thus letting .

Along the axis this is achieved by changing the integration limits for all integrals over by subtraction of the value of the coordinate. For the -coordinate we rotate the point of interest about the axis, so that it lies on the -axis. The rotation is defined by the angle . The integration limits of the integrals over are changed by subtracting this angle. Finally, the magnetic flux density should be rotated back about the axis with the angle .

The evaluated expressions for the components of the tensor field are given below. Some integrals have been evaluated using Rubi, the Rule-based Integrator Rich2018 . The integral expressions of the components of the tensor field can most easily be written by introducing a few helper-functions. These are defined as follows:

 A(r,x,θ,z) = √r2−2xrcosθ+x2+z2 B(x,θ,z) = x2(cos2θ−1)−z2 C(r,x) = 4rx(r+x)2 D(θ) = cosθ2 E(r,x,z) = 2√rx(r+x)2+z2 F± = A(r,x,θ,zs)±r√x2+z2s

The integral expressions of the components of the tensor field given below are given as indefinite integrals to compact the notation. As mentioned previously, the point of interest is defined as . Thus, the appropriate integration limits should be inserted in , respectively. Subscript indicates a constant value on a surface, e.g. is the value of on a surface, i.e. or . Incomplete elliptic integrals of the first, second and third kind are denoted , and , respectively. Using these functions, the integrals can be expressed as given below, respecting the singularities in the solutions given in the appendix.

 A(rs,θ′,z′) = 1rs12z′x2√(rs+x)2+z′2sgn(sinθ′2)((rs−x)2Π{D(θ′),C(rs,x),E(rs,x,z′)} B(r′,θs,z′) = atan(z′(r−xcosθs)cscθsxA(r,x,θs,z′))cosθs−atanh(z′A(r,x,θs,z′))sinθs C(r′,θ′,zs) = D(rs,θ′,z′) = 14x2((r2s+x2)(ln{A(rs,x,θ′,z′)−z′}−ln{A(rs,x,θ′,z′)+z′})−2z′A(rs,x,θ′,z′)) E(rs,θ′,z′) = −sgn(sinθ′/2)rsx√(rs+x)2+z′2(((rs+x)2+z′2)E{D(θ′),E(rs,x,z′)} −(r2s+x2+z′2)F{D(θ′),E(rs,x,z′)}) F(r′,θs,z′) = −ln{r−xcosθs+A(r,x,θs,z′)} G(rs,θ′,z′) = H(r′,θs,z′) = ∫rx(1−cos2θs)+z′2cosθsA(r,x,θs,z′)B(x,θs,z′)dz′ I(rs,θ′,z′) = −1rssgn(sinθ′/2)z′2x2(rs+x)√(rs+x)2+z′2((rs−x)(r2s+x2)Π{D(θ′),C(rs,x),E(rs,x,z′)} J(rs,θ′,z′) = −1rs√r2s−2rsxcosθ′+x2+z′2x K(r′,θ′,zs) = ∫cosθ′ln(r−xcosθ′+A(r,x,θ′,zs))dθ′ −∫cosθ′(2x2rcos2θ′−cosθ′(x3−xz2s)−r(x2+z2s))A(r,x,θ′,zs)B(x,θ′,zs)dθ′ +∫x(rxcosθ′−x2−z2s)A(r,x,θ′,zs)B(x,θ′,zs)dθ′ L(r′,θ′,zs) = +√x2+z2s(−(r2+x2+z2s)ln{r(r−xcosθ′+A(r,x,θ′,zs)} +(r2−2rxcosθ′+x2+z2s)ln{r−xcosθ′+A(r,x,θ′,zs)} −2rA(r,x,θ′,zs)−(r2+x2+z2s)ln2+2rxcosθ′−(r2+x2+z2s)))

The integrals and cannot fully be evaluated analytically and have to partially be computed numerically. We implement this by using the QuadPack numerical integration library QuadPack , specifically the qags function. The remaining integrals can be expressed analytically. This means that the following components of the tensor field can be evaluated analytically, and . The remaining components have to be evaluated numerically, which are and . Interestingly, as the tensor field is symmetric, as is the demagnetization tensor, it is curious that and have analytical expressions while and have to be evaluated numerically. However, when evaluated, the expressions have the same numerical value, and the tensor is thus symmetric.

## 3 Verification

The magnetic field as computed using the tensor field formulation is verified against a finite element method computation of the magnetic flux density from a cylindrical tile. The magnetic field is computed using the finite element framework Comsol. The equation solved in the FEM framework is the magnetic scalar potential equation

 −∇⋅(μ0∇Vm)=−∇⋅(μ0M) . (13)

The magnetic field is then calculated as . For the FEM model, a highly refined finite element mesh is used to ensure a high precision in the computed field, and a sufficiently large volume is modelled to make the boundary conditions not influence the generated magnetic field.

In the first example we consider a cylindrical tile similar to the one shown in Fig. 1. It has the following geometry: mm, mm, , , mm and mm. This is similar to specifying the center of the tile with and tile dimensions of . The magnetization is specified by the vector
T. The components of the magnetic field is calculated along a line from the point to , thus passing almost through the center of the tile and is shown in Fig. 2. As can be seen from the figure, there is an excellent agreement between the model framework presented here and the FEM model. This example is part of the verification examples for MagTense and is available online MagTense . Figure 2: The components of the magnetic field along a line from the point (x,y,z)=(2,−1,−3) to (8,5,3) as calculated using the components of the tensor field given above in the MagTense framework and as calculated using the FEM framework Comsol, for the cylindrical tile geometry specified in the text.

As a second example, we consider a cylindrical tile with the following geometry: mm, mm, , , mm, mm, with the cylinder central axis offset to mm, mm. The magnetization is specified by the vector
T, which corresponds to a direction specified by the spherical azimuthal angle and the spherical polar angle and a magnitude of . The norm of the magnetic field is calculated along each of the Cartesian axes centered on the tile and is shown in Fig. 3. As can be seen from the figure, there is again an excellent agreement between the model framework presented here and the FEM model. The individual components of the field, which are not shown in the figure due to brevity, show an equally excellent agreement. This example is part of the verification examples for MagTense and is available online MagTense . Figure 3: The norm of the magnetic field as calculated using the components of the tensor field given above in MagTense and as calculated using the FEM framework Comsol, for the cylindrical tile geometry specified in the text along each of the Cartesian axes through the center of the tile.

Regarding the computational speed of calculating the magnetic field, in the FEM model this crucially depends on the number of elements used, as do the precision. For example one the FEM solution has a fine mesh with tetrahedral elements. On an Intel Xeon W-3235 with 128 GB of memory the computation time for this FEM model was 47 seconds (disregarding the meshing time, which was about the same time). For MagTense the computation time was 0.126 seconds. For example two, an almost equal mesh with also tetrahedral element was used. Here the computational time was 43 seconds, while the MagTense computation time was 0.39 seconds, as three time as many points are evaluated in example two compared to example one. The computational speed of MagTense is so fast because in MagTense the field is only evaluated in a set of specified points, whereas it in a FEM model has to be evaluated in the entire simulation volume.

## 4 Conclusion

We have calculated the magnetic flux density components formulated as a tensor field for a cylindrical tile geometry. The components of the tensor field involved a number of integrals, which meant that not all components could be evaluated analytically. In the end, five components and have to be evaluated numerically, while and could be evaluated purely analytically. As the tensor field is symmetric, only six of these components needs to be evaluated in actual computations. There was an excellent agreement between the magnetic flux density of a cylindrical tile calculated using the tensor field and the flux density computed using an finite element approach.

## Acknowledgements

This work was financed partly by the Energy Technology Development and Demonstration Program (EUDP) under the Danish Energy Agency, project no. 64016-0058, partly by the Danish Research Council for Independent Research, Technology and Production Sciences projects no. 7017-00034 and 8022-00038 and partly by the Poul Due Jensen Foundation project on Browns paradox in permanent magnets, project no. 2018-016.

## Appendix A Singularities in the solutions

Some of the above given integrals have singularities for certain values of the variables. This puts limitations on the possible coordinates at which the tensor field can be evaluated. It is obvious that the following cannot be violated:

 r > 0 (14) x ≠ 0 (15) A(r,x,θ′,z) > 0 (16) B(x,θ′,z) ≠ 0 (17) r−xcosθ′+A(r,x,θ′,z) > 0 (18) A(r,x,θ,z)−z > 0 (19) A(r,x,θ,z)+z > 0 (20) F± ≠ ±1 ⇒A(r,x,θ,z)±r≠±√x2+z2, (21)

where is an integer.

#### a.0.1 The A function

It is a requirement that must be real and greater than zero and thus . In the following, we will consider the three variables, , and for any real value of . First, we note that is symmetric in and so the results derived from one of these variables will apply equally to the other.

We consider the above given inequality in the limit where it is an equality and solve for :

 x=rcosθ±√r2cos2θ−(r2+z2). (22)

We see immediately that can only be real valued when for being any integer including zero and conclude that the inequality holds for any value of and therefore also with this requirement.

Considering , we get that , which is to say that for the inequality always holds. However, in the case where we get:

 cosθ≤