On the equivalence of demagnetization tensors as discrete cell size approaches zero in three-dimensional space

Liang, Hao; Yan, Xinqiang

doi:10.1063/5.0226603

The calculation of the demagnetization field is crucial in various disciplines, including magnetic resonance imaging and micromagnetics. A standard method involves discretizing the spatial domain into finite difference cells and using demagnetization tensors to compute the field. Different demagnetization tensors can result in contributions from adjacent cells that do not approach zero, nor do their differences, even as the cell size decreases. This work demonstrates that in three-dimensional space, a specific set of magnetization tensors produces the same total demagnetization field as the Cauchy principal value when the cell size approaches zero. Additionally, we provide a lower bound for the convergence speed, validated through numerical experiments.

I. INTRODUCTION

Calculation of magnetic fields is important in various fields, including micromagnetics and magnetic resonance imaging (MRI). In the field of micromagnetics, numerical simulations often solve the Landau–Lifshitz–Gilbert equation along with its associated partial differential equations.¹ The primary variable is the magnetization vector field, $M$ ⁠. Given the magnetization field, various fields such as exchange, anisotropy, Zeeman, and demagnetization fields are computed,² which, in turn, influence magnetization. The most computationally demanding part is calculating the demagnetization field given the magnetization.

In MRI, the strong magnetic fields produced by scanners create varying magnetizations within different tissue types in the body. These magnetizations result in non-uniformly distributed demagnetization fields. The presence of such magnetic fields inhomogeneities is a major challenge in acquiring high-quality images, especially at higher field strengths where demagnetization fields are also of greater amplitude. Accurately computing these fields is essential for understanding how they vary across individuals with different body sizes, shapes, and anatomical structures.

In the above-mentioned disciplines, a common calculation involves determining the demagnetizing field given the magnetization. A common approach is to discretize the space into a regular grid of cubic cells. Within each cell, the magnetization and other fields are considered constant. To calculate the demagnetization field, an integral over the entire magnetic domain, in finite difference discretization, translates to a summation over all cells as follows:³

