Table Of Content

A Simple Model for Long-Range Interacting Pendula Owen Myers, Adrian Del Maestro, and Junru Wu Materials Science Program, University of Vermont, Burlington, Vermont 05405, USA. and Department of Physics, University of Vermont, Burlington, Vermont 05405, USA.∗ Jeffrey S. Marshall School of Engineering, University of Vermont, Burlington, Vermont 05405, USA. (Dated: January 30, 2015) 5 WeshowthattheHamiltonianmeanfield(HMF)modeldescribestheequilibriumbehavior 1 of a system of long pendula with flat bobs that are coupled through long-range interactions 0 2 (charged or self gravitating). We solve for the canonical partition function in the coordinate frame of the pendula angles. The Hamiltonian in the angles coordinate frame looks similar n a to the form of the HMF model but with the inclusion of an index dependent phase in the J interaction term. We also show interesting non-equilibrium behavior of the pendula angles, 9 namely that a quasistationary clustered state can exist when pendula angles are initially 2 ordered by their index. ] h p I. INTRODUCTION these characteristics compel cautious use of the - usualtoolsofstatisticalmechanics, theyarealso s s Systems with long-range interactions are a the source of many interesting dynamical and a l source of unique problems in the field of statisti- statistical features. Depending on the system of c cal mechanics and thermodynamics. This is due interest, such features include canonical and mi- . s c toseveralpropertiesoflong-rangesystemswhich crocanonical ensemble inequivalence and related i fall outside of the conditions normally needing negative specific heat [3], quasistationary states s y to be satisfied when applying the methodolo- (different than metastable states which lie on h gies of thermodynamics. Simply from the words local extrema of equilibrium potentials) whose p [ “long-range” the first infringement can be de- lifetimes increase with the number of particles duced, that long-range systems are not addi- [4], aninterestingdependenceofthelargestLya- 2 v tive. If two systems with short-range interac- punov exponent on particle number in a long- 6 tions are brought together to form a larger sys- range Fermi-Pasta-Ulam model [5], and sponta- 1 tem then the energy difference between the con- neous creation of macroscopic structures in non- 1 4 glomeratesystemandthesumofitsconstituents equilibrium states [6]. In some cases, long-range 0 is the new potential energy from the boundary interactions can greatly simplify problems. For . 1 between them. In the thermodynamic limit, the instance,meanfieldmodelsdependononeoftwo 0 potential energy of the boundary is small com- premises: (i)interactionsareshort-rangebutthe 5 1 pared to the bulk and can be neglected, mak- system is embedded in a space of infinite dimen- v: ing short-range systems additive. In the case of sion so that all bodies in the system are nearest i long-range interactions, one particle will feel a neighbors, or (ii) interactions are infinitely long. X significant potential created by every other par- For some time, the primary motivation r a ticle, so the additional potential energy of two for the study of long-range interactions was systems added together does not scale as the to understand galaxies, galaxy clusters and boundary but in a more complicated way that the general thermodynamic properties of self- dependsonthespecificnatureoftheinteractions gravitatingsystems. Asidefrommeanfieldmod- [1]. Directly related to the lack of additivity is els, interest has further built since the obser- the fact that systems with long-range interac- vation of modified scattering lengths in Bose- tions are not extensive because their energy di- Einstein Condensates (BEC) through the use of vergesinthethermodynamiclimit[2]. Although Feshbach resonances [7]. Using this technique, a BEC can be made to be almost non-interacting ∗ [email protected] by tuning the scattering length to zero. One 2 could even tune the scattering length to a neg- The index dependence in the Hamiltonian, that ative value, making the BEC collapse. More recently, O’dell et al. [8] has shown that it may be possible to produce an attractive 1/r potential between atoms in a BEC by applying an “ex- tremelyoffresonant”electromagneticfield. This has opened the possibility of creating table-top methods which physically model aspects of cos- mological behavior on a laboratory scale, as well as the possible development of entirely new dynamics in BEC. The challenges in understanding long-range systems drive the development of solvable mod- Figure 1. N pendulum system with parallel planes els that could help better explain some of the of rotation. The ith pendulum angle at some time t aforementioned phenomena. Campa et al. [9] is θ (t). i have recently published a collection of important solvable models. One particularly signifi- willbedescribedindetailinthenextsection,isa cant model, which is important to this work, is consequence of the pendula pivots being slightly the Hamiltonian Mean Field (HMF) XY spin offset from one another and appears as a phase model [6], often written in the form inthecosinetermoftheHMFmodel. Itinspires the investigation of non-equilibrium “repulsive” behaviorintheanglecoordinateframewherewe (cid:88)N p2 γ (cid:88)N find an interesting quasistationary state when H = i + [1−cos(θ −θ )], (1) i j the angles of the pendula are initially ordered 2 2N i=1 i,j=1 according to their indices. We find the clustered positionsintheusualHMFcoordinateframe(bi- where θ is the angular position of the ith parti- i clusters), butintheanglecoordinateframeclus- cle (spin), as shown in Fig. 1, and p is its con- i tering is only found for the initially ordered an- jugate angular momentum. The HMF model is gles and, unlike the biclusters, these are clearly generally used to describe two different classes quasistationary states. A quasistationary state of systems: 1) a mean field XY classical spin is defined as a dynamical state that can persist model, and2)aonedimensionalperiodicsystem for a length of time which goes to infinity as of itinerant particles with long-range interac- the thermodynamic limit is approached [9]. In tions. Though the connection between the HMF addition to discussing the clustered angle states modelandthesecondclassofsystemsmentioned exhibited by the system, we also solve for the could be thought of as contrived given the sim- canonicalpartitionfunctioninthependuluman- plifications under which the model is realized, it gle coordinate frame, finding that in equilibrium has been shown that the model produces useful withaheatbath,theprobabilitydistributionsof insights into how non-neutral plasmas and self theanglescanbedescribedbytheoriginalHMF gravitating systems behave [6]. model. This finding is similar to the work done In this paper, we study the dynamics of an by[10]onamodelsometimescalledtheHMFα- array of N pendula with long-range interacting model. IntheHMFαmodel,a1/rα dependence ij bobs. Byconsideringlongpendulawithflatbobs between the classical spins is introduced [9, 11– undergoingsmalloscillationsandhavingparallel 14], wherer isthedistancebetweentheith and ij planes of rotation, we produce a model related jth spins on a lattice. Though the physical mo- to the HMF model through a coordinate trans- tivations behind studying these various models formation. The transformation introduces a de- can be very different, it is interesting that their pendence on the indices of the particle labels. A equilibrium behavior is the same or nearly the cartoonofthephysicalpictureisshowninFig.1. same. We believe that the work in this paper 3 furthersuggeststhattheHMFmodeluniversally terms of φ , the positions can be rewritten as i describes an entire class of long-range interact- 2πi ing systems in equilibrium. x = +φ . (4) i i N B. Density Approximation II. THE MODEL We have not yet explicitly stated the phys- A. Coordinates ical mechanism through which the bobs interact. Connecting the interactions with specific In Fig. 1, we show an array of pendula ro- physical motivations should be discussed with tating in the same plane with bobs that inter- some discretion because the development of the act through a long-range potential. If we con- model leaves these motivations up to some free- siderthecasewhereallthependulaonlyundergo dom of interpretation. Imagine that the bobs all small oscillations, we may write the horizontal carry some charge. We will not distinguish be- location of the ith particle, x , as x = id+(cid:96)θ , tween particles in any other way than their in- i i i where d is the distance between the pivots of dices, so in the case where all particles carry the neighboring pendula and (cid:96) is the length of each same charge, repulsive behavior is expected. On pendulum. The small θ regime makes the prob- the other hand one could make the bobs attract lem one dimensional in x. We choose periodic one another, which could be thought of as the boundary conditions and rescale the system by self-gravitating case. To be solvable, the model 2π/Nd so that requires some simplifications. For the sake of brevity we will speak of the particle charge or 2π x → x (2) mass density as the “density”. Nd The approximation that we invoke is simi- making the total system length a dimensionless lar to that used when justifying the HMF model 2π where N is the number of particles in one (Eq. (1)) to describe free particles in a one- period. We will refer to a periodic space with dimensional ring [6, 15]. The distribution of the length 2π as a unit circle. The position of the bobs is such the mass density, ρ(x), is given by ith particle (bob) is now N (cid:88) 1 2πi 2π (cid:96) ρ(x) = δ(x−x )− . (5) i x = + θ . (3) 2π i i N N d i=1 For reasonable choices of (cid:96) and d ((cid:96)/d << N), The constant 1/2π subtracted from the delta the second term on the RHS is suitably small functionisnecessarytoproduceameaningfulex- such that the Hamiltonian can be written with pression for the potential Φ and corresponds to terms that are quadratic in θ. However, we are the inclusion of a neutralizing (of opposite sign) primarily interested in a regime where (cid:96)/d → homogeneous background density. Restricting ∞ as the thermodynamic limit is approached. the problem further to that of solving Poisson’s Physically this corresponds to the small oscil- equation for a one-dimensional potential phys- lations of very long pendula with suspension ically amounts to choosing large and flat bob points that are close together compared to their geometries oriented with their smallest axis par- lengths. In order to simplify the calculations allel to the x axis. Writing the delta function that follow, we define φ to be the last term as a cosine Fourier series, Poisson’s equation be- i on the RHS of Eq. (3), namely φ ≡ 2π(cid:96)θ /Nd. comes i i Given the choice of large (cid:96)/d, φ can take any i N ∞ value in the range [0,2π). This is only true be- ∇2Φ(x) = γ (cid:88)(cid:88)cos[n(x−x )]. (6) i cause (cid:96)/d is large, not because the θ s are. In π i i=1n=1 4 The parameter γ contains the particle (bob) ble, respectively. The force that the jth particle charge or mass and becomes the interaction experiences when φ and all φ are zero is given j i strength in the Hamiltonian. We can see that by the zeroth-order term in the Fourier series can- N (cid:26) (cid:20) (cid:21) (cid:27) celed the constant neutralizing background that γ (cid:88) 2π(j −i) −∇Φ(x ) = − sin +c . was superficially added. j π N 1 i=1 Themostimportantsimplificationinthispa- (8) per is truncating the sum of the Fourier coeffi- Thesum(cid:80) sin[2π(j −i)/N]equalszeroforany i cients used to represent the delta function after j, so c must be zero. Integrating once more to 1 the n = 1 coefficient. Antoni et al. defend the obtain the potential yields truncation by asserting that the “large scale col- lective properties” do not greatly change when γ (cid:88)N (cid:20) (cid:18) 2πi (cid:19)(cid:21) Φ(x) = c −cos x− −φ . (9) higherordertermsofthesum(includinginterac- 2 i π N tions at the smaller length scales) are included, i=1 and discuss the consequences of the approxima- To determine c we stipulate that if all φ = 0, 2 i tion in some detail [6]. The simplification also then Φ(0) = 0. Inserting Eq. (3) (or Eq. (4)) for warrants a brief discussion of the way that it x yields i could be physically interpreted. The truncation ofthesumisequivalenttosmearingouttheden- γ (cid:88)N (cid:20) (cid:18)2πi(cid:19)(cid:21) Φ(0) = c −cos . (10) sity of each particle over the system so that it is 2 π N peaked at its given location, x , but also having i=1 i a negative density peak on the opposite side of The sum over the cosine is zero, therefore c = 0 2 the unit circle. This could be thought of as dou- andwecannowwritethepotentialenergyofthe bling the number of particles and enforcing that jth particle as each particle has a negative partner that always remains on the opposing side of the unit circle. Afterthisdoubling,thenownebulousmassesare N (cid:20) (cid:21) γ (cid:88) 2π(j −i) dispersed such that a pair’s density is described Φ(xj) = − cos +φj −φi . π N by a cosine function with the positive peak cen- i=1 (11) tered at x . i C. Solving Poisson’s Equation D. The Hamiltonian Integrating Poisson’s equation once, we obtain: The Hamiltonian can be written as N γ (cid:88) H = H0+HI, (12) ∇Φ(x) = {sin(x−x )+c }. (7) i 1 π i=1 where H is the kinetic energy piece 0 In order to determine the constant c from the integration,thephysicalpictureshould1 beexam- H = (cid:88)N p2i , (13) 0 ined. A sensible requirement is that when all of 2 i=1 the bobs are hanging at their equilibrium positions,directlybelowtheirpivot(allφ = θ = 0), and i i the net force experienced by any bob is zero. N (cid:20) (cid:21) γ (cid:88) 2π(i−j) This is a valid requirement if the bobs are at- H = − cos +φ −φ I i j 2N N tractive or repulsive, the only difference being i,j=1 that the configuration would be unstable or sta- (14) 5 is the interaction piece, so librium with a heat bath by solving the partition function in the φ coordinate frame. In the i (cid:88)N p2 γ (cid:88)N (cid:20)2π(i−j) (cid:21) process of simplifying the Hamiltonian to solve H = i − cos +φ −φ . 2 2N N i j for the partition function, we will find expres- i=1 i,j=1 sions of the form cosφ and sinφ which we talk (15) i i about as the horizontal and vertical components The mass of the bobs has been set to unity, γ is of a magnetization m(cid:126) = (cosφ ,sinφ ). It could theinteractionstrength, afactorof1/2accounts i i i easily be stated that in the spin analogy, the for the double counting, and the 1/π coefficient φ are orientations of the spins, but we should in the potential energy has been absorbed into i make a more concrete connection between this γ. The factor of 1/N is a rescaling of the poten- idea and the original presentation of the model. tial energy that ensures that as the thermody- We would like to remind the reader that even namic limit is approached, the potential energy though the angles θ of the pendula are small, of the system does not diverge. The 1/N scal- i the long suspensions of the bobs ((cid:96)) allow φ to ing is known as the Kac prescription [16]. The i cover the entire system which, rescaled, has di- Kac serves to keep both the energy and entropy mensionless length 2π. The system is also peri- of a system proportional to the number of parti- odic, so the bobs can be thought of as moving cles in the system, an important prerequisite for on a unit circle where the position of the ith bob phase transitions [9]. is x = 2πi/N + φ . In order to think of φ i i i as the spin orientations, we start by considering each bob as living on its own individual unit cir- E. Relationship to the HMF model and the Spin Interpretation cle. An example of these unit circles is shown in Fig. 2, a visual aid to the following. Imagine Due to the simple bijective relationship be- stacking horizontal circles in the vertical direc- tween x and φ one can simply solve the equation and rotating each by an angle 2π/N from i i tions of motion for the HMF model and find the one below. The projection of these circles thedynamicsforφ viathecoordinatetransform onto the horizontal plane would be the system i x → φ . Previously it was mentioned that the viewed in x, i.e. the HMF model. If we twist i i HMF model is used to describe free particles on the stack so there is no rotation between adja- a ring with long-range repulsion or attraction, cent circles and then project onto the horizontal as well as describing a classical XY spin model. plane it creates the picture viewed in φ, where The θ played the role of either the position of the pivot points are all aligned. The reason for i the ith particle on the ring or the orientation of this artificial construction of stacked circles is theith spin. Therefore, itisinterestingtospecu- partly to pictorially depict the transformation late about the type of spin system the model de- between x and φ and partly to show how m(cid:126)i (as scribes in the φ picture. Thus far, the rescaled defined) is just the orientation of the ith spin i angle φ = 2π(cid:96)θ /Nd describes the distance of in the φ picture. Said differently, each circle in i i a pendulum bob from the point directly below the φ picture represents a spin with an orienta- its pivot, but it could also be interpreted as the tion in the horizontal plane determined by φi; orientation of spin. In the spin interpretation an infinite-range classical mean field spin model of Eq. (15), the potential energy of the ith and described by Eq. (15). jth spin pair depend on both their relative orientation as well as the difference between their indices. In the following discussion it will some- III. EQUILIBRIUM timesbeconvenienttospeakaboutφ inthespin i language. Inthissection,wesolveforthecanonicalpar- We will prove that in the φ picture, the sys- tition function, in the φ coordinate frame using tem in equilibrium with a heat bath is equiva- the Hamiltonian in Eq. (15) and show that, in lent to the HMF model (the x picture) in equi- equilibrium, the HMF model describes the an- 6 (cid:40) (cid:20) (cid:21) −γ (cid:88) 2π(i−j) H = cos [cosφ cosφ i i j 2N N i,j +sinφ sinφ ] i j (cid:41) (cid:20) (cid:21) 2π(i−j) −sin [sinφ cosφ −cosφ sinφ ] . i j i j N (16) The coefficients in the Hamiltonian cos[2π(i−j)/N] and sin[2π(i−j)/N] should be thought of as matrices with components C ij and S respectively. The Hamiltonian can now ij be written in the form γ (cid:88) (cosφ C cosφ +sinφ C sinφ i ij j i ij j 2N i,j −sinφ S cosφ +cosφ S sinφ ), (17) i ij j i ij j which is suggestive because it can be regarded as the matrix equation Figure 2. Example of a system of N = 8 particles when viewed as individual spins in the a) x coordinate frame, and b) the φ coordinate frame. a) In the   cosφ xcoordinateframethedirectionoftheith spingiven (cid:34) 1 tbeyrnthateivaenlgy,lexxiicoisueldxpbreestsheoduagshtxiof=as2πthi/eNpo+siφtiio.nAol-f HI = 2γN (cosφ1,cosφ2,...,cosφN)Ccos...φ2 ith particle on the unit circle, shown at the bottom cosφ N of the figure as the projection of all positions onto   the horizontal plane. The black circles on the rings sinφ1 in the figure mark the location of the pendulum piv- sinφ2 ots at 2πi/N in x. b) Twisting the column of rings +(sinφ1,sinφ2,...,sinφN)C ..   .  in a) such that the pivots are aligned transforms the sinφ system into the φ coordinate frame. In this picture, N i thedirectionoftheith spinisjustgivenbytheangle cosφ  1 φi. cosφ2 −(sinφ1,sinφ2,...,sinφN)S ..   .  cosφ N   sinφ 1 (cid:35) sinφ2 +(cosφ1,cosφ2,...,cosφN)S ..  . gles of long pendula with long-range interacting  .  bobs. In order to solve the configurational piece sinφ N of the partition function the Hamiltonian must (18) be modified. Using the cosine and sine sum and difference identities twice, the potential interac- Itishelpfultoconsidertheparticlespositionson tion piece of the Hamiltonian H can be written the unit circle with respect to their pivot (φ) as I as magnetizations. Defining 7 trix. The matrices C and S can be rewritten as m(cid:126)i ≡ (cosφi,sinφi), (19) C = U†DCU and S = U†DSU, where DC and DS arediagonalmatriceswithdiagonalelements and with mT = (m ,m ,...,m ) where µ 0,µ 1,µ N−1,µ thataretheeigenvaluesofC andS,respectively, µ holds the place of x or y, the Hamiltonian whichwedenoteasλC andλS. Fromhereonwe i i becomes label the indices i from 0 to N −1. It is worth pointing out that due to S being antisymmetric, γ HI = 2N(mTxCmx+mTyCmy U must be complex. Equation (22) becomes −mTSm +mTSm ). (20) y x x y γ (cid:16) H = mTU†DCUm +mTU†DCUm I 2N x x y y A closer examination of the structure of the (cid:17) −mTU†DSUm +mTU†DSUm . (22) coefficient matrices C and S indicates that they y x x y take the special form of circulant matrices, and We will move back to the index notation using thus can be simultaneously diagonalized by a the following relations: unitary matrix U. A circulant matrix has the form DC,S = λC,Sδ , (23) i ij   a a a ... a a 1 2 3 N−1 N  aN a1 a2 ... aN−2 aN−1 N aN−1 aN a1 ... aN−3 aN−2 X ≡ (cid:88)U mx, (24)  ... ... ... ... ... ... , (21) j k=1 jk k    a3 a4 a5 ... a1 a2  and a a a ... a a 2 3 4 N 1 N a special kind of Toeplitz matrix, where each Y ≡ (cid:88)U my., (25) j jk k subsequent row is a cyclic permutation of the k=1 row above or below it. Any matrix A with ele- where X is not to be confused with x. Using the ments a that can be written in terms of some ij Kronecker delta, we set all i = j since these are function f(i−j) is a circulant matrix. Because theonlynonzeroterms. TheHamiltonianisnow a circulant matrix is a normal matrix it can be given by diagonalized by a unitary matrix. We show that C and S are simultaneously diagonalizable by N−1 −γ (cid:88) (cid:16) showing that they commute, i.e. [C,S] = 0 H = (cid:107)X (cid:107)2λC +(cid:107)Y (cid:107)2λC i 2N j j j j where [C,S] = CS − SC. Starting with the j=0 second term, −SC = STCT = (CS)T which (cid:17) −Y∗X λS +X∗Y λS , (26) is found by arguing that C is symmetric since j j j j j j cosine is an even function and does not change with (cid:107)X(cid:107) = XX∗ and (cid:107)Y(cid:107) = YY∗. The in- under the exchange of i and j, whereas S is odd clusion of the eigenvalues λC and λS simplifies because sine is an odd function and does change the Hamiltonian further. We will now solve for sign under exchange of i and j. The commuta- λC and λS. Looking at the form of a circulant tion becomes [C,S] = CS + (CS)T. Also, an matrix shown in Eq. (21) reminds us that the odd function multiplied by an even function re- elements of a circulant matrix can be defined sults in an odd function so the entire matrix CS with a single label. We write the single labeled is odd. Therefore (CS)T = −CS bringing us to elements of the cosine and sine matrices respec- the final expression [C,S] = CS −CS = 0. We tively as have shown that C and S can be simultaneously diagonalized by U. The matrix U is known for 2πl c = cos , (27) circulant matrices and is called the Fourier Ma- l N 8 and and 2πl 1 (cid:18)2πik(cid:19) s = sin , (28) b ≡ √ sin . (35) l N ik N N where l = 0,1,2,...,N − 1. The eigenvalues, This was done to write the absolute squares λA, of a N × N circulant matrix A can be in Eq. (33) in terms of the squares of a and ik written in terms of the single label elements b . By noticing that a = a and b = ik 1k (N−1)k 1k a . The jth eigenvalue of A is known to be l √ −b we write the configurational partition λA = (cid:80)N−1a exp(2πilj/N), where i is −1 (N−1)k j l=0 l function as (not an index) and l = 0,1,2,...,N −1. There- fore, (cid:90) λC = N(cid:88)−1cos(cid:18)2πl(cid:19)e2πilj/N, (29) ZI = A dNφeβ2γ((cid:80)k[a1kmxk−b1kmyk])2 j N l=0 ×eβ2γ((cid:80)k[b1kmxk+a1kmyk])2, (36) and where β = 1/k T. B N−1 (cid:18) (cid:19) (cid:88) 2πl λS = sin e2πilj/N. (30) The Hubbard-Stratonovich transformation is j N l=0 now applied twice, once to each quadratic quan- tity in the partition function. The integration Writing cosine and sine in their exponential variables introduced through this transforma- forms gives tion are z and z with subscripts for first and 1 2 N−1 second quadratic quantities, respectively. After 1 (cid:88) (cid:104) (cid:105) λC = ei2πl(j+1)/N +ei2πl(j−1)/N , (31) after switching the order of integration, we find j 2 l=0 ZI = A (cid:90) ∞ dz1dz2e−(z12+z22)/2βγ(cid:89) 2πβγ −i N(cid:88)−1(cid:104) (cid:105) −∞ k λS = ei2πl(j+1)/N +ei2πl(j−1)/N , (cid:90) π j 2 × dφ e(z1a1k+z2b1k)cosφk+(z2a1k−z1b1k)sinφk. k l=0 −π (32) (37) The above representations of the eigenvalues show that C and S each have only two non-zero Theintegrationcanbeperformedusingtheiden- eigenvalues corresponding to j = 1,N −1 given tity by λC = λC = N/2 and λS = (λS )∗ = iN/21. The HNa−m1iltonian simplifi1es greatNly−1to (cid:90) π dφeξcosφ+ηsinφ = 2πI (cid:16)(cid:112)ξ2+η2(cid:17) (38) 0 −π where H = γ −(cid:0)(cid:107)X +iY (cid:107)2+(cid:107)X −iY (cid:107)2(cid:1). I 2 1 1 N−1 N−1 ξ2+η2 = (z a +z b )2+(z a −z b )2 (39) 1 1k 2 1k 2 1k 1 1k (33) which simplifies when a and b are included to The representation of H in Eq. (33) must I befurthermodifiedbeforethepartitionfunction (cid:20) 1 (cid:18)2πk(cid:19) 1 (cid:18)2πk(cid:19)(cid:21)2 z √ cos +z √ sin canbefound. WedothisbysplittingtheFourier 1 N N 2 N N matrix U into its real and imaginary compo- (cid:20) 1 (cid:18)2πk(cid:19) 1 (cid:18)2πk(cid:19)(cid:21)2 nents, a and b , given by + z √ cos −z √ sin ik ik 2 1 N N N N 1 (cid:18)2πik(cid:19) 1 a ≡ √ cos , (34) = (z2+z2) (40) ik N N N 1 2 9 It is convenient to make a change to polar coor- 1.0 (cid:112) dinates by introducing z = z2+z2, following ) 1 2 z which the partition function can be written as ( 00.8 I √ / A (cid:90) ∞ (cid:89) (cid:18) z (cid:19) z) Z = dze−z/2βγ 2πI √ . (41) (0.6 I βγ 0 N I1 −∞ k d0.4 Equation (41) is recognized to be an intermedi- n a ate step of the solution to the canonical parti- 0.2 β tion function for the HMF model. From here we / z jump to the main results, the details of which 0.0 0 1 2 3 4 5 are included in the HMF literature [6, 9, 15] . z The integration over z in Eq. (41) can be preformedusingthesaddlepointapproximation. Figure 3. The solid (black) curve is the fraction of The rescaled free energy per particle follows as modifiedBesselfunctionsI (z)/I (z),dashed(green) 1 0 β (cid:20)−z2 (cid:21) is z/β for β =1, solid (red) line is z/β for beta=2, −βF = − +inf +ln2πI (z) (42) dotted (blue) is z/β for β =4, all as a function of z. 0 2 z 2β (blue) line is z/β for β = 4. The values β = 1,2,4 correspond to the pre-phase transition, critical, and The expression above permits a convenient path post-phase transition values in that order. to finding the phase transition. Solving for the minimum values of z in order to satisfy the last term in Eq. (42) results in the equation ter is the total magnetization M(cid:126) = 1 (cid:80)N m(cid:126) N j=1 i z I1(z) where m(cid:126)i was defined to be (cosφi,sinφi). − = 0, (43) β I (z) Showing that the canonical partition func- 0 tion in the φ coordinate frame model and the which can be solved self consistently for z and HMF model are equivalent necessitates a more represented graphically for different values of β detailed discussion of the equilibrum behaviors as in Fig. 3. The reader will see that after β is in the φ frame. Campa et al. [9], in their re- increased passed the critical value (β = 2) there view of the HMF model, rigorously show ensem- are two well-defined solutions. ble equivalence between the canonical and mi- The Hubbard-Stratonovich transformation crocanoicalensembleoftheHMFmodel. Inlight decouplesspin-spin(squaredtermsintheHamil- of this fact, a large N microcanonical simulation tonian) contributions to the partition function should be able to produce equilibrium behavior at the price of needing to create a linear inter- like the phase transition mentioned above. The action between each spin with an auxiliary field temperature in a numerical simulation of a sys- z [17]. Again, a more detailed procedure can be tem with many particles would be “set” through found in [6, 9, 15] where discussion of the inter- a choice of the initial momenta distribution. In nal energy in the equilibrium state is followed this type of simulation, it is common practice by non-equilibrium behavior of the system pre- to compute the order parameter and free energy paredinmicrocanonicalensembles. Herewewill [1, 6, 18], begging the question: does a large mi- simply touch on the most important point of the crocaonical simulation of Eq. (15) approximate equilibrium behavior, being that for β < 2 the the expected equilibrium behavior? Also, since system is paramagnetic but for β ≥ 2 a pitch- the index-dependent model in equilibrium with fork bifurcation occurs resulting in two stable aheatbathcanbedescribedbytheHMFmodel, solutions. At this point there is a discontinuity would the dynamics of such a simulation quali- in the free energy, a second order phase transi- tatively resemble those in the HMF model? The tion occurs and the system can maintain finite answer to both of these questions is no if one magnetization. In this case, the order parame- were to find the equations of motions in φ for i 10 somelargeN andthencomparethemtoanHMF IV. NON-EQUILIBRIUM RESULTS model or the x coordinate frame. As stated, this discrepancy may appear to detract from our re- For a system of pendula, it is interesting to sult. Indeed, it uncovers a conceptual omission study an initial configuration where all pendula in the model, but it is one whose rectification aresettorandomsmalldisplacementsfromφ = i gives insight into the models ensemble equiva- 0. Specifically we initialize the ith pendulum lencepropertyofthemodel,orlackthereof. The angle, φ , randomly in the range [−π/N,π/N). i omission was in the arbitrary scaling of x which In x the indices are ordered in x such that we will now rectify. x < x < x < ... < x and the ith bob 1 2 3 N is randomly distributed in the range [2πi/N − π/N,2πi/N+π/N). It is possible to make some general statements about the dynamics of this configurationinxusingtheequationsofmotion. Expressing the Hamiltonian with terms that are We introduce the parameter L which gener- quadratic in φ yields alizes the scaling in Eq. (2) to 2πL x → x, (44) (cid:40) (cid:20) (cid:21) Nd γ (cid:88) 2π(i−j) H = sin (φ −φ ) I i j 2N N making the position of the ith particle ij (cid:32) (cid:33)(cid:41) (cid:20)2π(i−j)(cid:21) φ2 φ2 2πLi 2πL(cid:96) −cos 1− j − i +φ φ . xi = + θi. (45) N 2 2 i j N N d (46) and changing the definition of φ to φ ≡ i i 2πL(cid:96)θ /Nd. It can be shown that the intro- i Withthisexpression,theequationsofmotionfor duction of L only changes the final result of the ith particle can be written as the partition function by a constant factor of L due to the enlarged limits of integration. Nu- merically, we find is that if L >> 1, then the ∂H φ¨ = p˙ = , (47) i i simulations in φ closely reproduce the dynam- ∂φ i ics of HMF model simulations (dynamics in x). Therefore, for large L the mirocanonical simula- from which we obtain tions can approximate equilibrium and the an- swerstothepreviousquestions-doesalargemi- (cid:40) (cid:20) (cid:21) crocaonical simulation of Eq. (15) approximate −γ (cid:88) 2π(i−j) φ¨ = cos (φ −φ ) theexpectedequilibriumbehavior,andsincethe i 2N N j i j index-dependent model in equilibrium with a (cid:41) (cid:20) (cid:21) heat bath can be described by the HMF model, 2π(i−j) +sin . (48) would the dynamics of such a simulation quali- N tativelyresemblethoseintheHMFmodel? -becomes yes. Alternatively, the coordinate frame In the above equation, the last term and the inequivelenceismostextremeforsmallL. These cos[2π(i−j)/N]φ term sum to zero, leading to i numerical results were found using initial conditions that are randomly distributed φ about i (cid:20) (cid:21) −γ (cid:88) 2π(i−j) the domain [−Lπ,Lπ). It should be stated that φ¨ = cos φ . (49) i j 2N N for the rest of this paper we work with L = 1 j becuase we are inetersed in cases where the φ coordinate frame is markedly differnet than the Using the difference formula, we write the accel- x coordinate frame. eration as