A numerical simulation of a classical dipolar spin system
and a demonstration with magnetic compasses
Takeshi Nishimatsu1,2,3 , Umesh V. Waghmare4 , Yoshiyuki Kawazoe1 ,Benjamin Burton5 , Kazutaka Nagao6 , and Yoshihiko Saito7
We report a simple, safe and stimulating pedagogical demonstration with magnetic compasses
that facilitates intuitive understanding of dipolar interactions and their consequences in magnetic
materials. We determine and discuss special states of classical dipolar spin systems corresponding
to minima, maxima, and saddle-points of their energy surface. While our results for finite systems
readily show how one can expect interesting effects in nano-scale ferroelectric and ferromagnetic
materials, our analysis of macroscopic systems reveals low energy states of magnetic ordering andthe transition states on the paths of switching between them. A computer program developed fornumerical analysis of these model systems, based on the concepts of force constant matrix andan iterative energy-minimization scheme, can be readily used along with a lab-demonstration withmagnetic compasses.
I.

INTRODUCTION

Long-range dipole-dipole interaction plays importantroles in ferromagnets and ferroelectrics, especially in theformation of domain structures1 . For example, in finite ferroelectric systems, surface charges associated withuniform spontaneous polarization P cause depolarization
field which has an opposite direction to P , due to the inherently long-range interaction. In general, for reducing
the overall energy in the presence of depolarization field,
finite systems exhibit domain structures and stabilize.
But the dipole-dipole interaction, due to its long-range
and strongly anisotropic nature, is not an easy problem
for beginners in physics and material science.
It is also known that the stable states of dipolar
spin systems are sensitive to periodicity, dimensionality, boundary conditions (finite or infinite), size, aspectratio2 , symmetry, and defects. Recent developments
in nanotechnology introduce a new paradigm for ferromagnetic and ferroelectric finite-size particles and wires:e.g. an x-ray photo-lithography3 , a mask technique withmolecular bean epitaxy (MBE)4 , and a solution based
method5,6 . The role of dipole-dipole interaction has become more important than ever in such dipolar nanosystems.

1

Institute for Materials Research (IMR), Tohoku University,
Sendai 980-8577, Japan
2 Department of Physics and Astronomy, Rutgers University, 136
Frelinghuysen Road, Piscataway, NJ 08544-8019
3 Japan Society for the Promotion of Science (JSPS) Postdoctoral
Fellow for Research Abroad
4 Theoretical Sciences Unit, Jawaharlal Nehru Centre for Advanced Scientific Research (JNCASR), Jakkur, Bangalore, 560
064, India
5 National Institute of Standards and Technology (NIST), 100 Bureau Drive STOP 8520, Gaithersburg, MD 20899-8520
6 Research Institute of Electrical Communication (RIEC), Tohoku
University, Sendai 980-8577, Japan
7 Osaka Science Museum, 4-2-1 Nakanoshima, Kita-ku, Osaka 5300005, Japan

Microscopic dipolar spin systems have been studiedin crystals and crystalline surfaces7 , and in interfaces.
Recently, however, Wang et al. used a lithographic
method8 to synthesize a two-dimensional nanoscale frustrated dipole system of 80,000 elongated permalloy(Ni0.81 Fe0.19 ) islands on a native-oxide covered siliconsurface.
Here, we provide a macroscopic model for classicaldipolar spin system: an array of magnetic-compasses, inwhich compasses (classical spins) interact via magnetic
dipole-dipole interaction. This system is safe, simple, inexpensive, rearrangeable, and can be easily manipulatedwith hands and/or magnets.
In the following section (Sec. II), we introduce setup oflaboratory demonstration with magnetic compasses anddescribe some interesting results obtained with it. In
Section III, We discuss briefly the methods used in computational analysis and its results for ground and otherspecial states of various cases of dipolar spin systems. We
also point out the differences between infinite and finite
dipolar systems. Finally, we summarize in Section IV.

II.

MAGNETIC COMPASS DEMONSTRATION

If one puts two magnetic compasses on a table, well
separated from each other, both point north as pictured
in Fig. 1(a); i.e. earth’s magnetic field (24–66 µT)
dominates the dipole-dipole interaction. As one brings
them closer together, however the dipole-dipole interaction takes over below a certain distance. Compass
needles, then, align along the line between them, i.e.
< S N>< S N> in Fig. 1(b). As in Fig. 1(c), if one puts
three, four, or more compasses in a row, one also findsthat they align according to < S N>< S N >< S N > · · · .
Now, if we arrange the compasses in a square array on aboard and give a jiggling to the board, the compass needles relax into the states pictured in Fig. 2. As discussed
in section III, such systems have minimum, maximum,and saddle-point states. Figs. 2(a) and (b) are minimum
states.

