Nonlinear Excitations in Magnetic Lattices with Long-Range Interactions
Miguel Molerón, C. Chong, Alejandro J. Martínez, Mason A. Porter, P. G. Kevrekidis, Chiara Daraio
NNonlinear Excitations in Magnetic Lattices with Long-Range Interactions
Miguel Moler´on, C. Chong ∗ , Alejandro J. Mart´ınez, Mason A. Porter, P. G. Kevrekidis, and Chiara Daraio Institute of Geophysics, Department of Earth Sciences, ETH Zurich, 8092 Zurich, Switzerland Department of Mathematics, Bowdoin College, Brunswick, Maine 04011, USA Oxford Centre for Industrial and Applied Mathematics,Mathematical Institute, University of Oxford, Oxford OX2 6GG, UK Department of Mathematics, University of California, Los Angeles, CA 90095, USA Department of Mathematics and Statistics, University of Massachusetts, Amherst, MA, 01003, USA Division of Engineering and Applied Science California Institute of Technology Pasadena, CA 91125, USA (Dated: January 30, 2018)We study — experimentally, theoretically, and numerically — nonlinear excitations in lattices ofmagnets with long-range interactions. We examine breather solutions, which are spatially localizedand periodic in time, in a chain with algebraically-decaying interactions. It was established twodecades ago [S. Flach,
Phys. Rev. E , R4116 (1998)] that lattices with long-range interactionscan have breather solutions in which the spatial decay of the tails has a crossover from exponential toalgebraic decay. In this Letter, we revisit this problem in the setting of a chain of repelling magnetswith a mass defect and verify, both numerically and experimentally, the existence of breathers withsuch a crossover. Introduction.
There has been considerable progress inunderstanding localization in nonlinear lattices over thepast three decades [1]. A prototypical example are spa-tially localized and temporally periodic discrete breathers(or just “breathers”) [2]. The span of systems in whichbreathers have been studied is broad and diverse; they in-clude optical waveguide arrays and photorefractive crys-tals [3], micromechanical cantilever arrays [4], Josephson-junction ladders [5, 6], layered antiferromagnetic crys-tals [7, 8], halide-bridged transition metal complexes [9],dynamical models of the DNA double strand [10], Bose–Einstein condensates in optical lattices [11], and manyothers. Many of these studies concern models withcoupling between elements only in the form of nearest-neighbor interactions. However, there has been a greatdeal of theoretical and computational work in latticeswith interactions beyond nearest neighbors. For exam-ple, some models of polymers [12], quantum systems [13],and optical waveguide arrays [14, 15] have included inter-actions beyond nearest neighbors; see also [16, 17]. Dy-namical lattices with long-range interactions (e.g., withall-to-all coupling) have been used as models for energyand charge transport in biological molecules [18], andstudies of such long-range models have explored phenom-ena such as equilibrium relaxation [19], thermostatistics[20], chaos [21, 22], and energy thresholds [23, 24]. Os-cillators of numerous varieties have also been coupled vialong-range interactions on lattices (and more general net-work structures) [25, 26]. In fact, until recently, theywere often assumed to be a fundamental ingredient forthe formation of so-called “chimera states” [27–29].Long-range interactions can have a significant effecton nonlinear excitations and yield phenomena that arerather different from those that result from only nearest-neighbor coupling. For example, stationary solitarywaves with a nontrivial phase can arise both in dis-crete nonlinear Schr¨odinger (DNLS) equations with next-nearest-neighbor (NNN) interactions [16, 30] and in NNN discrete Klein–Gordon (KG) [31] equations, and bistabil-ity of solitary waves is possible in DNLS equations withlong-range interactions [32, 33]. Finally, most relevantfor the present paper, breathers in KG and Fermi–Pasta–Ulam–Tsingou (FPUT) lattices with long-range interac-tions can exhibit a crossover from exponential decay (atshort distances from the breather center) to algebraic de-cay (at long distances) if the interactions decay signifi-cantly slowly (specifically, algebraically slowly) [24]. Avariety of new studies continue to elucidate fascinatingconsequences of long-range interactions. For example,recent studies have revealed the emergence of travelingdiscrete breathers without tails in nonlinear lattices withsuitable long-range interactions [34] and the emergenceof a linear spectral gap, which enables the emergence ofa low-frequency breather [35], in nonlinear lattices withother long-range interactions. Although there are manytheoretical and computational studies of lattice systemswith long-range interactions, we are not aware of any ex-perimental realizations of breathers in such systems.In this paper, we use experiments, theory, and nu-merical computations to study a strongly nonlinear lat-tice with long-range interactions that decay algebraically.Specifically, we consider a one-dimensional (1D) chain ofrepelling magnets with a single mass defect. This sys-tem allows us to realize fundamental structures, such assolitary waves, in a tabletop setup with real-time spatio-temporal resolution [36, 37]. Moreover, the use of mag-netic interactions allows exciting applications. They havealready been used as a passive mechanism to couple nodesof a lattice for unidirectional wave-guiding [38]; and ithas been suggested that magnetic interactions can beused to design novel devices for frequency conversion [39]and shock absorption [36]. In our study, we focus onbreathers in a magnetic chain and demonstrate that thereis a crossover from exponential decay to algebraic decayin the spatial profile of these breathers. Our numericalfindings are consistent with theoretical predictions that a r X i v : . [ n li n . PS ] J a n were developed almost twenty years ago in [24], and theyagree quantitatively with our experiments. (a) n = 12 FIG. 1. (a) Picture of our experimental setup. The latticeconsists of 25 magnetic particles deposited on an air-bearingtable. The right boundary ( n = 12) is fixed, and the leftboundary ( n = −
12) is driven harmonically with an electro-dynamic transducer. The magnetic particles are composed ofa disc magnet (type Supermagnete S-03-01-N, with magneti-zation grade N48, a diameter of 3 mm, and a height of 1 mm).The inset shows a magnified view of the magnetic particlesembedded in a 3D-printed support: (left) normal particle and(right) defect particle. (b) Relationship between the force F and center-to-center separation distance d between two par-ticles. The plus signs represent experimental data, and thesolid curves represents a curve of the form F = Ad p . In theinset, we show a plot of log( F ) versus log( d ) that we use forfitting the exponent p and the magnetic coefficient A . Experimental setup.
In Fig. 1(a), we show a pictureof our experimental setup. We situate an array of discmagnets over a 150 mm ×
300 mm rectangular air-bearingtable from IBS Precision Engineering (to reduce surfacefriction) and between two Teflon rectangular rods (to re-strict the particle motion to one dimension). As shownin the inset of Fig. 1(a), we insert each magnet into a3D-printed support. We glue a glass slide below the 3D-printed support to obtain a desired amount of levitation.The magnets are axially magnetized, and they have thesame orientation, so each magnet repels its neighbors.The average mass of the non-defect particles in the 25-particle chain is M = 0 .
45 g (with a standard devia-tion of s = 0 . m = 0 .
20 g. To excite the chain harmonically, weglue the left boundary to an aluminum bar attached toan electrodynamic transducer (Beyma 5MP60/N). The measured total harmonic distortion of this transducer isbelow 10% in the amplitude range (between 0 and 4 cm)under consideration.We measure the motion using a digital image correla-tion (DIC) software from Correlated Solutions (VIC 2D).We use a camera (of model GS3-U3-41C6C-C from PointGray) to record the particles’ motion at a frame rate of200 fps. To help track the particles, we glue speckle pat-terns to the top of the 3D-printed support [see Fig. 1(a)].We postprocess the video files with the VIC softwareto extract particle displacements and velocities. As in[36, 38], we assume that the relationship between the re-pelling force and distance has the form F = Ad p , where F is the force and d is the center-to-center separation dis-tance between two particles. We estimate the magneticcoefficient A and exponent p by measuring the repellingforce at 22 separation distances [represented by plus signsin Fig. 1(b)]. We measure the repelling force by fixingone magnet to a load cell (of type OMEGA LCL-113G)and approaching another magnet using a high-precisiontranslation stage. Using a least squares fitting routine forlog( F ) versus log( d ) with our experimental data [see theinset in Fig. 1(b)] yields A ≈ . × − N/m p and p ≈ − . Theoretical Setup.
Our experimental setup motivatesthe following model (which assumes that each node, rep-resenting a magnet, is coupled to every node in a chain): M n ¨ u n = (cid:88) j =1 A ( jδ + u n − u n − j ) p (1) − A ( jδ + u n + j − u n ) p − η ˙ u n , where u n = u n ( t ) ∈ R is the displacement of the n thmagnet from its equilibrium position, the mass of the n th magnet is M n , the magnetic coefficient is A , and thenonlinearity exponent is p (Fig. 1(b) shows the spatialdecay in the force with respect to the center-to-centerdistance between particles). This model assumes thateach magnet, including its magnetic properties, is iden-tical. The equilibrium separation distance between twoadjacent magnets in an infinite lattice is δ . In a finitelattice, the equilibrium separation distance will dependon the lattice location, see the Supplementary Materialfor details. We model damping effects with a dashpotterm η ˙ u n , where we empirically estimate the dampingfactor η (see our discussion below). We apply a harmonicboundary drive u left ( t ) = a sin(2 πf b t ), where a denotesthe drive amplitude and f b denotes its frequency. Ourinitial theoretical considerations involve a Hamiltonianlattice, so we take a = η = 0. Later, when we compareour numerical results to experiments, we also considernonzero values of the drive amplitude and damping fac-tor.In a homogeneous chain (where all masses are identi-cal, so M n = M ) the linearization of (1) has plane-wave Frequency [Hz]2 3 4 5 6 PS D -4 -2 index F r equen cy [ H z ] P o w e r S pe c t r a l D en s i t y FIG. 2. Experimental power spectral density (PSD) for a ho-mogeneous chain of 25 magnets (dashed red curve) and for amagnetic chain with a mass defect located at site n d = − n d = − η = a = 0). The blue disk in theinset represents the numerical cutoff frequency, and the reddiamond shows the numerical defect mode. solutions u n = exp( ikn + iωt ), where ω ( k ) = K ∞ (cid:88) j =1 j s [1 − cos( jk )] (2)= K [ ζ ( s ) − Re { e ik φ ( e ik , s, } ] , where s = 1 − p , the linear stiffness is K = − Apδ p − /M , the Riemann zeta function is ζ ( s ), and φ ( z, s, a ) is the Hurwitz–Lerch transcendent function[40]. This dispersion curve is nonanalytic in thewavenumber k , because its κ th derivative (where κ isthe integer satisfying s − ≤ κ < s ) with respect to k is discontinuous at k = 0. Below we discuss the conse-quences of this nonanalyticity. The dispersion curve isanalytic at the upper band edge (i.e., at k = π ).Because we are interested in solutions that decay spa-tially to 0 at infinity, it is natural to seek breather fre-quencies that lie above the edge of the spectrum ω ( π )(to avoid resonances with linear modes). Equation (1)with M n = M is not an appropriate model for seekingsmall-amplitude (bright) breather solutions, because oneneeds the plane waves to have a modulational instability,which is not possible in a homogeneous magnetic chain[2]. Hence, to obtain breathers, we break the uniformityof the chain by introducing a light-mass defect, motivatedby the analysis of [41] for nonlinear lattices with nearest-neighbor interactions. This creates a defect mode thatlies above the edge of the linear spectrum, from whichbreathers can bifurcate. Breathers in nearest-neighborFPUT-like lattices with defects have been studied ex-tensively both theoretically [41] and experimentally [42].To find breathers in a magnetic chain, one can alter- -20 0 20 n -15 -10 -5 ˙ u n [ m / s ] (a) Frequency [Hz] n c (b) || FIG. 3. (a) Semi-log plot of a breather solution (black curvewith markers), with a frequency of f b ≈ .
54 Hz, of Eq. (1)with η = a = 0 for a magnetic chain with a defect particlein the center ( n d = 0). The vertical axis gives the absolutevalue of the velocity, and the horizontal axis gives the nodeindex. For comparison, we show a breather solution of thesame frequency for a lattice with only nearest-neighbor inter-actions (red dash-dotted curve). The vertical dashed line isthe predicted value of the crossover value n c from Eq. (5).(b) Numerically computed crossover point (black markers)and prediction based on Eq. (5) (curve). natively use a lattice with spatial heterogeneity (e.g., adimer) [43–45] or one with an on-site potential [46, 47]or local resonators [48, 49].A chain with a single mass defect is the starting pointfor our model with long-range interactions. We reducethe mass of the n d th node (but without modifying itsmagnetic properties) by adjusting the support in whichthe magnet is embedded [see Fig. 1(b)]. Consequently, M n = (cid:26) m , n = n d M , otherwise , (3)where n d is the index of the mass defect with mass m We start by numerically comput-ing time-periodic solutions of the Hamiltonian variant ofEq. (1) (i.e., with a = η = 0) and N = 65 nodes. Thevalues that we use for the magnetic potential parametersare A ≈ . × − N/m p and p ≈ − . M = 0 . 45 g; the mass of the defect node is m = 0 . 20 g.We numerically compute the linear spectrum and obtaina defect mode with frequency f d ≈ . n c , and the decay becomes algebraic ratherthan exponential. This feature was first observed abouttwo decades ago in a KG lattice with a cubic potential(i.e., in the φ model) [24], which has long-range interac-tions with coefficients with algebraic decay (in particu-lar, they have a power-law decay O (1 /n s ) with respect tonode n ). The linearization of Eq. (1) also has interactioncoefficients with power-law decay O (1 /n s ).The algebraic decay of the breather far away from itscenter arises as follows. Its amplitude is small away fromits center, so we can linearize the equations of motion.Additionally, because the breather is temporally peri-odic, we can express the time dependence of the solu-tion as a Fourier series u n ( t ) = (cid:80) j ˆ u n ( j ) e ijω b t , where ω b = 2 πf b is the breather’s angular frequency. One com-putes the Fourier coefficients using Green’s functions [24]to obtain ˆ u n ( j ) = (cid:90) π cos( kj )( jω b ) − ω ( k ) dk , (4)where ω ( k ) is given by the dispersion relation in Eq. (2).Now it is clear why it is important to highlight thenonanalytic nature of ω ( k ): the Fourier coefficients inEq. (4) with discontinuities in the κ th derivative yieldFourier series that converge algebraically. This impliesthat u n ∼ /n s for large n [24]. One can make simi-lar arguments to explain the exponential decay near thecenter. (See [24] for details.)Assuming that the proportionality constants of the ex-ponential decay and the algebraic decay are roughly thesame, there is a crossover point between the two types ofdecay that satisfies e νn c = n sc , where ν is the exponentialdecay rate of the breather near the center. This yieldsthe following prediction for the crossover site n c [24]:log n c n c = ν − p . (5)For the solution in Fig. 3(a), the predicted crossover is n c = 10, which is roughly where the decay propertieschange in the numerical solution [see Fig. 3(a)]. Be-cause we made several assumptions to derive Eq. (5), wealso compute the crossover point from the numerically-obtained breather solutions. We calculate this point nu-merically by determining the first node at which the devi-ation of the solution from the best-fit line in the semi-logscale exceeds 1% of the solution amplitude. In the exam-ple in Fig. 3(a), this yields a crossover point of n c = 12.Equation (5) predicts that the crossover location dependson the solution’s exponential decay rate ν , which in turndepends on the breather frequency f b . In Fig. 3(b), weshow a comparison of observed numerical crossovers andEq. (5) for various breather frequencies. Experimental Results. We now turn our attention tothe experimental realization of breathers in a nonlinearlattice with long-range interactions. For our experiments,we consider a chain of N = 25 magnets (including the boundaries) with a defect magnet at site n d = − 8. Weexperimentally probe the linear spectrum by performinga frequency sweep. To do this, we excite the chain at 33frequencies between 2 and 6 Hz and extract the resultingsteady-state displacement amplitudes at the excitationfrequency in different locations. The red dashed line inFig. 2 represents the power spectral density (PSD) ofparticles − η = a = 0) of Eq. (1) (which wascomputed numerically, as shown in the inset of Fig. 2)agrees with the experimentally-observed passband cutofffrequency f ≈ . 50 Hz and defect-mode frequency f d ≈ . 66 Hz.To further evaluate our model, we initialize the exper-imental chain using the displacements that correspondto the theoretically-predicted Hamiltonian breather withfrequency f b ≈ . 46 Hz. The nodes oscillate initially withthe predicted frequency [see Fig. 4(a)]. In this particularexperiment, we do not add energy to the system. Thus,as the oscillation amplitude decreases due to damping,the dynamics gradually becomes more linear and the os-cillation frequency approaches the sole linear defect-modefrequency f d ≈ . 66 Hz. We use this experiment to em-pirically determine the damping parameter η = 0 . 10 g/sto match the temporal amplitude decay of the defect par-ticle. [See the inset in Fig. 4(a).] We conduct an anal-ogous numerical experiment using Eq. (1) with dampingbut no driving (specifically, η = 0 . 10 g/s and a = 0),which matches the observed experimental data; see thesolid red disks in Fig. 4(a).Our final experiment probes the decay properties ofthe breather. To allow the experimental system to reacha steady state (which allows us to more closely exam-ine these properties), we again harmonically excite theleft boundary magnet, so the displacement of the bound-ary magnet is u left = a sin(2 πf b ). We thereby treat theboundary as a “core” of the breather, so we do not usea defect particle in these experiments. We seek time-periodic solutions of Eq. (1) that account for both theboundary excitation and damping effects. We use theparameter values η = 0 . 10 g/s and a = 3 . . 05 mm/s, below theamount of noise in the experiments. This value corre-sponds to the mean velocity amplitudes of particles 9–24,whose motion can be attributed primarily to ambient vi-brations. Thus, for the drive (breather) frequency f b = 6Hz, we observe only exponential decay.However, for a drive frequency of f b = 11 Hz, the tran-sition to algebraic decay occurs close to the core of thebreather, so there appears to be a glimpse of the asso-ciated decay prior to reaching the level at which ambi-ent noise vibrations overwhelm the algebraic tail. Notethat the crossover approaches the core of the breatheras the breather frequency increases [see Fig. 3(b)]. InFigs. 4(b,c), we show the tails of the breather in semi-log and log-log plots. For f b = 6 Hz, the experimen- u d e f ec t [ mm ] [s] M ea n n -8 -6 -4 -2 l og | ˙ u n | semi-log scale f = 11 Hz (experiment)f = 11 Hz (model)best-fit power lawf = 6 Hz (experiment)f = 6 Hz (model)best-fit exponential (b) log n -6 -4 -2 l og | ˙ u n | log-log scale (c) FIG. 4. (a) Experiment initialized with a Hamiltonian breather solution of Eq. (1) with frequency f b ≈ . 46 Hz. We showthe mean oscillation frequency of the defect particle for every 1 . 28 s for the experiment (blue markers with error bars) andmodel with damping (with η = 0 . 10 g/s) but no driving (red disks). The error bars indicate the standard deviation over 5experimental realizations. Note that the node oscillates initially at the predicted frequency. The frequency approaches the soledefect frequency of the linear system, as the damping cause displacements to approach 0. In the inset, we show an example ofdefect-particle decay from an experiment. (b) Semi-log plot of the experimental data for drive frequencies of f b = 6 Hz (openred circles with error bars) and f b = 11 Hz (open blue squares with error bars). The chain is homogeneous (there is no defectparticle), because the boundary drive is acting like the defect particle (which we label as n = 0). We show our predictionsfrom the damped, driven model (filled markers) as well as the best fit to exponential (yellow curve) and power-law (blue curve)decay. The experimental data for f b = 6 Hz follows a roughly linear trend in the semi-log plot, suggesting that its decay isexponential. (c) Same as panel (b) but as a log-log plot. The experimental data for f b = 11 Hz follows a roughly linear trendin the log-log plot, suggesting that its decay is algebraic. Panels (b) and (c) share the legend that we show in (b). The insetin panel (c) shows a similar result for a chain of length N = 29 (which has a smaller equilibrium distance). In this case, morenodes have an amplitude that is comparable to the amount of noise. tal data (open red circles with error bars) has a roughlylinear trend in the semi-log plot, suggesting that its de-cay is exponential. The experimental data follows themodel prediction (solid yellow circles) up to the pointat which it reaches the noise level (the horizontal graydashed line). We fit (using a least squares procedure)the model solution with an exponential curve of the form αe − βn (solid yellow curve), and we obtain α ≈ . β ≈ . f b = 11 Hz, the experimental data(open blue squares with error bars) has a roughly lineartrend in a log-log plot, suggesting its algebraic decay. Theexperimental data follows our model’s prediction (closedlight blue squares) until reaching the noise level (hori-zontal gray dashed line). We fit the model solution witha power-law curve of the form αn − β (solid blue curve),and we obtain α ≈ . 579 and β ≈ . f b = 9 Hz (red) and f b = 13 Hz (blue) for achain with N = 29 nodes. Because the lattice is confinedto a length of L ≈ . / N = 25 chain. This increasesthe linear stiffness and hence increases the passband cut-off. Consequently, we need higher frequencies to avoidresonance with the linear modes. Discussion and Conclusions. We studied a lattice ofmagnets with long-range interactions, and we obtainedquantitative agreement between theory, numerics, andexperiment. Specifically, using a combination of experi-ments, computation, and analysis, we explored the pre- diction of [24], made about twenty years ago, that thetail of a breather solution of this nonlinear lattice ex-hibits a transition from exponential to algebraic decay.As far as we are aware, our work represents the first ex-perimental realization of breathers in a nonlinear latticewith long-range interactions.The study of long-range interaction systems is an in-creasingly important topic in numerous and wide-rangingareas of physics. These include dipolar Bose–Einsteincondensates (BECs) [50], where the recent formation ofquantum droplets and their bound states [51] suggeststhat interesting types of long-range interactions can alsoarise in the study of BECs in optical lattices. Long-range interactions also play important roles in the studyof coupled phase oscillators in diverse physical settings[26], heat transport in oscillator chains coupled to ther-mal reservoirs [52, 53], and more. Moreover, our experi-mental setup of magnet lattices has the potential to en-able systematic, well-controlled studies of phenomena in-volving these mixed (exponential and algebraic) decayingbreathers. For instance, it would be especially interestingto examine what happens when such breathers interactand how the decay properties (and interactions betweenbreathers) depending on lattice dimensionality. Acknowledgements. This material is based upon worksupported by the National Science Foundation underGrant No. DMS-1615037 (awarded to CC) and EFRIGrant No. 1741565 (awarded to CD). AJM acknowledgessupport from CONICYT (BCH72130485/2013). PGKgratefully acknowledges support from the US-AFOSR viaFA9550-17-1-0114. [1] P. G. Kevrekidis, IMA J. Appl. Math. , 389 (2011).[2] S. Flach and A. Gorbach, Phys. Rep. , 1 (2008).[3] F. Lederer, G. I. Stegeman, D. N. Christodoulides, G. As-santo, M. Segev, and Y. Silberberg, Phys. Rep. , 1(2008).[4] M. Sato, B. E. Hubbard, and A. J. Sievers, Rev. Mod.Phys. , 137 (2006).[5] P. Binder, D. Abraimov, A. V. Ustinov, S. Flach, andY. Zolotaryuk, Phys. Rev. Lett. , 745 (2000).[6] E. Tr´ıas, J. J. Mazo, and T. P. Orlando, Phys. Rev. Lett. , 741 (2000).[7] L. Q. English, M. Sato, and A. J. Sievers, Phys. Rev. B , 024403 (2003).[8] U. T. Schwarz, L. Q. English, and A. J. Sievers, Phys.Rev. Lett. , 223 (1999).[9] B. I. Swanson, J. A. Brozik, S. P. Love, G. F. Strouse,A. P. Shreve, A. R. Bishop, W.-Z. Wang, and M. I.Salkola, Phys. Rev. Lett. , 3288 (1999).[10] M. Peyrard, Nonlinearity , R1 (2004).[11] O. Morsch and M. Oberthaler, Rev. Mod. Phys. , 179(2006).[12] D. Hennig, The European Physical Journal B — Con-densed Matter and Complex Systems , 419 (2001).[13] A. G. Choudhury and A. Roy Chowdhury, PhysicaScripta , 129 (1996).[14] P. G. Kevrekidis, B. A. Malomed, A. Saxena, A. R.Bishop, and D. J. Frantzeskakis, Physica D , 87(2003).[15] N. K. Efremidis and D. N. Christodoulides, Phys. Rev.E , 056607 (2002).[16] P. G. Kevrekidis, Physics Letters A , 3688 (2009).[17] P. G. Kevrekidis, Journal of Optics , 044013 (2013).[18] S. F. Mingaleev, P. L. Christiansen, Y. B. Gaididei,M. Johansson, and K. Ø. Rasmussen, Journal of Bio-logical Physics , 41 (1999).[19] G. Miloshevich, J.-P. Nguenang, T. Dauxois, R. Khome-riki, and S. Ruffo, Phys. Rev. E , 032927 (2015).[20] H. Christodoulidi, C. Tsallis, and T. Bountis, EPL (Eu-rophysics Letters) , 40006 (2014).[21] G. M. Zaslavsky, M. Edelman, and V. E. Tarasov, Chaos , 043124 (2007).[22] N. Korabel and G. M. Zaslavsky, Physica A , 223(2007).[23] M. Kastner, Nonlinearity , 1923 (2004).[24] S. Flach, Phys. Rev. E , R4116 (1998).[25] M. A. Porter and J. P. Gleeson, in Frontiers in Ap-plied Dynamical Systems: Reviews and Tutorials , Vol. 4(Springer-Verlag, 2016).[26] A. Arenas, A. D´ıaz-Guilera, J. Kurths, Y. Moreno, andC. Zhou, Physics Reports , 93 (2008).[27] J. Xie, E. Knobloch, and H.-C. Kao, Phys. Rev. E ,022919 (2014).[28] J. Xie, E. Knobloch, and H.-C. Kao, Phys. Rev. E ,042921 (2015).[29] M. J. Panaggio and D. M. Abrams, Nonlinearity , R67(2015).[30] C. Chong, R. Carretero-Gon´azlez, B. A. Malomed, andP. G. Kevrekidis, Physica D: Nonlinear Phenomena ,1205 (2011).[31] V. Koukouloyannis, P. G. Kevrekidis, J. Cuevas, andV. Rothos, Physica D: Nonlinear Phenomena , 16 (2013).[32] Y. B. Gaididei, S. F. Mingaleev, P. L. Christiansen, andK. Ø. Rasmussen, Phys. Rev. E , 6141 (1997).[33] K. Ø. Rasmussen, P. Christiansen, M. Johansson, Y. Gai-didei, and S. Mingaleev, Physica D: Nonlinear Phenom-ena , 134 (1998).[34] Y. Doi and K. Yoshimura, Phys. Rev. Lett. , 014101(2016).[35] Y. Yamaguchi and Y. Doi, arXiv:1710.08675 (2017).[36] M. Moler´on, A. Leonard, and C. Daraio, Journal of Ap-plied Physics , 184901 (2014).[37] A. Mehrem, N. Jim´enez, L. J. Salmer´on-Contreras,X. Garc´ıa-Andr´es, L. M. Garc´ıa-Raffi, R. Pic´o, and V. J.S´anchez-Morcillo, Phys. Rev. E , 012208 (2017).[38] N. Nadkarni, A. F. Arrieta, C. Chong, D. M. Kochmann,and C. Daraio, Phys. Rev. Lett. , 244501 (2016).[39] M. Serra-Garcia, M. Moler´on, and C. Daraio,arXiv:1704.07226v1 (2017).[40] National Institute of Standards and Technology (2015),available at http://dlmf.nist.gov/ ; release 1.0.10(2015-08-07).[41] G. Theocharis, M. Kavousanakis, P. G. Kevrekidis,C. Daraio, M. A. Porter, and I. G. Kevrekidis, Phys.Rev. E , 066601 (2009).[42] N. Boechler, G. Theocharis, and C. Daraio, Nat. Mater. , 665 (2011).[43] G. Theocharis, N. Boechler, P. G. Kevrekidis, S. Job,M. A. Porter, and C. Daraio, Phys. Rev. E , 056604(2010).[44] N. Boechler, G. Theocharis, S. Job, P. G. Kevrekidis,M. A. Porter, and C. Daraio, Phys. Rev. Lett. ,244302 (2010).[45] G. Huang and B. Hu, Phys. Rev. B , 5746 (1998).[46] G. James, Math. Models Meth. Appl. Sci. , 2335(2011).[47] G. James, P. G. Kevrekidis, and J. Cuevas, Physica D , 39 (2013).[48] L. Liu, G. James, P. G. Kevrekidis, and A. Vainchtein,Nonlinearity , 3496 (2016).[49] L. Liu, G. James, P. G. Kevrekidis, and A. Vainchtein,Physica D , 27 (2016).[50] T. Lahaye, C. Menotti, L. Santos, M. Lewenstein, andT. Pfau, Reports on Progress in Physics , 126401(2009).[51] I. Ferrier-Barbut, H. Kadau, M. Schmitt, M. Wenzel,and T. Pfau, Phys. Rev. Lett. , 215301 (2016).[52] C. Olivares and C. Anteneodo, Phys. Rev. E , 042117(2016).[53] S. Iubini, P. D. Cintio, S. Lepri, R. Livi, and L. Casetti,arXiv:1712.07979 (2017). Supplementary Material I. EQUATIONS OF MOTION FOR FINITE CHAIN n-30 -20 -10 0 10 20 30 u n [ mm ] -0.4-0.200.20.40.60.8 (a) n-30 -20 -10 0 10 20 30 ˙ u n [ m / s ] -0.2-0.100.10.20.30.4 (b) FIG. 5. Breather solution of Eq. (7) with η = a = 0. We show the displacements of the solution in (a) and their velocities in(b). For experimental and numerical approaches, we consider a chain of N (where N is odd) magnets that we arrangeas a lattice confined within a distance L ∈ R with fixed boundary conditions (i.e., u − N +12 = u N +12 = 0). Under theseconditions, the equilibrium distance between magnets n − n depends on n . The N + 1 equilibrium distances δ ,n (with n ∈ {− N − , . . . , N +12 } ) satisfy L = N +12 (cid:88) n = − N − δ ,n and the following N equations:0 = n − (cid:88) j = − N +12 n (cid:88) i = j +1 δ ,i p − N +12 (cid:88) j = n +1 (cid:32) j (cid:88) i = n +1 δ ,i (cid:33) p , n ∈ {− N − , . . . , N − } . (6)Thus, for a finite chain, we obtain the following N equations of motion: M n ¨ u n = n − (cid:88) j = − N +12 A n (cid:88) i = j +1 [ δ ,i ] + u n − u j p − N +12 (cid:88) j = n +1 A (cid:32) j (cid:88) i = n +1 [ δ ,i ] + u j − u n (cid:33) p − η ˙ u n , n ∈ {− N − , . . . , N − } , (7) u − N +12 ( t ) = a sin(2 πf b t ) ,u N +12 ( t ) = 0 . For an infinite lattice (i.e. in the limit N → ∞ ) the equilibrium distances are constant with respect to lattice site.This is easily verified by substituting δ ,n = δ into Eq. (6), n − (cid:88) j = −∞ (( n − j ) δ ) p − ∞ (cid:88) j = n +1 (( j − n ) δ ) p = (8) ∞ (cid:88) j =1 − n (( j + n )) p − ∞ (cid:88) j = n +1 (( j − n )) p = (9) ∞ (cid:88) k =1 ( k ) p − ∞ (cid:88) (cid:96) =1 ( (cid:96) ) p = 0 (10)where new indices where defined k = j + n and (cid:96) = j − n . Substituting δ ,n = δ into Eq. (7) and redefining indicesonce again leads to Eq. (1) in the main text, which is valid for an infinite lattice.We find time-periodic solutions of Eq. (7) with period T by numerically computing roots x of the map f ( x ) = x − ˜ x ( T ), where x is the initial value of Eq. (7) and ˜ x ( T ) is the solution at time T of Eq. (7) with initial value x . See [2] for details. We numerically integrate Eq. (7) with an adaptive-size Runge–Kutta method. We use thelinearization of (7) to determine our initial guess for the Newton iterations. We show a numerical solution with abreather frequency f b ≈ ..