Self-assembly of three-dimensional ensembles of magnetic particles with laterally shifted dipoles

We consider a model of colloidal spherical particles carrying a permanent dipole moment which is laterally shifted out of the particles' geometrical centres, i.e. the dipole vector is oriented perpendicular to the radius vector of the particles. Varying the shift $\delta$ from the centre, we analyze ground state structures for two, three and four hard spheres, using a simulated annealing procedure. We also compare to earlier ground state results. We then consider a bulk system at finite temperatures and different densities. Using Molecular Dynamics simulations, we examine the equilibrium self-assembly properties for several shifts. Our results show that the shift of the dipole moment has a crucial impact on both, the ground state configurations as well as the self-assembled structures at finite temperatures.


I. INTRODUCTION
Recent advances in particle synthesization and the permanent need for novel materials meeting more and more specialized requirements encourage the search for novel types of functionalized particles.Promising candidates in this area are colloids with directional interactions.These interactions are the key for the controlled self-assembly of colloidal particles into specific structures.Recent research on this topic resulted in complex colloidal particles characterised by complex shape [1][2][3], anisotropic internal symmetry [4,5] or surface charges [6,7].Understanding the interactioninduced behaviour of such particles is crucial for optimizing their application e.g. in material science biomedicine or sensors.Yet, not only the application of functionalized particles is of interest but also their capability to serve as model systems to study fundamental concepts of physics such as self-organization [8,9], chirality [2,[10][11][12], synchronisation [8,[13][14][15], critical phenomena [16,17], entropic effects [18,19] and active motion [20,21], to name a few.A paradigm example is the model of dipolar hard and soft spheres which is a well-established model to examine and understand the properties of magnetic colloidal particles immersed in a solvent, also called ferrofluids.From numerous studies of the phase behavior of dipolar liquids (e.g.[22,23]), and especially of their structural properties (e.g.[24,25]), it is known that the dipoles assemble into chains, rings and branched structures at sufficiently large dipolar strengths and low densities.
Here, we focus on (permanent) ferromagnetic colloids with anisotropic symmetry, i.e. the magnetic moment within the colloidal particle is not located at the geometric centre of the particle.A first theoretical description for decentrally located dipoles was introduced by Holm et al. [26][27][28].In their model, spherical particles carry a dipole moment which is shifted out of the particle * klapp@physik.tu-berlin.decentre and is oriented parallel to the raduis vector of the particle.The model describes very well the cluster formation of particles carrying a magnetic cap [29].Yet, it is insufficient to mimic the self-assembly of so-called Patchy colloids [9], that is, silicon balls carrying magnetic cubes beneath their surfaces.Furthermore, the model does not reproduce the zig-zag chained structures formed by magnetic Janus particles in an external field, i.e. particles where one hemisphere of silica spheres are covered with a magnetic coat [8,30].The concept of shifting the dipole was later extended by fixing the amount of shift and varying the orientation of the dipole moment vector within the particle [31] which was proven to be more convenient for patchy collids.In the present contribution we consider a model in which the dipole moment is laterally shifted such that the radius vector and the dipole moment vector are oriented perpendicular.The same model was also proposed in Ref. [31].However, here we fix the orientation of the dipole moment and vary the amount of shift.Thereby, we do not only aim at modeling synthesized particles mentioned above.Rather, we are interested in understanding the impact of successively shifting the dipole on the self-assembly of such particles.To this end, we perform ground state calculations for a small number of dipolar hard spheres and conduct Molecular Dynamics (MD) simulations to study the bulk at finite temperature in three dimensions.Very recently, Novak et al. have considered the same model [32], however, they restricted their study to systems where the particles are fixed in a plane with freely rotating dipoles.Thus, they considered a quasi-twodimensional (q2D) system.Besides, the authors examined the system at one fixed density.Here, we examine a three dimensional system of such particles at zero temperature and conduct MD simulations of the bulk at several thermodynamic state points.Thereby, we aim at a quantitive characterization in which the shift is the stateparameter in the system.The remainder of the paper is structured as follows.In section II, we present the model and the equations of motion, sec.III refers to the computaitonal methods and sections IV and V include the results for the ground state calculations and for the structural analysis of the bulk systems, respectively.We close the paper by a summary and outlook.

II. MODEL AND EQUATIONS OF MOTION
Our model consists of N spherical particles carrying a permanent dipole moment µ i , (i = 1, .., N ), which is laterally shifted with respect to the center of particles.A sketch is given in Fig. 1.In the body-fixed reference frame (in the following denoted by the subscript b), the location of µ i is specified by the shift vector , and its orientation is given by the vector µ b i = µ (0, 0, 1), with d and µ being constant for all particles.Hence, d i and µ i are oriented perpendicular to one another.Thus, our particles differ from those considered in Ref. [27] where d i and µ i are arranged parallel and hence µ i is shifted radially.In the laboratory reference frame, r i is the position vector of the particle centre while the position vector of µ i is given by r ′ i = r i + d i , where d i now denotes the shift vector in the laboratory frame.For d = 0, r ′ i coincides with r i yielding conventional dipolar systems with centered dipoles.The total pair potential between two particles i and j consists of a short-range repulsive potential, u short (r ij ), and the dipole-dipole potential, yielding Here, determines the distance between the dipoles.We employ two different types of repulsive interactions.First, for the finite temperature MD simulations discussed in Sec.V, the repulsive potential is modeled via the shifted soft sphere (SS) potential defined as The parameters for potential depth and steepness, ǫ and n, respectively, are specified in Sec.V.At the cut-off distance r c , the shifted potential given in Eq. ( 3) and its first derivative continuously go to zero such that corrections due to the cut-off are not required.Finally, σ represents the diameter of the particles.Second, for the ground-state calculations presented in Sec.IV, we set u short (r ij ) equal to the hard sphere (HS) potential defined as We now derive the equations of motion of the particles in the absence of a solvent.Each particle i experiences the total force F i = F i short + F i DD at its centre of mass.The force is due to steric interactions with all other particles j and is the dipolar force due to the dipole-dipole potential u DD (µ i , µ j , r ′ ij ) given in Eq. ( 1).Note that although the force F i DD acts at r ′ i , the same force also acts on the center r i due to the rigidity of the particle.Moreover, the finite shift d i generates a torque T i d = d i ×F i DD acting at r i , which supplements the torque T i µ = µ i ×G i DD stemming from angle dependent dipolar forces [33].Here, Thus, the total torque on the particle centre is given by T i = T i µ + T i d .For d = 0, i.e. r ′ i = r i , the additional torque T i d vanishes and the forces and torques reduce to the expressions familiar for centered dipoles (e.g.[34]).We also note that our treatment of the forces and torques in a system of shifted dipoles is equivalent to the virtual sites method introduced by Weber et al. [27].Having derived the forces and torques, the (Newtonian) equations of motion are given by for translation (with m beeing the mass of the particles), and FIG. 1. Sketch of a dipolar sphere with a laterally shifted dipole moment.Also shown are the axes of the body-fixed coordinate system.
for rotation [34].In Eq. ( 9), ω i is the angular velocity and I is the moment of inertia.Further, quantities in the body fixed frame can be transformed to the laboratory frame via a rotation matrix given in [34].In Eq. ( 9), the quantity Qi is the time derivative of the quaternion Q i = (q i 0 , q i 1 , q i 2 , q i 3 ) which we employ to describe the orientation of the particle (specified in [34]).The matrix W is defined as (see [34]) while the quaternion Ω i b = (0, ω ix b , ω iy b , ω iz b ) corresponds to the x, y and z components of the angular velocity.It can be shown that Eq. ( 10) is equivalent to the expression ṡi = ω i × s i known for the rotation of linear molecules (e.g.[34]), where s i is the unit vector of the particle orientation and ṡi its time derivative.

III. COMPUTER SIMULATIONS A. Molecular Dynamics simulation
In our MD simulations, we constrain the temperature to a constant value T by using a Gaussian isokinetic thermostat [34].Hence, the linear and angular momenta of the particles are rescaled by the factors X trans = T /T trans and X rot = T /T rot , respectively, where 2 (with ṙi = |dr i /dt|) and are the translational and rotational kinetic temperatures of the system.Further, k B is the Boltzmann constant.We solve the corresponding isokinetic equations for translational and rotational motion with a Leapfrog algorithm, following the schemes suggested in Refs.[34] and [35].To account for the long range dipolar interaction u DD , we apply the three-dimensional Ewald summation technique [33].Specifically, we use a cubic simulation box with side length L x = L y = L z = L and employ periodic boundary conditions in a conducting surrounding.The parameter α which determines the convergence of the real space part of the Ewald sum is chosen to be α = 6.0/L which is large enough to consider only the central box with n = 0 in the real space.For the Fourier part of the Ewald sum we consider wave vectors k up to (k) 2 = 54, giving a total number of wave vectors n k = 1500.In the MD simulations, we use the following reduced units: The simulations were carried out with N = 864 particles and with a time step of ∆t * = 0.0025.Typical simulations lasted for 3 × 10 6 steps.

B. Simulated annealing
To investigate ground state configurations of small clusters of particles interacting via the pair potential u(i, j) = u HS + u DD [see Eqs. ( 1), ( 2) and ( 4)], we employ a simulated annealing procedure which involves a Monte Carlo simulation using the Metropolis algorithm [34].Within this method, we choose initial states with comparable dipolar and thermal energies, i.e. u DD /k B T ≈ 1.Here, u DD is the dipolar energy of two hard spheres in contact with central dipoles having head to tail orientation.We then lower the temperature stepwise to zero.At each temperature, 10 6 trial moves are performed while conducting the usual Metropolis scheme involving translational and rotational trial moves.We realize an acceptance ratio of 60% by regularly adjusting the absolute value of the translational displacement during the simulation.New orientational configurations are generated by rotating the particles with a constant angle of dφ = π/18 around one of the three axes of the laboratory fixed frame.In order to ensure that we reach the state with lowest energy, we start several simulations for each set of parameters and choose those results with the lowest energy as the minimum energy state.

IV. GROUND STATE CONSIDERATIONS A. Analytical expression for the pair energy
As a first step towards understanding the impact of the lateral shift, we consider two hard spheres with shifted dipoles (see Eq. ( 2) with u short = u HS ).Specifically, we derive an analytical expression for the pair energies as function of the relative shift δ = |d|/(2|R|), where |R| = σ/2 is the particle radius.A similar derivation (leading to the same result) was very recently presented in [32].Here, we include the derivation as a background for later investigations of N > 2 particles.The basis of the derivation is the coordinate system shown in Fig. 2. Note that this is a two-dimensional system (x-z-plane) where the orientations of the dipoles along the y-axis, i.e. out-of-plane orientations, are neglected.This assumption is confirmed by simulation studies of q2D dipolar systems showing that out-of-plane fluctuations vanish for decreasing temperatures [? ].In Fig. 2, the angles α and β describe the orientations of the shift vectors d i and d j with respect to the z-axis.As the shift vector and the dipole vector have a fixed orthogonal orientation to each other, the orientations of µ i and µ j with respect to the z-axis follow as α + π/2 and β + π/2.With these definitions of the angles, the results for our lateral shift can be directly compared to those for the radial shift given in Ref. [27].Clearly, the distance |r ′ ij | varies with α, β and δ.Finally, we obtain for the dipolar potential u DD in terms of the parameters α, β and δ.This expression is equivalent to that of Ref. [32] (as can be seen after some rewriting.) We now aim at finding the minimum energy states, E G , of two dipolar hard spheres as a function of δ.To this end, we minimize Eq. ( 12) with respect to α and β and compare the results with simulated annealing calculations in three dimensions, as described in Sec.III B. As the plot in Fig. 3 clearly shows, the analytically gained results perfectly fit the numerical ones.Furthermore, it can be seen that the ground state energy E G (δ) (which agrees with that calculated in [32]) is a quantity which decreases with increasing shift.Initially, E G (δ) changes slowly and is comparable to that of nonshifted dipoles suggesting that in this region, shifting the dipole moments out of the centres does not have a significant effect on the system.Upon further increase in δ, E G (δ) starts to rapidly decrease.This is a result of the fact that shifting the dipoles out of the centres enables them to reduce their distance compared to the case with zero shift.This effect becomes more and more pronounced with growing δ as the dipolar potential of Eq. (1) follows a power law of the dipolar distance.
When the results shown in Fig. 3 are compared to the corresponding results of radially shifted dipoles of Ref. [26], a qualitative agreement of the function E G (δ) can be seen.Yet, in the case of lateral shifts, the reduction of energy sets in earlier, i.e. for smaller shifts δ than those of radial shifts for which the energy starts to decrease only at δ ≈ 0.25.Further light on this issue is shed by inspecting the ground state configurations presented in the next section.

B. Ground state pair configurations
The ground state configurations pertaining to a given shift δ are determined by those values for the angles α and β that minimize Eq. (12).In Fig. 4, the angles α and β, as well as the cosine of the enclosed angle ∆(δ) = ∡(µ 1 , µ 2 ) between the dipoles in their ground state arrangements at different shifts are shown.For δ = 0, cos(∆ = 0) = 1 holds.We also note that for δ = 0, the two sets α = β = π/2 and α = β = 3π/2 both describe the ground state configuration (see Fig. 2) of nonshifted dipoles, which is the parallel head-to-tail ori- entation.Here, we choose the latter set of initial values, α = β = 3π/2, as a starting point for our examination.Shifting the dipoles out of the centres, the parallel orientation of the dipoles is gradually abandoned in favour of reducing the dipolar distance.In detail, upon increasing δ from zero, α is reduced until it reaches the value π (see inset of Fig. 4a).Correspondingly, β grows with increasing shift towards 2π, as shown in the inset of Fig. 4. In other words, with increasing shift, the upper particle in Fig. 4 (b) rotates clockwise while the lower one rotates counterclockwise and α and β evolve in a completely symmetric manner.Thereby, ∆ increases and cos(∆) decreases, reflecting that the dipoles more and more deviate from their parallel orientation.At the value δ ≈ 0.13, cos(∆) passes the zero line where ∆ ≈ π/2 and the dipoles attain a perpendicular orientation.Finally, cos(∆) reaches its lowest value cos(∆) = −1 (and thus ∆ = π) at δ = 0.2.This corresponds to an antiparallel configuration of µ 1 and µ 2 relative to each other, and to a perpendicular orientation of each of the dipoles relative to the connecting line between the particle centres.For all higher shifts, the antiparallel orientation is kept and only the dipolar distance is further reduced.Interestingly, the value of δ = 0.2 does not point any significance in the energy plot of Fig. 3 but is highly significant for the preferred orientation of the dipoles.Thus we conclude that δ = 0.2 represents the border between the two regimes of parallel (small shifts) and anti-parallel (high shifts) orientations (consistent with [32]).Compared to radially shifted dipoles of Ref. [26], the main difference in the ground state structures is that radially shifted dipoles keep their parallel head-to-tail orientation for small shifts.For large shifts, the two radially shifted dipoles also attain an antiparallel oriented relative to each other whereas at the same time, each dipole is orientated along the connecting line between the centres of the particles.This demonstrates that not only the location but also the orientation of the dipole vector within the particle plays a crucial role for the ground states of the particles as also confirmed in Ref. [37] in which the authors study the influence of shape and geometric anisotropy of the particles on their interaction.

C. Triplet configurations and beyond
The principle impact of shifting the dipoles out of the particles' centres, namely, the decrease of the ground state energy and a preferred non-parallel orientation of the dipoles with increasing shift, becomes even more pronounced in systems of three and four dipolar hard spheres.For a detailed investigation, we have performed simulated annealing calculations to determine the ground state configurations of three-dimensional systems with three and four hard spheres for different shifts.We first consider the three-particle case.In Fig. 5(a)-(d), we sketch the obtained configurations.Starting from the chainlike headto-tail orientation known for nonshifted dipoles, the particles first organize into slightly curved chainlike geometries [Fig.5(b)].This occurs for very small shifts up to δ ≈ 0.03.Our simulation results show that the corresponding ground state energies for this curved chain configuration is indeed slightly lower (see Table I) than those for the corresponding structure proposed in Ref. [32] which the authors call a "zipper".In a "zipper" con-figuration the dipoles have head-to-tail orientation and are organized in a staggered manner.When the shift takes values above δ ≈ 0.03, the two particles at the ends of the chain approach each other in such a way that they form a planar triangular arrangement.This behavior remains for all higher shifts [Fig.5(c) and (d)], in agreement with the results of Ref. [32].In terms of the orientations of the dipoles within the particles, in the case of chainlike geometries, the dipoles show head-to-tail orientation.Within the triangular geometries, there are two qualitatively different types of dipolar orientations.The first type is likewise triangular with all the pair angles (i.e. the angles between the three dipolar pairs) attaining the value of 120 • for shifts up to δ ≈ 0.38, as confirmed by the plot in Fig. 5(e).The second type is a rectangular orientation in which two of the dipoles form an antiparallel pair and the third one joins the pair in a perpendicular manner [Fig.5(d)].Correspondingly, two of the three pair angles have a value of 90 • and the third one of 180 • , as shown in Fig. 5(e).We note that at higher shifts δ 0.4, we find again a difference to the results in [32].The authors propose a configuration containing an antiparallel pair which is joined by the third particle via a head to tail orientation with one of the dipoles of the antiparallel pair.To clarify this issue, we have derived an analytical expression for the rectangular configuration in Fig. 5

(d). It is given by
Evaluating this energy, we find that the rectangular configuration is energetically slightly more favourable than that of Ref. [32]. Figure 6 shows the results for the absolute values of u rect (δ), the results for the absolute values of Eq. ( 7) of Ref. [32], u ap+p (δ), and the difference |u rect (δ)|− |u ap+p (δ)|, which is positive for all values considered.
Finally, in the case of four particles, the nonshifted ground state configuration is a ring geometry with rectangular, cyclic orientation of the dipoles, as it is known from other ground state studies [27] [see Fig. 7(a)].This configuration remains for small shifts where only the dipolar distances are reduced while the orientations are maintained.Upon increasing the shift, opposing dipoles within the rectangular geometry more and more approach each other and form two pairs of antiparallel dipoles which are perpendicular oriented to each other.This is accompanied by a change from the planar rectangular towards a tetrahedral configuration as shown in Fig. 7(c) and (d).Thus, in the four particle system, we observe for the first time a cross-over from planar to three-dimensional configurations.

V. SYSTEMS AT FINITE TEMPERATURE A. Preliminary considerations
In this section, we investigate finite temperature systems with soft-sphere repulsive interactions, which seem more realistic for the real colloidal particles mentioned in the introduction.To this end, we set in Eq. ( 2) the parameters ǫ = 50 and n = 38.Due to the that the magnitude the ground energy E G (δ) is an increasing function of the shift (see previous discussions), also the dipolar coupling strength λ, which is defined as the ratio of the half ground state energy and the thermal energy, λ(δ) = |E G (δ)|/2k B T , becomes an increasing function of the shift.This yields an irreversible agglomeration of the particles, which cannot be counteracted by the soft-core potential.For the present choices for ǫ and n this situation occurs if the shift exceeds the value of δ = 0.33.We examined higher shifts than δ = 0.33 by appropiate choices for ǫ and n but did not gain any new insights of the system beyond those already observed for smaller shifts.Therefore, instead of adjusting λ(δ), e.g by appropriate reduction of µ * with increasing shifts, or instead of enhancing the soft-sphere potential values ǫ and n, we limit the shift at δ limit = 0.33 in order to prevent agglomeration.In this way the structural properties of the system can be directly related to the amount of shift which hence is the parameter of interest in our examinations.We consider a strongly coupled system with µ * = 3 with the densities ρ * = 0.07, ρ * = 0.1 and ρ * = 0.2 and at the two temperatures T * = 1.0 and T * = 1.35, respectively.This yields coupling strengths ranging from λ(δ = 0) = µ 2 /(k B T σ 3 ) = 9 to λ(δ = 0.33) = 72 for T * = 1.0, and λ(δ = 0) ≈ 6.67 to λ(δ = 0.33) ≈ 53.33 for T * = 1.35.For a thorough investigation of the equilibrium properties of the shifted system, we performed MD simulations and calculated various structural properties, as described in the next section.

B. Results
For a first overview, we present in Fig. 8 representative MD simulation snapshots illustrating typical self assembling structres.Specifically, we consider systems at T * = 1.0 and ρ * = 0.1 for δ = 0, δ = 0.21 and δ = 0.33 .Qualitatively, the structures appearing for the considered values of δ can be divided into four groups.These are chains (A), staggered chains (B), rings built by staggered chains (C) and small clusters (D) of the types presented in Figs.4(d), 5(c) and 7(c).Structures of type (A) can consist of a few (2 − 5) as well as of many (more than 10) particles, i.e., the chains can be short or long.Structures of types (B) and (C) always consist of more than 10 particles [Fig.8(d), (e)].In accordance with the ground state configurations (see Figs. 5 and 7), the structures found in the finite temperature systems for different shifts pass from chainlike geometries to circular close-packed clusters upon the increase of δ.Accordingly, structures of the first group are formed for zero and small shifts in the range δ = 0.01 − δ ≈ 0.1 [Fig.8(a) and (d)].In this shift region, the overall chainlike structure with head-to-tail orientation as formed by nonshifted dipoles is maintained.Yet, the shift causes more and more curved structures compared to the nonshifted particles.As is generally known for dipolar systems, the chain length, i.e. the number of particles within a chain, has a polydisperse distribution [38].This holds also for the shifted system (see also the discussion of the cluster analysis in Sec.V B 2).For intermediate shifts, e.g.δ = 0.24, Fig. 8(b) and (e), the particles within the chains become staggered and we observe coexistence of structures of the types (B), (C) and (D).Structures of group (D) are consistent with ground the state configurations of this and higher shifts.Although groups (B) and (C) are not observed for zero temperature, they can be understood as a modification of chains, as they appear for small δ, and of rings which occur at zero temperature.If δ takes values near 0.33, all large aggregates (B) and (C) vanish and only small clusters (D) remain, as shown in Figs.8(c) and (f).The same structural behaviour at the different shift regions is observerd for the other state points considered.Thus we conclude that the described self-assembly of the particles at different shifts is a quite general behaviour which results from the increasing dipolar coupling strength for increasing shifts.The latter causes more and more close-packed structures as we already confirmed in the case of hard spheres.

Radial distribution function
a first quantitative measure of the structure formation, we consider the radial distribution function for several shifts.The plots in Fig. 9 show g(r) for δ = 0 and δ = 0.33 for T * = 1.0 and T * = 1.35.The g(r) at zero shift is dominated by first and second neighbour correlations.This is a typical feature of strongly coupled dipolar systems [40? ] and reflects the formation of chain-like structures.When we successively increase the shift, the second peak exists up to a value of δ ≈ 0.25.Beyond this value, only nearest neighbour correlations at r/σ = 1 are present in the system signifying the presence of only small and close-packed clusters (D), as seen in the snap shots of Fig. 8(c).Noticeably, the results for the higher temperature T * = 1.35 completely coincide with those of T * = 1.0 in the high shift region [Fig.9(b) and (d)].This is because for sufficiently high shifts, the increase of the dipolar coupling strength is already enhanced and thus, the increase of temperature does not affect the self-assembly.

Cluster analysis
To further characterize the aggregates, we perform a cluster analysis.In particular, we are interested in the cluster size distribution for several shifts, the mean cluster size and the mean cluster magnetization as a function of δ.The basis of this analysis are distance and energy criteria.Specifically, all particles with a distance lower than r c = 1.3 σ and binding energy u c = i,i ′ >i u ii ′ DD < 0 are regarded as being clustered.Here, u ii ′ DD denotes the dipolar energy [see Eq.Specifically, values of M c (S) near to one reflect a high degree of head to tail orientation, while vanishing values of this quantity indicate antiparallel or triangular orientation.
Therefore, the mean cluser magnetization gives inside into the organization of the dipoles within the formed structures and thus allows to evaluate if a given assembly is chainlike [types (A) and (B)] or closed [types (C) and (D)].Note that the total magnetization, which is usually calculated by summing over all particles, has vanishing values as the system is globally isotropic at the state points considered here.Finally, the mean cluster size is obtained from (a) Normalised cluster size distribution.The results for h(S) for different characteristic shifts, namely for δ = 0.1 (small shift), δ = 0.16 (intermediate shift) and δ = 0.27 (high shift) are presented in Fig. 10.The figures 10(a) and (d) show that mostly large aggregates, that can contain up to 25 − 30 particles, are formed.On the other hand, Figs.10(c) and (f) indicate the formation of only small assemblies with 3 − 4 particles.However, in Fig. 10(b) and (e), although there is a preferential emergence of small assemblies, large aggregates of up to 20 particles are present in a non-negligible number and secondary peaks at e.g. S = 15 (for T * = 1.0) and S = 13 (for T * = 1.35) are visible.Evidently, for this and comparable shifts, small and large assemblies can coexist.One also finds that for higher temperature, large aggregates are less often formed than for the smaller temperature.This is indicated by the fact that the peaks in Figs.10(e (c) Mean cluster magnetization.In order to evaluate the types of the occuring structures for a given shift, we determine M as a function of the shift and plot the results in Figs.11(b) and (d).For zero and initial shifts, M takes the value ≈ 0.7, reflecting predominantly parallel orientation of the dipoles within their aggregates.From this and from the cluster size distribution [Fig.10(a),(d)] we conclude that for small shifts (up to δ ≈ 0.1), mainly short and long polar chains of type (A) or (B) are formed.If the shift is further increased, M decreases, indicating that polar chains occur less often.Instead, the aggregates become more and more closed structures of the types (C) or (D) with increasing shifts.Hence, the decrease of M implies the coexistence of types (B), (C) and (D) [see Figs.8(b) and (e)].At the high shift end, M drops down to vanishing values indicating only pairwise antiparallel or triangular arrangements of the dipoles within the clusters, which is also consisent with the results shown in Fig. 10(c) and (f).The fact that the mean cluster magnetization has vanishing values at large δ also suggests that the clusters poorly interact.
Note that for all values of δ, the according aggregates are isotropically oriented such that the total magnetization is zero for all shifts (not shown here).
(b) Mean cluster size.
Finally, we examine the influence of the shift on the mean cluster size and plot in Figs.11(a) and (c) S as a function of the shift.Starting at δ = 0, the mean cluster size grows to its maximum with about 17 particles for T * = 1.0 and about 13 particles for T * = 1.35.The maximum is reached at δ ≈ 0.05, respectively.This increase can be understood by the effective increase of the dipolar coupling strength λ (see preceding discussion) such that shifts result in the formation of longer chains of type (A).If δ exceeds this value, S starts to gradually decrease because with increasing shift, smaller aggregates are formed more frequently (see Fig. 10).Finally, S attains the value of about 3 particles in the high shift end, which is a highly representative value for both temperatures considered [Figs.10(c) and (f)].Significant differences between the results of the two temperatures can be seen only for shifts smaller than δ ≈ 0.1 where mainly chainlike aggregates are formed.Here, the increase of temperature, which involves the decrease of the coupling strength from λ = 9 to λ ≈ 6.67, causes the formation of chains with less particles.Moreover, for these values of δ, shifting the dipoles does not impose fundamentally different self-assembly patterns compared to nonshifted dipoles.Therefore, small shifts can be regarded as perturbation of the nonshifted system.On the other hand, high shifts impose significantly different structures: the particles exclusively form structures of type (D) that correspond to ground state configurations of two, three and four hard spheres [see Figs.Finally, for intermediate shifts, where large aggregates as well as small clusters are formed, the decrease of S (and at the same time of M ) can be interpreted as a transition region in which large aggregates gradually dissolve into small clusters until no large structures appear at all.Within this region, the competition between energy minimization and entropy maximization results in the coexistence of both, small and large aggregates.With increasing shift (i.e., effectively increasing λ(δ)), the particles accomplish to form structures equivalent to ground state configurations.
To summarize, in the bulk systems at the finite temperatures and densities considered here, we can qualitatively distinguish between three shift regions (small, intermediate and high) each of which is characerized by it's own structural characteristics.By contrast, in the ground states of two particles, we determined only a small and a high shift region (see Fig. 4(a) and the related discussion).The intermediate shift region, observed for the bulk systems is not detected for zero  temperature.This is consistent with the fact that the corresponding structures of types (B) and (C) are not observed in the ground state calculations.

VI. CONCLUSIONS AND OUTLOOK
In this paper, we investigated a model of spherical particles with laterally shifted dipole moments which is inspired by real micrometer sized particles that carry a magnetic component on or right beneath their surfaces (e.g.[2,8,9]).Aiming at understanding the principle impact of the shift of the dipole moment on the threedimensional system, we determined the ground state structures of two, three and four dipolar hard spheres.It turns out that shifting the dipole fundamentally affects ground state energies and configurations, as well as self-assembly patterns in finite temperature systems.For these, we could determine three regions of shift, being small, intermediate and high.In each region, the self-assembly of the particles is fundamentally different.The system passes from a state which is similar to that of nonshifted dipoles to a clustered structural state.
Further, it is an interesting observation that the asymmetry of the particles, caused by the off-centred location of the dipole moment, is overcome for small shifts insofar as the behaviour of the small shift region can be recognized as a perturbation of the nonshifted system.On the other hand, if the shift is too high, the system compensates the off-centred location of the dipole by building symmetric aggregates.
So far, we examined the equilibrium properties of systems of shifted dipoles.Further investigations should clarify the interaction between the aggregates in the different shift regions.Moreover, it would be desirable to have a full phase diagram as it is known for centred dipolar soft spheres [41].
In view of the severe effects of the shift on the equilibrium properties, one expects new types of pattern formation if the system is out of equilibrium.An interesting case are systems of shifted dipoles exposed to several types of external magnetic fields.The case of a constant field was examined in Ref. [32] demonstrating that shifted dipoles form staggered chains for appropriate values for the field strength and the shift.Even more exciting phenomena are provoked if the field is time-dependent, e.g.precessing or rotating, driving the particles accordingly into tubular [8] or crystalline structures [42] as a result of synchronization effects in the systems.The immediate interest particularly lays in the question to which extent the model of laterally shifted permanent dipoles can be used to theoretically describe phenomena observed in real systems such as Janus particles [8,42].Computer simulations in these directions are on the way.

FIG. 2 .
FIG.2.Sketch of two dipolar hard spheres i and j and the orientations of their shift and dipole vectors in the x-z-plane.

FIG. 3 .
FIG. 3. Ground state pair energy EG normalized by the corresponding ground state energy E 0 G = −2µ 2 /σ 3 of centred dipoles.The results are obtained by simulated annealing (circles) and by minimization of Eq. (12) (solid line).

FIG. 4 .
FIG. 4. (a) Simulated annealing results for the enclosed angle between the two dipoles in their ground states.The solid line indicates the zero line.The inset shows α and β as defined in Fig. 1 and gained by minimizing Eq. (12).(b)-(d) Ground state configurations of two dipoles.

FIG. 5 . 4 FIG. 6 .
FIG. 5. (Color online) (a)-(d) Ground state configurations of three particles.(e) Simulation results for the three pair angles occuring between each of the three dipolar pairs in the ground state configurations.Each symbol represents one pair angle, respectively.For rectangular dipolar configurations, e.g at δ = 0.45 [see (d)], two of the angles approach 90 • and the remaining one, accordingly, 180 • .

FIG. 8 .
FIG. 8. (Color online) Snapshots for δ = 0 (a), δ = 0.21 (b) and δ = 0.33 (c) with revealing structures of the group (A) (a), the groups (B), (C) and (D) (b) and only (D) (c).In each snapshot, some randomly chosen clusters are coloured for a better visibility.Particles of the same color besides magenta belong to the same cluster.Magenta colored clusters represent small single clusters (D).(d)-(f) Magnification of randomly chosen clusters of the snap shots in the left column.
(1)] between all pairs i, i ′ within the critical distance r c .The detected clusters were collected in a histogram in which the number of clusters with size S, n(S), is counted and normalized by the total number of clusters, N c = S≥2 n(S), such thath(S) = n(S) N c ,gives the normalised cluster size distribution.Only S ≥ 2 enters to the sum, i.e., single particles are disregarded.Based on the function n(S), the mean cluster magnetization is calculated byM = S≥2 n(S) • M c (S) N c , where M c (S) = S i=1 µ i /(µ • S)gives the normalized magnetization of a cluster with size S.The quantity M c (S) is a measure of parallel allignment of the dipole vectors within the individual clusters.
) and (f) are enhanced compared to those in Figs.10(b) and (c).

FIG. 11 .
FIG. 11. (Color online) Mean cluster size S and mean cluster magnetization M as a function of the shift at two temperatures T * = 1.0 ((a),(b)) and T * = 1.35 ((c),(d)).Colors are the same as in Fig. 9.