2

FIG. 1: (a) Two car-windshield magnetic compasses on atable, sufficiently far apart they both point north. (b) sufficiently close together, the dipole-dipole interaction dominates. (c) Five compasses in a row.

FIG. 3: (a) A triangular array of about 1000 magnetic compasses, on a rotating table top, that is on display in the Osaka Science Museum. (b) Magnetic domains are visible. Observers can make sure that the effects of earth’s magnetic
field is negligible by rotating the table on which compasses
are placed. Domain structures do not change when the table
is rotated.

FIG. 2: Observed stable states in the 4 × 4 square magneticcompass array. In section III, it will be explained that (a)and (b) are minimum states and (b) is the ground state.

One of the authors, YS, invented this demonstration9
and displayed a triangular array with about 1000 magnetic compasses in the Osaka Science Museum as depicted in Fig 3(a). The domain structure clearly visible (Fig 3(b)) in this finite and symmetric triangulararray has been attracting both children and adults, andstimulating more questions pertaining to the nature ofmagnetism in crystals.10
These experiments are easily performed when one uses
magnetic compasses that are designed for use on automobile windshields, as shown in Figs. 1, 2, and 3. Eachcompass has a transparent 28 mm diameter bubble-likeplastic container filled with oil, a suction cup for attach-ing it to a windshield, and a ferrite permanent magnet
arrow that floats freely in the oil. Dimensions of the ferrite magnet are 8 × 6 × 3 mm. The 8 × 3 mm arrowsides are magnetic poles with surface magnetic fields of≈ 50 mT, that decay to 1.5 mT at the surface of theplastic container. For compasses of this strength, the influence of earth’s magnetic field on them can be neglectedin comparison with that from the interaction between two
compasses, if they are placed within 50–60 mm of eachother. For compasses arranged on a square array suchthat the nearest neighbor compasses touch each other,
dipole-dipole interactions dominate over the influence ofearth’s magnetic field for the first three nearest neighborcompasses.

III.

COMPUTATIONAL ANALYSIS
A.

Force constant matrix

Consider a 2-dimensional finite system of N dipoles.
The dipoles are placed at coordinates rI = (xI , yI ),I = 1, · · · , N , their momenta are µI = (µIx , µIy ), andtheir interactions are exclusively dipole-dipole in charac-

3ter. The total interaction energy V is
1 XX
V ({µI }) =µIα ΦIα,Jβ µJβ ,2(a) maximum(b) saddle-point(c) saddle-point+2.3863-1.8694-1.9580(d) minimum(e) minimum(f) minimum-1.9758-1.9840-2.1670

(1)

Iα Jβ

where α and β are independent orthogonal directions (xand y), and ΦIα,Jβ is the force constant matrix∂∂1δαβ − 3(drIJ )α (drIJ )β=,∂rIα ∂rJβ |rI − rJ ||rIJ |3(2)
where rIJ = rI − rJ , a hat indicates a vector of unit
magnitude, e.g. rb = r/|r|, and ΦIα,Jβ = 0 for I = J.
We fix the amplitude of each moment µI to unity,X(3)µ2Iα = µ2Ix + µ2Iy = 1 ,ΦIα,Jβ = −αto simulate the magnetic-compass system. For finding
extrema (minima, maxima, and saddle-points) of the total interaction energy (1), subject to the constraint (3),
we use an iterative minimization technique: First, only
once at the beginning, we prepare the force constant matrix as in (2). From a spin configuration µIα(k) of the k-th
iteration,
we evaluate the force exerted on each dipole,P− Jβ ΦIα,Jβ µJβ(k) , then, we update from µIα(k) toµIα(k+1) :XµIα(k+1) 0 = µIα(k) − γ
ΦIα,Jβ µJβ(k) ,(4)Jβwith normalization to satisfy the fixed-amplitude restriction (3),where for I 6= J,
cos(θJ − θI ) + 3 sin(θJ − φ) sin(φ − θI )∂2V=,
∂θI ∂θJ|rIJ |3(8)
and for I = J,
X − cos(θJ − θI ) + 3 cos(θJ − φ) cos(φ − θI )
∂2V
=
.
2
∂θI
|rIJ |3
J6=IµI(k+1) 0µI(k+1) =.
|µI(k+1) 0 |(5)