H {(r)}^{e x c l u d e c e l l a t r} = - \sum_{r^{'} \neq r} N (r^{'} - r) \cdot M (r^{'}),

(1)

where

N

is the demagnetization tensor. The cell at

r^{'} = r

is excluded from the summation for convenient comparison among different methods. Due to the approximation of the discretization, there are field errors, especially in the case of rapid changes in magnetization. However, the impact of discretization (e.g., cell size) on field errors has not yet been studied based on a comprehensive theoretical analysis. This article aims to fill this gap and provides a comprehensive theoretical explanation of these errors.

Typically, in micromagnetics, the demagnetization tensor is defined based on the interaction energy between cells, assuming each cubic cell is uniformly magnetized.^4–6 This approach is referred to as the uniformly magnetized cube (UMC) method, and the demagnetization tensor of this approach is denoted as $N_{c}$ ⁠. Another approach treats each cell as a point dipole located at its center. This approach is used in both micromagnetics^7,8 and MRI^9–11 owing to its simplicity of calculation for the demagnetization tensor. This approach is referred to as the dipole method, and the demagnetization tensor of this approach is denoted as $N_{d}$ ⁠. Besides the UMC and dipole methods, the cell at $r^{'}$ can be treated as a uniformly magnetized cube; however, the cell at $r$ is treated as a point dipole.¹² This method can benefit from averaging $r^{'}$ over a cell and keeping the calculation relatively simple. We refer to this method as the uniformly magnetized cube-dipole (UMCD) method and its demagnetization tensor as $N_{c d}$ ⁠.

The equivalence of these methods is non-trivial. The demagnetization tensor is proportional to the volume of the small cell (⁠

h^{3}

⁠) and also

- \nabla \nabla (1 / | r |) \sim (r^{2} - 3 z^{2}) / r^{5}

⁠, which is about

h^{- 3}

for cells adjacent to the cell where the field is being evaluated. This means that the demagnetization field contribution from a single nearby cell is NOT small but finite. This situation differs significantly from integration on a smooth function, where each piece is small, and the sum of small pieces is a finite value. In our situation, the contribution of a single cell adjacent to the cell at

r

does not approach zero, nor does the difference between different methods, even as the cell size approaches zero. For example, while

N_{d}

is a good approximation of

N_{c}

for

r^{'}

far from

r

⁠, this is not the case for cells adjacent to the cell at

r

⁠. For example, the component

{(N_{d})}_{z z}

of the tensor

N_{d}

is given by

{(N_{d})}_{z z} (r) = \frac{1}{4 π} \frac{i^{2} + j^{2} - 2 k^{2}}{ρ^{5}},

(2)

where

i = x / h

⁠,

j = y / h

⁠, and

k = z / h

⁠. It should be noted that

M

and

H

have the same dimension, making the demagnetization tensor unitless in the metric system. A single cell at

i = j = 0

and

k = 1

with magnetization of

M_{z}

contributes

2 M_{z} / (4 π)

to the demagnetization field

H_{z}

⁠. Similarly, using other methods, this single cell can contribute finite but different values to the demagnetization field. These differences can exceed ten percent, independent of

h

⁠,^13,14 as illustrated in Fig. 1. Therefore, the conditions under which these methods yield the same result for the demagnetization field require investigation.

FIG. 1.

View large Download slide

Upper panel: demagnetization tensor, $N_{z z}$ ⁠, for the UMC method, dipole method, variant of dipole method with $b = 0.2 ı$ ⁠, the asymptotic expansion for the UMC method up to $n = 6$ ⁠, i.e., $O ({(h / r)}^{7})$ ⁠, and UMCD method. In the plots, $x - x^{'} = y - y^{'} = 0$ ⁠, and $k$ is defined by $z^{'} - z = k h$ ⁠. Lower panel: the relative differences between these methods and the UMC method. Specifically, the relative difference is defined as $| N_{z z} / {(N_{c})}_{z z} - 1 |$ for demagnetization tensor $N$ ⁠.

To our knowledge, the relationship between these methods is not sufficiently discussed. Numerous literature implicitly assumed that the UMC, dipole, and UMCD can yield the correct results, i.e., converge to the same result.^4–11 However, other works hint that these methods cannot converge to the same result, due to significant discrepancy between the demagnetization tensors for individual cells.^13–16 Additionally, some works treat the result of the dipole method as equivalent to that of the Cauchy principal value.^9–11

One main goal of this article is to prove that all these methods consistently produce results at the limit of $h$ approaching zero in three-dimensional space using cubic cells. The equivalence is not on a cell-by-cell basis but on the total field for sufficiently smooth magnetization.

This article is structured as follows: Sec. II reviews the demagnetization tensors, presents the proof of our statement, provides a lower bound on the convergence speed, and outlines the implementation using FFT. Section III demonstrates the validation of the statement through numerical experiments. Section IV discusses the extension and limitation of the current work. Section V offers concluding remarks.

II. THEORETICAL ANALYSIS

A. Cauchy principal value

The demagnetization field for open boundary conditions is given by^2,3

H = - \nabla ψ (r),

(3)

ψ (r) = - \frac{1}{4 π} \int \frac{\nabla^{'} \cdot M (r^{'})}{| r - r^{'} |} d^{3} r^{'},

(4)

where

\nabla

is applied to the functions with respect to

r

⁠, and

\nabla^{'}

is applied to the functions with respect to

r^{'}

⁠,

ψ

is the scalar potential. Although the function

1 / | r - r^{'} |

with respect to

r^{'}

diverges at

r

⁠, the divergence is slow, and the integration of Eq. exists as an improper integral. After some mathematical processing, Eq. can be formally written as

H (r) = - \frac{1}{4 π} \int M (r^{'}) \cdot \nabla \nabla^{'} \frac{1}{| r - r^{'} |} d^{3} r^{'} .

(5)

However, the term

\nabla \nabla^{'} (1 / | r - r^{'} |)

as a function of

r^{'}

diverges rapidly at

r

⁠, and the result of Eq. depends on how the singularity is treated. One approach is to use the Cauchy principal value, which we will briefly review. We decompose

M (r^{'})

into two functions,

u^{- S} (r^{'}) M (r^{'})

and

u^{S} (r^{'}) M (r^{'})

⁠. Here,

u^{S} = 1

inside a sufficiently small sphere

S

⁠, and

u^{S} = 0

elsewhere; conversely,

u^{- S} = 1

outside the sphere and

u^{- S} = 0

elsewhere. This ensures

u^{S} (r) + u^{- S} (r) \equiv 1

⁠. Using this notation,

\int

always means the integral over the whole space, Eqs. and can be written as

H = \frac{1}{4 π} \nabla (\int \frac{\nabla^{'} \cdot u^{- S} (r^{'}) M (r^{'})}{| r - r^{'} |} d^{3} r^{'} + \int \frac{\nabla^{'} \cdot u^{S} (r^{'}) M (r^{'})}{| r - r^{'} |} d^{3} r^{'}) .

(6)

The surface magnetic charge has been implicitly included in Eq. .

Since

u^{- S} = 0

inside of the sphere, we can avoid the singularity issue in evaluating the first part in Eq. . The first part in Eq. can be simplified as follows using the method of integration by parts:

\frac{1}{4 π} \nabla \int (\nabla^{'} \cdot \frac{u^{- S} (r^{'}) M (r^{'})}{| r - r^{'} |} - u^{- S} (r^{'}) M (r^{'}) \cdot \nabla^{'} \frac{1}{| r - r^{'} |} d^{3} r^{'}) .

(7)

The first term in is simply zero by Gauss’s law, noting that

u^{- S} M

inside the sphere and around infinity is zero, respectively. Thus, Eq. becomes

\begin{aligned} - \frac{1}{4 π} \nabla \int u^{- S} (r^{'}) M (r^{'}) \cdot \nabla^{'} \frac{1}{| r^{'} - r |} d^{3} r^{'} \\ = - \frac{1}{4 π} \int_{- S} M (r^{'}) \cdot \nabla \nabla^{'} \frac{1}{| r - r^{'} |} d^{3} r^{'} . \end{aligned}

(8)

We identify this as the Cauchy principal value, denoted as

H_{\,p}

⁠. The integration of the second term of the integrands in is a well-known problem, and the result¹⁷ is

- \frac{1}{3} M

⁠. The macroscopic field

H

can be decomposed into two parts,

H = - \frac{1}{3} M + H_{p} .

(9)

H_{p}

is given by

H_{p} (r) = - \int_{- S} N_{p} (r - r^{'}) \cdot M (r^{'}) d^{3} r^{'},

(10)

where

\begin{aligned} N_{p} (r - r^{'}) & = \frac{1}{4 π} \nabla \nabla^{'} \frac{1}{| r^{'} - r |} \\ = - \frac{1}{4 π} \frac{3 (r^{'} - r) (r^{'} - r) - | r^{'} - r |^{2}}{| r^{'} - r |^{5}} \end{aligned}

(11)

is the so-called

3 \times 3

demagnetization tensor.

{(N_{\,p})}_{i j}

represents the field component in direction

i

at position

r

⁠, generated by the magnetization component in direction

j

of a volume element at position

r^{'}

⁠, where

i

and

j

are labels for axes. The integration of Eq. is operational because the integration domain excludes the singularity point of

N_{\,p}

⁠. It is important to note that the Cauchy principal value at

r

is distinct from the macroscopic magnetization by

- M / 3

⁠. More precisely, it represents the field observed in a spherical cavity generated from the remaining parts of magnetization.

In MRI, the magnetic field with the so-called sphere Lorentz correction, denoted as

B^{'}

is of interest,¹⁸ namely,

B^{'} = B - \frac{2}{3} μ_{0} M = μ_{0} (H + \frac{1}{3} M) .

(12)

Thus, we have

B^{'} = μ_{0} H_{p}

⁠, which is already discussed in the literature.^19,20

B. UMC and UMCD methods

As mentioned above, the UMC refers to the approach where each cubic cell is uniformly magnetized (however, the magnetization is not necessarily equal for different cells). In this method, the field is given by

H_{c} (r) = - \sum_{i, j, k}^{r_{i j k} \neq r} N_{c} (r_{i j k} - r) \cdot M (r_{i j k}),

(13)

r_{i j k} = (i e_{x} + j e_{y} + k e_{z}) h,

(14)

where

i

⁠,

j

⁠, and

k

are indices for 3D grid, and

h

represents the cell size. To make the article more concise, we express it as follows:

H_{c} (r) = - \sum_{r^{'} \neq r} N_{c} (r^{'} - r) \cdot M (r^{'}) .

(15)

The cell at

r^{'} = r

is excluded from the summation for convenient comparison with other methods. The value of the demagnetization tensor is defined as the average demagnetization field in a cell centered at

r

⁠, which is generated by another uniformly magnetized cubic cell (with a unit magnetization) centered at

r^{'}

⁠. Specifically, the demagnetization tensor is given by

N_{c} (r) = \frac{1}{4 π h^{3}} \int_{c e l l 0} d^{3} r^{'} \int_{c e l l r} d^{3} r^{″} \nabla^{'} \nabla^{″} \frac{1}{| r^{″} - r^{'} |} .

(16)

The demagnetization tensor

N_{c}

can be calculated using various methods, including numerical integration,^4,21 exact analytical formulas,^1,3,14,15 or analytical formulas combined with asymptotic expansions.¹⁶ As an example,

N_{z z}

is calculated and shown in Fig. 1.

For the UMCD method, the cell at

r^{'}

is treated as a uniformly magnetized cube, but the cell at

r

is treated as a point dipole. This results in a demagnetization tensor of

N_{c d} (r) = \frac{1}{4 π} \int_{c e l l r} d^{3} r^{″} {(\nabla^{'} \nabla^{″} \frac{1}{| r^{″} - r^{'} |}) |}_{r^{'} = 0} .

(17)

The analytical expression of the demagnetization tensor for the UMCD method^12,15 is simpler than that for the UMC.

N_{zz}

as an example is shown in Fig. 1. Similar to the UMC method, the cell at

r^{'} = r

is excluded from the summation for the UMCD method, for convenient comparison with other methods.

However, for both methods, by noticing that

Tr [\nabla \nabla^{'} (1 / | r^{'} - r |)] = \nabla \cdot (\nabla^{'} (1 / | r^{'} - r |)) = 4 π δ (r^{'} - r)

and considering the symmetries of

x

⁠,

y

⁠, and

z

axes, it is not difficult to demonstrate that demagnetization tensor for cell at

r^{'} = r

is

(1 / 3) δ_{a b}

⁠, where

δ (r - r^{'})

is the Dirac delta function,

δ_{a b}

Kronecker delta,

Tr

denotes the trace of a

3 \times 3

tensor. Therefore, the demagnetization field generated by the cell itself is

- (1 / 3) M

⁠. Thus,

H_{c}^{\, f u l l} \equiv - \sum_{r^{'}} N_{c} (r - r^{'}) \cdot M (r^{'}) = H_{c} (r) - \frac{1}{3} M

(18)

and

H_{c d}^{\, f u l l} \equiv - \sum_{r^{'}} N_{c d} (r - r^{'}) \cdot M (r^{'}) = H_{c d} (r) - \frac{1}{3} M .

(19)

C. Dipole method

As mentioned above, the dipole method replaces a uniformly magnetized cubic cell with a single dipole at the cell’s center. The demagnetization field is expressed by

H_{d} (r) = - \sum_{r^{'} \neq r} N_{d} (r^{'} - r) \cdot M (r^{'}),

(20)

where

N_{d}

represents the demagnetization tensor for the dipole method. This tensor quantifies the demagnetization field at

r

⁠, generated by a dipole with a moment of

h^{3}

at

r^{'}

⁠, specifically,

N_{d} (r) = - \frac{3}{4 π ρ^{5}} (\begin{matrix} i^{2} - ρ^{2} / 3 & i j & i k \\ i j & j^{2} - ρ^{2} / 3 & j k \\ i k & j k & k^{2} - ρ^{2} / 3 \end{matrix}),

(21)

ρ \equiv \sqrt{i^{2} + j^{2} + k^{2}},

(22)

r = (i e_{x} + j e_{y} + k e_{z}) h,

(23)

where

i

⁠,

j

⁠, and

k

are integral indices for 3D grids;

e_{x}

⁠,

e_{y}

⁠, and

e_{z}

are the unit vectors along

x

⁠,

y

⁠, and

z

directions, respectively.

Note that Eq. (21) does not explicitly depend on $h$ ⁠. Therefore, the contribution from a single cell adjacent to the cell at $r$ does NOT approach zero as $h \to 0$ ⁠. It is also noteworthy that since the interaction energy between two uniformly magnetized spheres equals that between two dipoles of the same moment,^3,22,23 the method can also be referred to as the uniformly magnetized sphere method. $r^{'} = r$ must be excluded from the summation because the demagnetization tensor for the dipole method is singular at this point.

D. Asymptotic expansions

The demagnetization tensors for both UMC and UMCD methods have an asymptotic expansion in the form of (see Appendix A)

\sum_{n} N^{(n)} (r) = \sum_{n = 2, 4, \dots} \sum_{s, t} c_{n, s, t} {(\frac{x}{r})}^{s} {(\frac{y}{r})}^{t} {(\frac{z}{r})}^{n - s - t} {(\frac{h}{r})}^{n + 1},

(24)

where

c_{n, s, t}

are coefficients of mathematical constants. It can also be written in the form of

\sum_{n = 2, 4, \dots} \sum_{s, t} c_{n, s, t} {(\frac{i}{ρ})}^{s} {(\frac{j}{ρ})}^{t} {(\frac{k}{ρ})}^{n - s - t} {(\frac{1}{ρ})}^{n + 1},

(25)

where

i

⁠,

j

⁠,

k

⁠,

ρ

are defined by

x / h

⁠,

y / h

⁠,

z / h

⁠, and

| r | / h

⁠, respectively. Again,

h

does not explicitly appear in the formula. Evaluating

c_{n, i, j}

according to Appendix A for both UMC or UMCD methods, we find that: for

n = 2

⁠, this represents the dipole approximation; For

n = 4

⁠, these values are zero for cubic cells, which explains why the demagnetization tensor of dipole method converges to those of the UMC or UMCD fast for distant cells. For

n = 6

⁠, the asymptotic expansions for

N_{c}

⁠, such as the components

N_{x x}^{(6)}

and

N_{x y}^{(6)}

⁠, are given as follows:

\begin{aligned} N_{x x}^{(6)} = \frac{1}{4 π} \cdot \frac{7}{16 ρ^{13}} (2 i^{6} - j^{6} - k^{6} - 15 i^{4} (j^{2} + k^{2}) + 15 i^{2} (j^{4} + k^{4})) \end{aligned}

(26)

and

\begin{aligned} N_{x y}^{(6)} = \frac{1}{4 π} \cdot \frac{7}{16 ρ^{13}} i j (7 i^{4} - 19 i^{2} j^{2} + 7 j^{4} - 13 (i^{2} + j^{2}) k^{2} + 13 k^{4}) . \end{aligned}

(27)

For the UMCD method, the correction of

n = 6

to the demagnetization tensor is only half of that for the UMC method. The asymptotic expansion of the magnetization tensor for the UMC is shown in Fig. 1. The figure shows that the asymptotic expansions for the UMC method up to

n = 6

converge much faster than the dipole approximation for distant cells. However, higher-order corrections increase the error for the cell adjacent to the cell where the field is being evaluated.

E. Near and far field decomposition

We use an analysis method similar to that described in the literature¹⁷ to split the demagnetization field into near and far fields. The near magnetic field originates from a cubic volume

V

of size

L \times L \times L

⁠, while the far field encompasses contributions outside this cubic volume. We can always choose

L

to be macroscopically small but significantly larger than

h

⁠, ensuring that

V

contains exactly

2 N + 1

cubic cells, satisfying

h ≪ 1, L ≪ 1, a n d \frac{L}{h} = 2 N + 1 ≫ 1.

(28)

We will prove that these methods result in the same field under the limit

h \to 0, L \to 0, a n d \frac{L}{h} \to + \infty .

(29)

It is intuitive that these methods yield the same far field, given that

L / h \to \infty

under the limit defined in Eq. . Therefore, it is only necessary to prove these methods produce the same near field under the limit .

F. Near field for uniform magnetization

First, we prove that the near part of the demagnetization field produced by all these methods equals zero in the case of uniform magnetization, i.e., the magnetization is uniform in the near cells. We utilize a similar treatment as in literature.¹⁷ For the Cauchy principal value method, the near field can be expressed as

H_{\, p}^{V} (r) = - \int_{V - S} N_{\, p} (r - r^{'}) \cdot M_{0} (r^{'}) d^{3} r^{'} .

(30)

The superscript

V

indicates that only the contribution from the cubic volume

V

is counted. We can directly verify that

{(N_{\,p})}_{x x} (r - r^{'}) + {(N_{\,p})}_{y y} (r - r^{'}) + {(N_{\,p})}_{z z} (r - r^{'}) = δ (r^{'} - r)

according to , thus

\int_{V - S} ({(N_{\, p})}_{x x} + {(N_{\, p})}_{y y} + {(N_{\, p})}_{z z}) d^{3} r^{'} = 0.

(31)

Since

r

is exactly at the center of

V

⁠, there is the symmetry between

x

⁠,

y

⁠, and

z

axes. Thus, the integral of the three terms in are equal and must be equal to zero, namely,

\int_{V - S} {(N_{\, p})}_{x x} d^{3} r^{'} = \int_{V - S} {(N_{\, p})}_{y y} d^{3} r^{'} = \int_{V - S} {(N_{\, p})}_{z z} d^{3} r^{'} = 0.

(32)

For the dipole method, we need to evaluate the magnetic field from dipoles inside the cube, according to , namely,

H_{d}^{V} (r) = - \sum_{\begin{matrix} _{r^{'} \in V} \\ _{r^{'} \neq r} \end{matrix}} N_{d} (r - r^{'}) \cdot M (r^{'}) .

(33)

We can directly examine

{(N_{d})}_{x x} + {(N_{d})}_{y y} + {(N_{d})}_{z z} = 0

according to , thus,

\sum_{\begin{matrix} _{r^{'} \in V} \\ _{r^{'} \neq r} \end{matrix}}^{N} ({(N_{d})}_{x x} + {(N_{d})}_{y y} + {(N_{d})}_{z z}) = 0.

(34)

Since the symmetry between

x

⁠,

y

⁠, and

z

axes, the summations of each of three terms in are equal to zero, as follows:

\sum_{\begin{matrix} _{r^{'} \in V} \\ _{r^{'} \neq r} \end{matrix}} {(N_{d})}_{x x} = \sum_{\begin{matrix} _{r^{'} \in V} \\ _{r^{'} \neq r} \end{matrix}} {(N_{d})}_{y y} = \sum_{\begin{matrix} _{r^{'} \in V} \\ _{r^{'} \neq r} \end{matrix}} {(N_{d})}_{z z} = 0.

(35)

Due to the parity (mirror) symmetry, the off diagonal elements of the demagnetization tensor are also zero

\int_{V - S} {(N_{\, p})}_{x y} d^{3} r^{'} = \int_{V - S} {(N_{\, p})}_{y z} d^{3} r^{'} = \int_{V - S} {(N_{\, p})}_{x z} d^{3} r^{'} = 0,

(36)

\sum_{\begin{matrix} _{r^{'} \in V} \\ _{r^{'} \neq r} \end{matrix}} {(N_{d})}_{x y} = \sum_{\begin{matrix} _{r^{'} \in V} \\ _{r^{'} \neq r} \end{matrix}} {(N_{d})}_{y z} = \sum_{\begin{matrix} _{r^{'} \in V} \\ _{r^{'} \neq r} \end{matrix}} {(N_{d})}_{x z} = 0.

(37)

For the UMC, UMCD, and their asymptotic expansion methods, the symmetries also exist, allowing similar equations to Eqs. and to hold. Note that the dipole method can be seen as the leading order of UMC and UMCD’s asymptotic expansions. In summary, we have proven that, for all these methods, the near part of the demagnetization field is exactly zero, namely,

H_{\, p}^{V} (r) = H_{c}^{V} (r) = H_{c d}^{V} (r) = H_{d}^{V} (r) = 0.

(38)

Modifications to the demagnetization tensors, while preserving Eq. and symmetries of

x

⁠,

y

⁠, and

z

axes, can also result in . This guides us in creating variation of the dipole method. An example is replacing

ρ = \sqrt{i^{2} + j^{2} + k^{2}}

in the denominator in Eq. with

ρ^{'} = \sqrt{i^{2} + j^{2} + k^{2} + b}

⁠, namely,

N_{d} (r) = - ℜ \frac{3}{4 π ρ^{' 5}} (\begin{matrix} i^{2} - ρ^{2} / 3 & j k & i k \\ j k & j^{2} - ρ^{2} / 3 & j k \\ i k & j k & k^{2} - ρ^{2} / 3 \end{matrix}),

(39)

ρ = \sqrt{i^{2} + j^{2} + k^{2}}, ρ^{'} = \sqrt{i^{2} + j^{2} + k^{2} + b},

(40)

where

b

can be positive real number or imaginary number,

ℜ

⋆ denote the real part of ⋆. The additional term

b

softens the demagnetization tensor, possibly leading to a smoother result where

M

changes rapidly.

G. Near field for nonuniform magnetization

In the previous subsection, it has been established that the near demagnetization field is exactly zero for uniform magnetization. Then, we prove that the near demagnetization field at $r$ approaches zero as $L \to 0$ ⁠, provided that the magnetization is Hölder continuous at $r$ ⁠.

By subtracting a constant magnetization, $M_{0} = M (r)$ ⁠, from $M (r^{'})$ ⁠, the result remains unaffected. Consequently, with this subtraction, ${(M (r^{'}) - M_{0}) |}_{r^{'} = r} = 0$ ⁠, allowing us to potentially address the challenges arising from the discrepancy of the demagnetization tensor around $r$ ⁠.

For the principal integral method, we have

H_{\, p}^{V} (r) = - \int_{V - S} N_{\, p} (r - r^{'}) \cdot (M (r^{'}) - M_{0}) d^{3} r^{'} .

(41)

We assume

M (r^{'})

hold Hölder condition at

r

⁠, namely,

| M (r^{'}) - M_{0} | < K | r^{'} - r |^{α},

(42)

where

K

is a constant independent of

r^{'}

and

α > 0

⁠.

First, we can estimate the magnitude of the demagnetization field for the principal integral method. Since

N_{\,p}

diverges as

r^{- 3}

as

r \to 0

⁠, the integral in spherical coordinates is bounded by

H_{\, p}^{V} \sim \int_{V - S} K r^{α} (1 / r^{3}) r^{2} d r d θ d ϕ \sim K L^{α} .

(43)

In the estimation, we neglect constants on order of one. Equation approaches zero as

L \to 0

⁠. Then, we turn to the dipole method. Similarly, we can estimate the magnitude of the demagnetization field of the dipole method as

\begin{aligned} H_{d}^{V} (r) & \sim \sum_{i, j, k} K {(ρ h)}^{α} (1 / ρ^{3}) \\ \sim \sum_{I = 1}^{N} \sum_{max (| i |, | j |, | k |) = I} K {(I h)}^{α} (1 / I^{3}) \\ \sim \sum_{I = 1}^{N} K {(I h)}^{α} I^{2} / I^{3} \\ \sim \sum_{I = 1}^{N} K h^{α} I^{α - 1} \\ \sim K h^{α} (N^{α} - 1) \\ \sim K L^{α} . \end{aligned}

(44)

Thus,

H_{d}^{V} (r) \to 0

at limit of

L \to 0

⁠.

The differences between the dipole method and other methods are bounded by higher-order terms in the asymptotic expansions, possibly multiplied by a factor of order one. For these higher-order corrections, we replace $ρ^{- 3}$ with $ρ^{- (1 + n)}$ in Eq. (44), where $n \geq 4$ ⁠. After summing over all cells in $V$ ⁠, the higher-order correction is on the order of $K h^{α} \times max {1, N^{α - (n - 2)}}$ ⁠, where $N \sim L / h$ ⁠. Thus, the higher-order correction to the near field is smaller than or on the same order as that for the dipole method.

In summary, provided

M (r^{'})

holds Hölder condition at

r

⁠, all these methods result in a zero near demagnetization field as

L \to 0

⁠, specifically

lim_{L \to 0} H_{d}^{V} (r) = lim_{L \to 0} H_{c}^{V} (r) = lim_{L \to 0} H_{c d}^{V} (r) = lim_{L \to 0} H_{\, p}^{V} (r) = 0.

(45)

Note that to prove the demagnetization field is zero for a sufficiently small cubic volume centered at

r

⁠, Eq. is sufficient. It does not necessarily require Hölder continuity to hold everywhere in

V

⁠.

H. Lower bound on the convergence speed

We have proven that the near field is zero as

L \to 0

⁠, assuming Hölder condition holds at

r

⁠. All methods intuitively yield the same result for the far field; thus, a detailed analysis is not provided. Now, we give a lower bound on the convergence speed with the specific condition for magnetization, i.e., the Hölder condition holds at every cell. All these methods, including the exact Cauchy principal value method, yield values on the order of

K L^{α}

for the near fields, excluding

r^{'} = r

from summation. Consequently, the error in the near field is of the order

K L^{α}

⁠. For the analysis of the far field, we begin with the UMCD method, in which each cell has an error on the order of

K h^{α} / ρ^{3}

⁠. The total error for the far field is bounded by the order of

\sum_{i, j, k \notin V} \frac{K h^{α}}{ρ^{3}} \sim \sum_{I \notin V} \frac{K h^{α} I^{2}}{I^{3}} \sim K h^{α} \log (L^{m a x} / L),

(46)

where

L^{max}

is the maximum material length. Taking

L \sim h

⁠, the total error of the near and far field is on the order of

K h^{α} \log (L^{m a x} / h) .

(47)

For the other methods, the error has two parts: first, the error of the UMCD method, and second, the differences between these methods. For the far field, the difference between these methods is bounded by the asymptotic expansion of

n = 6

⁠, i.e.,

| max M | / ρ^{7}

for a cell. After summing over all cells, the total error is bounded by the order of

| max M | {(h / L)}^{4} .

(48)

The total error for the near field and the far field part of Eq. is minimized to the order of

h^{4 α / (4 + α)},

(49)

for

L \sim h^{4 / (4 + α)}

⁠. For

L \sim h^{4 / (4 + α)}

⁠, the error of the far field for the UMCD method as in Eq. is on the order of

K h^{α} \log (L^{max} / h)

⁠, which is smaller than the error described in Eq. and can be neglected. Given that this analysis is conservative and the errors are over-estimated, the actual convergence speed may exceed these estimates.

I. Implementation using FFT

Noticing the translation symmetry property of the demagnetization tensor, the summation of the demagnetization tensor weighted by magnetization is a convolution operation. According to the cyclic convolution theorem, cyclic convolution can be performed using fast Fourier transform (FFT) and inverse FFT (IFFT) operations. To avoid the side effects of cyclic convolution, zero padding on the original array

M

is needed. The detailed steps are as follows: (1) Evaluate the array

N

of dimension

2 N_{x} - 1 \times 2 N_{y} - 1 \times 2 N_{z} - 1

representing the demagnetization function, assuming

r

is at the center of the array. Then, circularly shift the center to the frontmost position. (2) Extend the dimension of

M

from

N_{x} \times N_{y} \times N_{z}

to

(2 N_{x} - 1) \times (2 N_{y} - 1) \times (2 N_{z} - 1)

by padding zeros on the end side. (3) Apply the cyclic convolution to

M

and

N

according to the cyclic convolution theorem,²⁴

M \otimes G = I F F T (F F T (M) \times F F T (N)),

(50)

where

\otimes

is the cyclic convolution operation,

\times

is the element-wise product operation. (4) Finally, clip the

N_{x} \times N_{y} \times N_{z}

elements of

M \otimes G

at the front. Except for rounding errors, the result would be identical to the direct summation. However, the time complexity is reduced from

O (N_{x}^{2} N_{y}^{2} N_{z}^{2})

to

O (N_{x} N_{y} N_{z} \log (N_{x} N_{y} N_{z}))

using FFT. Note the

FFT (N)

is independent of

M

⁠, thus it can be reused for different

M

s provided the dimension of the problem is unchanged.

III. NUMERICAL VALIDATION

To validate the conclusion that the three methods agree under certain conditions, we construct two problems. Problem I with magnetization

M = e_{z} M_{0} {\begin{matrix} | η_{x} |^{α} (1 + η_{y}^{2} + η_{z}^{2}) & \, f o r 0 \leq η_{x} \leq 1, | η_{y} | \leq 1, and | η_{z} | \leq 1, \\ 0 & o t h e r w i s e, \end{matrix}

(51)

and problem II with magnetization

M = e_{z} M_{0} {\begin{matrix} | η_{x} |^{α} + | η_{y} |^{α} + | η_{z} |^{α} & \, f o r 0 \leq η_{x} \leq 1, | η_{y} | \leq 1, a n d | η_{z} | \leq 1, \\ 0 & o t h e r w i s e, \end{matrix}

(52)

where

α > 0

⁠,

η_{x} = x / L_{0}

⁠,

η_{y} = y / L_{0}

⁠,

η_{z} = z / L_{0}

⁠,

L_{0}

is the length of the magnet,

M_{0}

is a constant with a unit of magnetization. The magnetizations are shown in Fig. 2. For problem I, the magnetization is discontinuous at the boundaries of the magnet at

x = L_{0}

⁠,

z = \pm L_{0}

⁠, and

y = \pm L_{0}

⁠. While magnetization is continuous at the boundary of the magnet at

x = 0

⁠, its derivative is not for

0 < α \leq 1

⁠. For problem I, at the center of each cell,

r

⁠, we have

| M (r^{'}) - M (r) | \leq 3 M_{0} | r^{'} - r |^{α} / L_{0}^{α}

⁠, thus it satisfies the Hölder continuous condition. For problem II, magnetization is discontinuous at the boundaries except for the point

x = y = z = 0

⁠. At

x = y = z = 0

⁠, magnetization is continuous; however, its derivative is not for

0 < α \leq 1

⁠. The Hölder condition still holds at

r = 0

⁠.

FIG. 2.

View large Download slide

Upper panel: magnetization for $z = 0$ described by Eq. (51) for different $α$ ’s. Lower panel: magnetization for $z = 0$ described by Eq. (52) for different $α$ ’s.

For problem I, the exact analytical result can be expressed using special functions, as shown in Appendix B, with the aid of symbolic processing software.²⁵ Thus, the result can be effectively calculated to arbitrary precision. We calculate the result to at least 20 digits, although the precision is truncated to about 16 digits when converting from a multi-precision representation to a double-precision floating-point format. For problem II, the exact result is zero for symmetry reasons.

On the numerical side, we employ several methods: the UMC method, the dipole method, and a variant of the dipole method described in Eq. (39) with $b = 0.2 ı$ (where $ı$ represents the imaginary unit), the asymptotic expansions for the UMC and UMCD methods. In the asymptotic expansion methods, the asymptotic expansions are used for both far and near cells. The coordinate system is set such that the origin is at the center of a cell. The cell size is adjusted so that at $x = L_{0}$ ⁠, $y = \pm L_{0}$ ⁠, and $z = \pm L_{0}$ ⁠, the boundaries of the magnet align with the boundaries of the cells. At $x = 0$ ⁠, the cells’ centers are on the magnet’s boundary. Thus, we minimize the error fluctuations arising from the discretization of the magnet boundary. Computations are performed using double-precision floating-point arithmetic. The demagnetization tensors for the UMC and UMCD methods are calculated using exact formulas for cells at short distances and asymptotic expansions for cells at long distances. These asymptotic expansions include terms up to $O ({(h / r)}^{7})$ ⁠. The demagnetization tensor has a maximum error on the order of $10^{- 12}$ around the crossover between analytical formulas and asymptotic expansions. The numerical errors observed in this experiment, which are larger than $10^{- 7}$ ⁠, significantly exceed those of double-precision floating-point arithmetic (⁠ $10^{- 16}$ ⁠) and the error of the demagnetization tensors. Therefore, all errors are attributed to discretization errors.

The numerical errors of $H_{z}$ for $α$ of 0.2, 0.6, 1, and 2 at $r = 0$ are illustrated in Figs. 3 and 4. The numerical errors of $H_{z}$ are modeled by $c \cdot {(h / L_{0})}^{κ + d {(\log (h / L_{0}))}^{- 2}}$ ⁠. The parameter $κ$ determines the convergence speed as $h \to 0$ ⁠. $κ$ as a function of $α$ are shown in Fig. 5. We observe that $κ > 0$ for $α > 0$ ⁠. Thus, we verify that these methods converge to the same value for $α > 0$ ⁠, which is consistent with the prediction of our analysis. All convergence speeds slow down as $α \to 0$ where magnetization becomes more non-smooth. In problem I, the convergence speeds for the UMC method UMCD methods are much faster and have almost identical convergence speeds. These two methods are significantly better than the other methods. This can be explained as follows: for the magnetization relevant to the term $η_{x}^{α} e_{z}$ ⁠, the volume magnetic charge is zero and only surface magnetic charges appear on $z = \pm L_{0}$ ⁠. For these two methods, the error can be attributed to the discretization error for the magnetic surface charge at the magnet surface of $z = \pm L_{0}$ ⁠. This surface is far from where the field is under evaluation; thus, the error is suppressed. For the magnetization relevant to $η_{x}^{α} (η_{y}^{2} + η_{z}^{2}) e_{z}$ ⁠, the discretization error from the near field is highly suppressed by the term $η_{y}^{2} + η_{z}^{2}$ ⁠. Thus, these two methods show an advantage over the other methods. We also note that the variation of the dipole method shows a better convergence speed than the dipole method. The asymptotic expansion does not have an advantage compared to other method. This is as expected since the asymptotic expansions are not good for near cells. The convergence speed is faster than our conservative estimate; thus, our results are not violated. For $α \geq 2$ ⁠, the $κ$ value is close to 2 for all methods, consistent with other studies.^26,27 In problem II, all these methods exhibit similarly slow convergence speeds. This issue can be attributed to the discontinuity at $x = 0$ ⁠, which worsens the results from all methods. Additionally, as $α$ approaches zero, $κ$ for all methods approaches zero, indicating non-convergence. This confirms the assertions of this work.

FIG. 3.

View large Download slide

Error in the demagnetization field of numerical calculations at $x = y = z = 0$ as a function of cell size using magnetization for problem I and for different values of $α$ ⁠.

FIG. 4.

View large Download slide

Error in the demagnetization field of numerical calculations at $x = y = z = 0$ as a function of cell size using magnetization for problem II and for different values of $α$ ⁠.

FIG. 5.

View large Download slide

The converge speed as a function of $α$ for problem I (upper panel) and problem II (lower panel).

IV. DISCUSSIONS

One might use a rectangular prism instead of a simple cube. For the uniformly magnetized prism method, the macroscopic demagnetization field can still be correctly calculated.²⁶ However, in the dipole method, using different cell sizes along different directions will disrupt Eqs. (34) and (35), leading to discrepancies between the two methods. They represent different physics in this case. It is worth noting again that this discrepancy cannot be eliminated as cell size approaches zero.

For sufficiently smooth magnetization, the convergence speed estimated by Eqs. (47) and (49) is not faster than $h$ ⁠. However, numerical experiments indicate a convergence speed of $h^{2}$ for sufficiently smooth magnetization. Therefore, the lower bound of the convergence speed is much worse than that observed in numerical experiments. To achieve a better lower bound estimate for the convergence speed, one may need methods beyond the simple near and far field splitting, and cell-by-cell analysis. Given the potential for error cancelation among cells, it is crucial to carefully analyze the correlations between errors of different cells.

The section of our proof that demonstrates the near field is zero for uniform magnetization closely resembles the approach described in the literature.¹⁷ However, it serves a different purpose: while the literature¹⁷ aims to establish the relationship between molecular polarizability and electric susceptibility, our goal is to prove that different magnetization tensors yield consistent results for both non-uniform and uniform magnetization.

The conclusion does not readily extend to two-dimensional materials in three-dimensional space. The magnetization of a two-dimensional thin material that is one cell thick cannot be considered a three-dimensional Hölder continuous function. Specifically, the symmetries between $x$ ⁠, $y$ ⁠, and $z$ axes are broken for two-dimensional materials.

V. CONCLUSION

Due to the discretization in numerical computation, the target object is typically divided into small cubic cells. The magnetic field at the target location is the sum of contributions from sources beyond the target cell. We proved that the UMC, UMCD, and dipole methods, their asymptotic expansions, and the Cauchy principal value method yield consistent results in three-dimensional space as the cell size approaches zero, namely,

lim_{h \to 0} H_{d} (r) = lim_{h \to 0} H_{c} (r) = lim_{h \to 0} H_{c d} (r) = H_{p} .

(53)

If the self-contribution of the cell to the demagnetization field is included in both the UMC and the UMCD methods, then

\begin{aligned} lim_{h \to 0} H_{d} (r) - \frac{1}{3} M & = lim_{h \to 0} H_{c}^{\, f u l l} (r) \\ = lim_{h \to 0} H_{c d}^{\, f u l l} (r) = H_{\, p} - \frac{1}{3} M = H . \end{aligned}

(54)

The agreement among these demagnetization tensors is significant. In the calculation of the demagnetization tensor for the UMC method, both numerical integrals and lengthy exact analytical formulas introduce complexity in implementation. Additionally, exact analytical formulas can suffer from numerical stability issues, such as catastrophic cancelation.^16,21 In contrast, the implementation of the demagnetization tensor in the dipole method is considerably simpler. For numerical calculations, the consistency among these methods allows us to employ the dipole method in situations where convergence speed is not a critical factor.

A comprehensive theoretical analysis is conducted to study the effect of discretization on field errors, offering valuable insights into the relationship between error amplitudes and cell size, addressing the scientific question. The theoretical analysis is presented in detail. The approach used is rigorous, comparing numerical solutions to analytical solutions to determine the errors. The findings in this manuscript are valuable to the micromagnetics and MRI communities, providing a better understanding of different methods for calculating demagnetization fields and how to use those methods more accurately.

We acknowledge that, while the conclusions are applicable to three-dimensional materials, further investigation is required for two-dimensional materials in three-dimensional space.

ACKNOWLEDGMENTS

This work was in part supported by National Institutes of Health (NIH) Grant No.R01 EB 031078. The content is solely the responsibility of the authors and does not necessarily represent the official views of the National Institutes of Health.

AUTHOR DECLARATIONS

Conflict of Interest

The authors have no conflicts to disclose.

Author Contributions

Hao Liang: Conceptualization (equal); Formal analysis (lead); Investigation (equal); Methodology (lead); Software (lead); Validation (equal); Writing – original draft (lead); Writing – review & editing (equal). Xinqiang Yan: Conceptualization (equal); Funding acquisition (lead); Investigation (equal); Project administration (lead); Resources (lead); Supervision (lead); Validation (equal); Writing – original draft (supporting); Writing – review & editing (equal).

DATA AVAILABILITY

The data that support the findings of this study are available from the corresponding author upon reasonable request.

APPENDIX A: ASYMPTOTIC EXPANSIONS FOR THE UMCD AND UMC METHODS

We start from the UMCD method. Using

\nabla = - \nabla^{'}

and substituting the following expansion:

\frac{1}{| r + r^{'} |} = \sum_{m = 0} a_{n} \frac{{(2 r \cdot r^{'} + | r^{'} |^{2})}^{m}}{| r |^{2 m + 1}},

(A1)

into Eq. , the demagnetization tensor element

{(N_{c d} (r))}_{i j}

can be expressed in Cartesian coordinates as

\begin{aligned} - \frac{1}{4 π} & \sum_{n = 0} a_{m} \int_{- h / 2}^{h / 2} \int_{- h / 2}^{h / 2} \int_{- h / 2}^{h / 2} d x^{'} d y^{'} d z^{'} \\ \times (\frac{\partial}{\partial x_{i}^{'}} \frac{\partial}{\partial x_{j}^{'}}) \frac{{(2 (x^{'} x + y^{'} y + z^{'} z) + x^{' 2} + y^{' 2} + z^{' 2})}^{m}}{r^{2 m + 1}}, \end{aligned}

(A2)

where

x_{1}^{'} = x^{'}

⁠,

x_{2}^{'} = y^{'}

⁠, and

x_{3}^{'} = z^{'}

⁠. By substituting expansion of

{(2 (x^{'} x + y^{'} y + z^{'} z) + x^{' 2} + y^{' 2} + z^{' 2})}^{m}

as

\sum_{p = 0}^{m} b_{m, p} {(x^{' 2} + y^{' 2} + z^{' 2})}^{p} {(x^{'} x + y^{'} y + z^{'} z)}^{m - p},

(A3)

into , we obtain a polynomial in

h

⁠,

x

⁠,

y

⁠,

z

⁠, and

1 / r

expressed as

\sum_{m} \sum_{p = 0}^{m} d_{m, p, s, t} \frac{h^{2 p + m - p + 3 - 2} x^{m - p - s - t} y^{s} z^{t}}{r^{2 m + 1}} .

(A4)

The sum of exponents on

x

⁠,

y

⁠, and

z

is

m - p

⁠. The total exponents on

h

is

2 p + m - p + 3 - 2

⁠, where

2 p

comes from

{(x^{' 2} + y^{' 2} + z^{' 2})}^{p}

⁠,

m - p

comes from

{(x^{'} x + y^{'} y + z^{'} z)}^{m - p}

⁠, 3 from

\int_{- h / 2}^{h / 2} \int_{- h / 2}^{h / 2} \int_{- h / 2}^{h / 2} d x^{'} d y^{'} d z^{'}

⁠, and

- 2

from

\frac{\partial}{\partial x_{i}^{'}} \frac{\partial}{\partial x_{j}^{'}}

⁠. By defining

n = m + p

⁠, we rewrite Eq. as

\sum_{n = 0} \sum_{p = 0}^{2 p \leq n} d_{n - p, p, s, t} \frac{h^{n + 1} x^{n - 2 p - s - t} y^{s} z^{t} {(x^{2} + y^{2} + z^{2})}^{p}}{r^{2 n + 1}} .

(A5)

For terms with odd

n

⁠, before integration, the total exponents on

x^{'}

⁠,

y^{'}

⁠, and

z^{'}

are

n

⁠, thus those terms must be zero due to parity symmetry. The term with

n = 0

is also zero because, for

n = 0

⁠,

m

must also be zero, making it obvious. Therefore, only terms with

n = 2, 4, 6, \dots

remain. After summing over

p

in Eq. , we obtain Eq. .

Now turn to the UMC method. By replacing

r^{'}

with

r^{'} - r^{″}

in Eq. , and substituting it into Eq. , the demagnetization tensor element

{(N_{c} (r))}_{i j}

can be expressed in Cartesian coordinates as

\begin{aligned} - \frac{1}{4 π h^{3}} \sum_{m = 0} a_{m} \int_{- h / 2}^{h / 2} \int_{- h / 2}^{h / 2} \int_{- h / 2}^{h / 2} d x^{″} d y^{″} d z^{″} \\ \times \int_{- h / 2}^{h / 2} \int_{- h / 2}^{h / 2} \int_{- h / 2}^{h / 2} d x^{'} d y^{'} d z^{'} (\frac{\partial}{\partial x_{i}^{'}} \frac{\partial}{\partial x_{j}^{'}}) \\ \times ((x^{'} - x^{″}) x + (y^{'} - y^{″}) y + (z^{'} - z^{″}) z \\ + {(x^{'} - x^{″})}^{2} + {(y^{'} - y^{″})}^{2} + {(z^{'} - z^{″})}^{2}))^{m} \frac{1}{r^{2 m + 1}} . \end{aligned}

(A6)

In the following, using the same strategy for the UMCD method, there is no major difficulty in obtaining Eq. for the UMC method, although the process is more complex and the specific coefficients are different. We omit the detailed process here.

APPENDIX B: SOLUTION FOR PROBLEM I

H_{z}

at

r = 0

for problem I is given by

\begin{aligned} \frac{M_{0}}{12 π {(α + 1)}^{2} (α + 2)} [- 4 π α^{2} - 18 α^{2} \log (2) - 24 α^{2} \log (2 - \sqrt{3}) + 36 α^{2} \log (\sqrt{3} - 1) - 9 α^{2} ψ^{(0)} (\frac{α + 1}{4}) \\ + 9 α^{2} ψ^{(0)} (\frac{α + 3}{4}) - 8 π α - 36 α - 54 α \log (2) - 72 α \log (2 - \sqrt{3}) + 108 α \log (\sqrt{3} - 1) \\ - 27 α ψ^{(0)} (\frac{α + 1}{4}) + 27 α ψ^{(0)} (\frac{α + 3}{4}) - 18 ψ^{(0)} (\frac{α + 1}{4}) + 18 ψ^{(0)} (\frac{α + 3}{4}) \\ - 24 \sqrt{2} α^{2} F_{1} (\frac{α + 1}{2}; - \frac{1}{2}, 1; \frac{α + 3}{2}; - \frac{1}{2}, - 1) - 36 \sqrt{2} α F_{1} (\frac{α + 1}{2}; - \frac{1}{2}, 1; \frac{α + 3}{2}; - \frac{1}{2}, - 1) \\ - 24 \sqrt{2} F_{1} (\frac{α + 1}{2}; - \frac{1}{2}, 1; \frac{α + 3}{2}; - \frac{1}{2}, - 1) + 18 \sqrt{2} α^{2}_{2} F_{1} (\frac{1}{2}, \frac{α + 1}{2}; \frac{α + 3}{2}; - \frac{1}{2}) \\ + 18 \sqrt{2} α_{2} F_{1} (\frac{1}{2}, \frac{α + 1}{2}; \frac{α + 3}{2}; - \frac{1}{2}) + 12 \sqrt{2}_{2} F_{1} (\frac{1}{2}, \frac{α + 1}{2}; \frac{α + 3}{2}; - \frac{1}{2}) \\ + 18 α^{2} Φ (- 1, 1, \frac{α + 3}{2}) + 54 α Φ (- 1, 1, \frac{α + 3}{2}) + 36 Φ (- 1, 1, \frac{α + 3}{2}) - 4 π - 72 - 36 \log (2) \\ - 48 \log (2 - \sqrt{3}) + 72 \log (\sqrt{3} - 1)], \end{aligned}

(B1)

where log is the natural logarithm,

ψ^{(0)}

is the polygamma function of order 0,

F_{1}

is the Appell series,

_{2} F_{1}

the hypergeometric function, and

Φ

is the Lerch transcendent.

REFERENCES

1.

Y.

Nakatani

,

Y.

Uesaka

, and

N.

Hayashi

, “

Direct solution of the Landau–Lifshitz–Gilbert equation for micromagnetics

,”

Jpn. J. Appl. Phys.

28

,

2485

(

1989

).

https://doi.org/10.1143/JJAP.28.2485

Google Scholar

Crossref

2.

C.

Abert

, “

Micromagnetics and spintronics: Models and numerical methods

,”

Eur. Phys. J. B

92

,

1

–

45

(

2019

).

https://doi.org/10.1140/epjb/e2019-90599-6

Google Scholar

Crossref

3.

A. J.

Newell

,

W.

Williams

, and

D. J.

Dunlop

, “

A generalization of the demagnetizing tensor for nonuniform magnetization

,”

J. Geophys. Res.: Solid Earth

98

,

9551

–

9555

, https://doi.org/10.1029/93JB00694 (

1993

).

https://doi.org/10.1029/93JB00694

Google Scholar

Crossref

4.

B.

Van de Wiele

,

F.

Olyslager

,

L.

Dupré

, and

D.

De Zutter

, “

On the accuracy of FFT based magnetostatic field evaluation schemes in micromagnetic hysteresis modeling

,”

J. Magn. Magn. Mater.

322

,

469

–

476

(

2010

).

https://doi.org/10.1016/j.jmmm.2009.09.077

Google Scholar

Crossref

5.

R.

Victora

and

P.-W.

Huang

, “

Simulation of heat-assisted magnetic recording using renormalized media cells

,”

IEEE Trans. Magn.

49

,

751

–

757

(

2013

).

https://doi.org/10.1109/TMAG.2012.2219300

Google Scholar

Crossref

6.

R.

Ferrero

and

A.

Manzin

, “

Adaptive geometric integration applied to a 3D micromagnetic solver

,”

J. Magn. Magn. Mater.

518

,

167409

(

2021

).

https://doi.org/10.1016/j.jmmm.2020.167409

Google Scholar

Crossref

7.

N.

Inami

,

Y.

Takeichi

,

C.

Mitsumata

,

K.

Iwano

,

T.

Ishikawa

,

S.-J.

Lee

,

H.

Yanagihara

,

E.

Kita

, and

K.

Ono

, “

Three-dimensional large-scale micromagnetics simulation using fast Fourier transformation

,”

IEEE Trans. Magn.

50

,

1

–

4

(

2013

).

https://doi.org/10.1109/TMAG.2013.2278221

Google Scholar

Crossref

8.

A. L.

Wysocki

and

V. P.

Antropov

, “

Micromagnetic simulations with periodic boundary conditions: Hard-soft nanocomposites

,”

J. Magn. Magn. Mater.

428

,

274

–

286

(

2017

).

https://doi.org/10.1016/j.jmmm.2016.11.128

Google Scholar

Crossref

9.

B.

Müller-Bierl

,

H.

Graf

,

U.

Lauer

,

G.

Steidle

, and

F.

Schick

, “

Numerical modeling of needle tip artifacts in MR gradient echo imaging

,”

Med. Phys.

31

,

579

–

587

(

2004

).

https://doi.org/10.1118/1.1640971

Google Scholar

Crossref

PubMed

10.

B.

Müller-Bierl

,

H.

Graf

,

G.

Steidle

, and

F.

Schick

, “

Compensation of magnetic field distortions from paramagnetic instruments by added diamagnetic material: Measurements and numerical simulations

,”

Med. Phys.

32

,

76

–

84

(

2005

).

https://doi.org/10.1118/1.1828674

Google Scholar

Crossref

PubMed

11.

Y.

Shang

,

S.

Theilenberg

,

M.

Terekhov

,

W.

Mattar

,

B.

Peng

,

S. R.

Jambawalikar

,

L. M.

Schreiber

, and

C.

Juchem

, “

High-resolution simulation of B

_{0}

field conditions in the human heart from segmented computed tomography images

,”

NMR Biomed.

35

,

e4739

(

2022

).

https://doi.org/10.1002/nbm.4739

Google Scholar

Crossref

PubMed

12.

M.

Jenkinson

,

J. L.

Wilson

, and

P.

Jezzard

, “

Perturbation method for magnetic field calculations of nonconductive objects

,”

Magn. Reson. Med.

52

,

471

–

477

(

2004

).

https://doi.org/10.1002/mrm.20194

Google Scholar

Crossref

PubMed

13.

E.

Della Torre

, “

Magnetization calculation of fine particles

,”

IEEE Trans. Magn.

22

,

484

–

489

(

1986

).

https://doi.org/10.1109/TMAG.1986.1064541

Google Scholar

Crossref

14.

M.

Schabes

and

A.

Aharoni

, “

Magnetostatic interaction fields for a three-dimensional array of ferromagnetic cubes

,”

IEEE Trans. Magn.

23

,

3882

–

3888

(

1987

).

https://doi.org/10.1109/TMAG.1987.1065775

Google Scholar

Crossref

15.

H.

Fukushima

,

Y.

Nakatani

, and

N.

Hayashi

, “

Volume average demagnetizing tensor of rectangular prisms

,”

IEEE Trans. Magn.

34

,

193

–

198

(

1998

).

https://doi.org/10.1109/20.650225

Google Scholar

Crossref

16.

M. J.

Donahue

, “Accurate computation of the demagnetization tensor” (2007), see https://math.nist.gov/~MDonahue/talks/hmm2007-MBO-03-accurate_demag.pdf

17.

J. D.

Jackson

,

Classical Electrodynamics

(

John Wiley & Sons

,

1999

), pp.

159

–

162

.

Google Scholar

18.

R.

Salomir

,

B. D.

De Senneville

, and

C. T.

Moonen

, “

A fast calculation method for magnetic field inhomogeneity due to an arbitrary distribution of bulk susceptibility

,”

Concepts Magn. Reson. Part B

19

,

26

–

34

(

2003

).

https://doi.org/10.1002/cmr.b.10083

Google Scholar

Crossref

19.

Y.

Wang

,

D.

Zhou

, and

P.

Spincemaille

, “What is the Lorentz sphere correction for the MRI measured field generated by tissue magnetic susceptibility: The spatial exclusivity of source and observer and the Cauchy principal value,” in Proceedings of the International Society for Magnetic Resonance in Medicine (ISMRM, 2015), Vol. 23, p. 1698.

20.

Y.

Wang

and

T.

Liu

, “

Quantitative susceptibility mapping (QSM): Decoding MRI data for a tissue magnetic biomarker

,”

Magn. Reson. Med.

73

,

82

–

101

(

2015

).

https://doi.org/10.1002/mrm.25358

Google Scholar

Crossref

PubMed

21.

D.

Chernyshenko

and

H.

Fangohr

, “

Computing the demagnetizing tensor for finite difference micromagnetic simulations via numerical integration

,”

J. Magn. Magn. Mater.

381

,

440

–

445

(

2015

).

https://doi.org/10.1016/j.jmmm.2015.01.013

Google Scholar

Crossref

22.

J. D.

Jackson

,

Classical Electrodynamics

(

John Wiley & Sons

,

1999

), pp.

198

–

199

.

Google Scholar

23.

B. F.

Edwards

,

D. M.

Riffe

,

J.-Y.

Ji

, and

W. A.

Booth

, “

Interactions between uniformly magnetized spheres

,”

Am. J. Phys.

85

,

130

–

134

(

2017

).

https://doi.org/10.1119/1.4973409

Google Scholar

Crossref

24.

W. H.

Press

,

W. T.

Vetterling

,

S. A.

Teukolsky

, and

B. P.

Flannery

,

Numerical Recipes

(

Citeseer

,

1988

).

Google Scholar

25.

Wolfram Research Inc.

, “Mathematica, Version 14.0,” Champaign, IL (2024).

26.

J. E.

Miltat

and

M. J.

Donahue

, “Numerical micromagnetics: Finite difference methods,” in Handbook of Magnetism and Advanced Magnetic Materials, Vol. 2 (NIST, 2007), pp. 742–764.

27.

C.

Abert

,

L.

Exl

,

G.

Selke

,

A.

Drews

, and

T.

Schrefl

, “

Numerical methods for the stray-field calculation: A comparison of recently developed algorithms

,”

J. Magn. Magn. Mater.

326

,

176

–

185

(

2013

).

https://doi.org/10.1016/j.jmmm.2012.08.041

Google Scholar

Crossref

On the equivalence of demagnetization tensors as discrete cell size approaches zero in three-dimensional space

I. INTRODUCTION

II. THEORETICAL ANALYSIS

A. Cauchy principal value

B. UMC and UMCD methods

C. Dipole method

D. Asymptotic expansions

E. Near and far field decomposition

F. Near field for uniform magnetization

G. Near field for nonuniform magnetization

H. Lower bound on the convergence speed

I. Implementation using FFT

III. NUMERICAL VALIDATION

IV. DISCUSSIONS

V. CONCLUSION

ACKNOWLEDGMENTS

AUTHOR DECLARATIONS

Conflict of Interest

Author Contributions

DATA AVAILABILITY

APPENDIX A: ASYMPTOTIC EXPANSIONS FOR THE UMCD AND UMC METHODS

APPENDIX B: SOLUTION FOR PROBLEM I

REFERENCES

Citing articles via

Submit your article

Sign up for alerts

Resources

Explore

pubs.aip.org

Connect with AIP Publishing

On the equivalence of demagnetization tensors as discrete cell size approaches zero in three-dimensional space

I. INTRODUCTION

II. THEORETICAL ANALYSIS

A. Cauchy principal value

B. UMC and UMCD methods

C. Dipole method

D. Asymptotic expansions

E. Near and far field decomposition

F. Near field for uniform magnetization

G. Near field for nonuniform magnetization

H. Lower bound on the convergence speed

I. Implementation using FFT

III. NUMERICAL VALIDATION

IV. DISCUSSIONS

V. CONCLUSION

ACKNOWLEDGMENTS

AUTHOR DECLARATIONS

Conflict of Interest

Author Contributions

DATA AVAILABILITY

APPENDIX A: ASYMPTOTIC EXPANSIONS FOR THE UMCD AND UMC METHODS

APPENDIX B: SOLUTION FOR PROBLEM I

REFERENCES

Citing articles via

Submit your article

Sign up for alerts

Related Content

Resources

Explore

pubs.aip.org

Connect with AIP Publishing

This Feature Is Available To Subscribers Only