In Eq. (4), γ is the adjustable dumping factor for steadyminimization. Larger positive γ may achieve faster convergence to the minima but may possibly fall into an
oscillation in the updating scheme. On the other hand, if
γ is set to a negative value, the scheme will converge toa maximum energy state. The initial guess of spin configuration µIα(1) determines which specific special (minimum, maximum or saddle-point) state the scheme will
converge to.
The total interaction energy V is also expressed in rotational angle θI of each dipole as
V ({θI }) =FIG. 4: Maximum state (a), saddle-point states (b) and (c),and minimum states (d)–(f) for the 4 × 4 square array. (f)is the ground state. Below each panel, averaged interactionenergy per dipole V /16 is indicated. In the idealized system,
the dipoles are “points” at each center of circle. Here, we
emphasize the directions of dipoles and their packings witharrows and circles.

(9)
Using this θI expression of the total interaction energy,
we can classify extrema into minimum, maximum, andsaddle-point states: At the “minimum ¯state”, all the2
V ¯
eigenvalues of the N × N matrix ∂θ∂I ∂θ
are
¯
J
0
0
θI =θI ,θJ =θI

positive. At the “maximum state” all the eigenvalues arenegative. At the “saddle-point state” some of eigenvalues
are positive and others are negative. We also define the“ground state” as the lowest-energy minimum state. In
following two subsections, using the minimization scheme(4) and (5), we search ground states of some square andtriangular dipole arrays which have high symmetries.

1 X cos(θJ − θI ) − 3 cos(θJ − φIJ ) cos(φIJ − θI ),2
|rIJ |3

B.

Square array of dipoles

I6=J

(6)
where φIJ is the rotational angle of rIJ , i.e. tan φIJ =yIJ /xIJ . Because the first θI derivertives of Eq. (6) are
zero for all I in the neighborhood of extremum V 0 =V ({θI0 }), the total interaction energy can be expressedwith deviation angles δθI of each dipole as¯
X X ∂2V ¯0¯δθI δθJ ,(7)
V =V +∂θI ∂θJ ¯θI =θ0 ,θJ =θ0
I

J

I

J

In Fig. 4, we show the calculated maximum, saddlepoint, and minimum states for the 4 × 4 square array.
It is very interesting that the state shown by Fig. 4(c)is a saddle-point state. On the demonstration with themagnetic compasses, the state as Fig. 4(c) behaves likea stable state. It may be thought that this phenomenonshows a reasonable state as follows. When the set of
4 × 4 dipoles is divided into by half of the top and bottom, the structure of Fig. 4(c) is made of two clusters of

4

k

-1.6
-1.8

y

V / N [µ2 / a3]

k
y
x

x

(a) k=(π/a, 0)(b) k=(0, π/a)

4 16 36 64

N = M2
144

256

400

finiteinfinite (≈ -2.5495)(a) square

-2.0

-1.6
-1.8
-2.0

-2.2

-2.2

-2.4

-2.4

-2.6

-2.6

-2.8

N = 3M(M-1)
6 18 36 60 90 168 270

396

finiteinfinite (≈ -2.7585)(b) triangular-2.8
2 4 6 8 10 12 14 16 18 20
M

2 3 4 5 6 7 8 9 10 11 12
M

FIG. 6: Ground-state averaged interaction energy per dipole
V /N for: (a) M × M = N finite square arrays; (b) finite
triangular arrays with M dipoles on a side and a total of
N = 3M (M −1) dipoles. The interaction energy per dipole in
the ground-state for the infinite square and triangular latticesystem is plotted as dashed holizontal asymptote lines.

yx(c) linear combination

FIG. 5: (a) Antiferromagnetic ground-state of the infinite
periodic dipole square lattice which has wave vector k =(π/a, 0), i.e. —(rI ) = (0, 1) exp[i πa xI ], where a is the distance
between nearest neighbor dipoles. —(rI ) are illustrated withsmall arrows at each site. The direction of the wave vector
k is indicated with a dashed arrow. (b) Degenerate antiferromagnetic ground-state which has wave vector k = (0, π/a),i.e. —(rI ) = (−1, 0) exp[i πa yI ]. (c) Canted antiferromagnetic
structure whichof (a) and (b), i.e.
√ is a linear combination√—(rI ) = (0, 1/ 2) exp[i πa xI ] + (−1/ 2, 0) exp[i πa yI ].

2 × 4 dipoles, which show the same type on the groundstate each other. If these two clusters got closer to fromfar away each other and the state shown by Fig. 4(c), it
would be easily imagined that the structure of Fig. 4(c)is stable with the line symmetry. It, however, is impossible that the magnets of each compass put precisely likeFig. 4(c), because an error of around 10 percentage of thesize occurs in structure of each compass. So, the two clusters on the demonstration with the magnetic compasses
are stable. If the magnets of each compass put precisely
as Fig. 4(c), the line symmetry would brake and theywould fall into the structure of Fig. 4(d). It is interestingthat the above-mentioned consideration on the saddlepoint state shown by Fig. 4(c) and the minimum point
shown by Fig. 4(d) resembles the spontaneously broken
symmetry on Higgs field in elementary particle physics
and creation of the universe in cosmology.
It is surprising that even this high-symmetry 4 × 4
square array exhibits three kinds of minimum states.
It is known that the ground-state of the infinite periodic square lattice of dipoles forms a continuously degenerate manifold of antiferromagnetic states11 , which
is a linear combination of two independent spin configurations as illustrated in Fig. 5 (In this report, we
use the word “array” for finite systems and “lattice” forinfinite periodic systems.). The interaction energy per
dipole for the infinite periodic square lattice is −2.5495(We express the interaction energy in units of µ2 /a3 thedipole moment squared divided by the nearest neighbordistance cubed.).7 The averaged interaction energy per

FIG. 7: Calculated ground-state of the N = 20 × 20 = 400
dipole square array.

dipole V /N for the M × M = N finite dipole square array approaches to −2.5495 as M increases, as illustrated
in Fig. 6(a). Figure 7 shows the 20 × 20 square array,in which the spin configuration in central region array
approaches the canted antiferromagnetic structure of theinfinite system (Fig. 5(c)). Exactly the same structure
as in Fig. 7 is achieved in the 20 × 20 magnetic compass
square array, but some by-hand one-by-one alignments
and a jiggling of the board on which compasses are placed
are required to reach the ground-state configuration. It
seems that friction between the plastic container and thearrow of the compasses preventing the system from falling
into its ground state. Compasses with greater accuracy
and stronger magnets may reach the ground state moreeasily.

C.

Triangular array of dipoles

We also calculated ground-states for the two different
types of finite triangular arrays shown in Fig. 8. Bothexhibit six-fold symmetry on their outer edges, but in
Fig. 8(a) 6-fold symmetry is broken by the central dipole,
while in Fig. 8(b) full six-fold symmetry is maintained.
The latter has a defect (a missing dipole) at the centre.

5

(a)

(b)

−2.3633−2.4081

FIG. 8: Ground states of finite triangular arrays with (a) N =19 and (b) N = 18 dipoles. The averaged interaction energy
per dipole V /N are −2.3633 and −2.4081, respectively.

FIG. 9: Calculated ground-state of the N = 396 triangulararray. Six ferromagnetic domains separated by domain walls
are evident.

This demonstrates how sensitive ground-states of finite
dipolar systems are to their symmetries, defects and size.
This has direct implications to nano-scale electric dipolarsystems, for example ferroelectrics.
The ground-state of an infinite periodic triangular lattice of dipoles is also continuously degenerated but ferromagnetic.7 In finite triangular arrays, the continuous
degeneracy breaks down, and domains form to reducethe depolarization field as depicted in Fig. 3(b), Fig. 8 ormore clearly in Fig. 9. The interaction energy per dipole
in the infinite triangular lattice is −2.7585.7 We evaluate
the averaged interaction energy per dipole V /N in the
ground-states of finite triangular arrays, e.g. Fig. 8(b)and Fig. 9, in which there are: M dipoles on a side; a totalof N = 3M (M − 1) dipoles per array; and 6-fold symmetry because they do not have central dipoles. Asymptotic
approach of the averaged interaction energy, to that of
the infinite system, is shown in Fig. 6(b). Again, a spin
configuration that is identical to that in Fig. 9 can beachieved in the magnetic-compass array with one-by-one
manipulation and a jiggling of the board on which compasses are placed; analogous to simulated annealing.

FIG. 10: On a numerical simulation, a crystal structure of the397 dipoles on a hexagon boundary condition is transformed
from the triangle lattice arrangement (a) to the rectangle oneby the quasi-static process. Though antiferromagnetism is
not observed by 92.6◦ (b), a phase transition is observed at
92.55◦ (c).

D.

Phase transition between ferromagnetism andantiferromagnetism

We show calculated phase transitions between ferromagnetism and antiferromagnetism. Figure 10 shows thephase transitions transformed from triangle array to thesquare array of the 397 dipoles on the case of a hexagon
boundary condition.
Figure 11 is the case of a diamond boundary condition.
Domains of antiferromagnetism grow up from a energetically unstable place, from the center toward the top left
side and the lower right one.

6
on the boundary condition. The angles are 97.7◦ on thecase of a diamond boundary condition and 92.55◦ on thecase of a hexagon boundary condition.
Phenomena like these phase transitions can be observed by the demonstration with tha magnetic compasses.9 10

IV.

SUMMARY

We introduced an attractive demonstration that uses
magnetic compasses to facilitate intuitive understandingof dipolar spin systems. Even with a few dozen cheap
magnetic compasses, safe and fascinating experiments
can be performed. Maximum, minimum and saddle-point
states of finite classical dipole arrays are analyzed by using the force constant matrix and an iterative energyminimization scheme. This is a good pedagogical tool
for teaching how the dipole-dipole interactions influence
the structure of special states of dipolar systems.

Acknowledgment

FIG. 11: On a numerical simulation, a crystal structure of
the 20 × 20 dipoles is transformed from the triangle latticearrangement (a) to the rectangle one by the quasi-static process. Domains of antiferromagnetism grow up toward the top
left side and the lower right one from the center by 97.7◦ (b).
A phase transition is observed at 97.6◦ in the top right sideand the lower left one (C).

The angles at which the phase transition occurs depend

1

C. Kittel, Introduction to Solid State Physics (Wiley, New
York, 2004), 8th ed, Chap. 12 and 16.
2
P. Politi, M. G. Pini, and R. L. Stamps, Phys. Rev. B 73,
020405 (2006).
3
M. Hehn, K. Ounadjela, J. P. Bucher, F. Rousseaux, D. Decanini, B. Bartenlian, and C. Chappert, Science 272, 1782
(1996).
4
C. Stamm, F. Marty, A. Vaterlaus, V. Weich, S. Egger,
U. Maier, U. Ramsperger, H. Fuhrmann, and D. Pescia,
Science 282, 449 (1998).
5
J. J. Urban, W. S. Yun, Q. Gu, and H. Park, J. Am. Chem.
Soc. 124, 1186 (2002).
6
W. S. Yun, J. J. Urban, Q. Gu, and H. Park, Nano Lett.

We thank International Frontier Center for Advanced
Materials (IFCAM) of IMR who supported UVW tovisit to Sendai and authors’ collaborative study in ferroelectrics. We would like to extend our gratitude toProfessor Gyo Takeda (Emeritus Professor of Universityof Tokyo and Tohoku University) and Professor Noboru
Takigawa (Department of Physics, Tohoku University)
who showed TN the demonstration with magnetic compasses invented by YS. We also extend our gratitude toProfessor Takayuki Hamaguchi (Nada Junior and SeniorHigh School) who gave suggestive photos and commentsof the compass arrays to the authors. We also wish to
thank an anonymous referee for a careful reading of ourmanuscript and a useful comment.
A part of this work was supported by a Grant-inAid for Encouragement of Scientists, 2006, No.18912015,
from Japan Society for the Promotion of Science.

2, 447 (2002).
V. M. Rozenbaum, V. M. Ogenko, and A. A. Chuı̆ko, Sov.
Phys. Usp. 34, 883 (1991).
8
R. F. Wang, C. Nisoli, R. S. Freitas, J. Li, W. McConville,
B. J. Cooley, M. S. Lund, N. Samarth, C. Leighton, V. H.
Crespi, et al., Nature 439, 303 (2006).
9
Y. Saito and K. Yasue, Frontier Perspectives 10, 28 (2001).
10
Y. Saito, Butsuri Kyouiku (Physics Education) 53, 103(2005), [in Japanese].
11
K. De’Bell, A. B. MacIsaac, I. N. Booth, and J. P. Whitehead, Phys. Rev. B 55, 15108 (1997).
7